3

I'm using the maptools package in R to merge some data to a shapefile.

When I open the original input shapefile in QGIS, QGIS automatically picks up on the CRS.

When I open the output shapefile from R, I need to select the CRS.

This happens even if I specify the CRS in my R code. Is there a solution to this?

# Load libraries ----------
library("maptools")

# Read shapefile ----------
LSOA <- readShapeSpatial("c:\\Temp\\LSOA11.shp")

# Set CRS ----------
# http://spatialreference.org/ref/epsg/27700/
OSGB36.CRS <- CRS("+proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 +y_0=-100000 +ellps=airy +datum=OSGB36 +units=m +no_defs")
proj4string(LSOA)<-OSGB36.CRS

# Check CRS ----------
proj4string(LSOA)

# Filter ----------
LSOA$Flag <- grepl("YORK",LSOA$LSOA11NM, ignore.case = TRUE)
LSOA <- LSOA[LSOA$Flag==TRUE,]

# Plot to check filtering ----------
plot(LSOA, col = 'light blue')

# Export shapefile ----------
writeSpatialShape(LSOA,fn="c:\\Temp\\LSOA_Fuel_Poverty")

1 Answer 1

5

Use the rgdal package. If the shapefiles have a projection defined rgdal::readOGR() is recommended. This package provides bindings to the GDAL library (Geospatial Data Abstraction Library) and access to projection/transformation operations from the PROJ.4 library. This supersedes the shapefile read/write functionality in maptools.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.