Projecting sp objects in R
You can use the spTransform()
methods in rgdal - using your example, you can transform the object to NAD83 for Kansas (26978):
library(rgdal)
library(maptools)
P4S.latlon <- CRS("+proj=longlat +datum=WGS84")
hrr.shp <- readShapePoly("HRR_Bdry", verbose=TRUE, proj4string=P4S.latlon)
plot(hrr.shp)
hrr.shp.2 <- spTransform(hrr.shp, CRS("+init=epsg:26978"))
plot(hrr.shp.2)
To save it in the new projection:
writePolyShape(hrr.shp.2, "HRR_Bdry_NAD83")
EDIT: Or, as per @Spacedman's suggestion (which writes a .prj file with the CRS info):
writeOGR(hrr.shp.2, dsn = getwd(), layer = "HRR_Bdry_NAD83", driver="ESRI Shapefile")
If one is not certain which CRS to project from, refer to the following post:
- Choosing the correct value for proj4string for shapefile reading in R?
And if one wants to define/assign a CRS when data doesn't have one, refer to:
- Assigning CRS to shapefile when it doesn't have one, in R?
Since the introduction of the sf-package (have a look at the vignettes sf1, sf2, sf3, sf4 and a migration guide here) you can use st_transform()
for re-reprojecting your vector data:
require(sf)
hrr_sf = st_read('HRR_Bdry.shp', stringsAsFactors = FALSE,
crs = 4326) # has +proj=longlat +datum=WGS84
plot(hrr_sf)
hrr_sf2 = st_transform(hrr_sf, "+init=epsg:26978") # 1st option sp::CRS() not working/ needed
hrr_sf2 = st_transform(hrr_sf, 26978) # 2nd option - EPSG code as an integer
plot(hrr_sf2)
# don't think about doing this:
hrr_sf3 = st_read('HRR_Bdry.shp', stringsAsFactors = FALSE,
crs = 26978)
# Output layer
st_write(hrr_sf2, dsn = getwd(), layer = "HRR_Bdry_NAD83", driver = "ESRI Shapefile")
sf will replace sp in the future and has due to its simplicity and speed imho several advantages compared to sp.