6.1 Introduction
Section 2.4 introduced coordinate reference systems (CRSs) and demonstrated their importance.This chapter goes further.It highlights issues that can arise when using inappropriate CRSs and how to transform data from one CRS to another.
As illustrated in Figure 2.1, there are two types of CRSs: geographic (‘lon/lat’, with units in degrees longitude and latitude) and projected (typically with units of meters from a datum).This has consequences.Many geometry operations in sf, for example, assume their inputs have a projected CRS, because the GEOS functions they are based on assume projected data.To deal with this issue sf provides the function st_is_longlat()
to check.In some cases the CRS is unknown, as shown below using the example of London introduced in Section 2.2:
london = data.frame(lon = -0.1, lat = 51.5) %>%
st_as_sf(coords = c("lon", "lat"))
st_is_longlat(london)
#> [1] NA
This shows that unless a CRS is manually specified or is loaded from a source that has CRS metadata, the CRS is NA
.A CRS can be added to sf
objects with st_set_crs()
as follows:24
london_geo = st_set_crs(london, 4326)
st_is_longlat(london_geo)
#> [1] TRUE
Datasets without a specified CRS can cause problems.An example is provided below, which creates a buffer of one unit around london
and london_geo
objects:
london_buff_no_crs = st_buffer(london, dist = 1)
london_buff = st_buffer(london_geo, dist = 1)
#> Warning in st_buffer.sfc(st_geometry(x), dist, nQuadSegs, endCapStyle =
#> endCapStyle, : st_buffer does not correctly buffer longitude/latitude data
#> dist is assumed to be in decimal degrees (arc_degrees).
Only the second operation generates a warning.The warning message is useful, telling us that the result may be of limited use because it is in units of latitude and longitude, rather than meters or some other suitable measure of distance assumed by st_buffer()
.The consequences of a failure to work on projected data are illustrated in Figure 6.1 (left panel):the buffer is elongated in the north-south direction because lines of longitude converge towards the Earth’s poles.
The distance between two lines of longitude, called meridians, is around 111 km at the equator (execute geosphere::distGeo(c(0, 0), c(1, 0))
to find the precise distance).This shrinks to zero at the poles.At the latitude of London, for example, meridians are less than 70 km apart (challenge: execute code that verifies this).Lines of latitude, by contrast, have constant distance from each other irrespective of latitude: they are always around 111 km apart, including at the equator and near the poles.This is illustrated in Figures 6.1 and 6.3.
Do not interpret the warning about the geographic (longitude/latitude
) CRS as “the CRS should not be set”: it almost always should be!It is better understood as a suggestion to reproject the data onto a projected CRS.This suggestion does not always need to be heeded: performing spatial and geometric operations makes little or no difference in some cases (e.g., spatial subsetting).But for operations involving distances such as buffering, the only way to ensure a good result is to create a projected copy of the data and run the operation on that.This is done in the code chunk below:
london_proj = data.frame(x = 530000, y = 180000) %>%
st_as_sf(coords = 1:2, crs = 27700)
The result is a new object that is identical to london
, but reprojected onto a suitable CRS (the British National Grid, which has an EPSG code of 27700 in this case) that has units of meters.We can verify that the CRS has changed using st_crs()
as follows (some of the output has been replaced by …
):
st_crs(london_proj)
#> Coordinate Reference System:
#> EPSG: 27700
#> proj4string: "+proj=tmerc +lat_0=49 +lon_0=-2 ... +units=m +no_defs"
Notable components of this CRS description include the EPSG code (EPSG: 27700
), the projection (transverse Mercator, +proj=tmerc
), the origin (+lat_0=49 +lon_0=-2
) and units (+units=m
).25The fact that the units of the CRS are meters (rather than degrees) tells us that this is a projected CRS: st_is_longlat(london_proj)
now returns FALSE
and geometry operations on london_proj
will work without a warning, meaning buffers can be produced from it using proper units of distance.As pointed out above, moving one degree means moving a bit more than 111 km at the equator (to be precise: 111,320 meters).This is used as the new buffer distance:
london_proj_buff = st_buffer(london_proj, 111320)
The result in Figure 6.1 (right panel) shows that buffers based on a projected CRS are not distorted:every part of the buffer’s border is equidistant to London.
Figure 6.1: Buffers around London with a geographic (left) and projected (right) CRS. The gray outline represents the UK coastline.
The importance of CRSs (primarily whether they are projected or geographic) has been demonstrated using the example of London.The subsequent sections go into more depth, exploring which CRS to use and the details of reprojecting vector and raster objects.