1

I am supporting an package that serves {sf} objects using EPSG:4326 as default projection; I have little control over the setup my users have (some rely on defaults, others compile from source).

I have run into issue where EPSG:4326 as interpreted in context of PROJ 4.9.3 (the standard when using Windows binaries) looks differently from the same projection in context of higher versions of PROJ.

The PROJ 4.9.3 version of a crs object is

Coordinate Reference System:
  User input: EPSG:4326 
  wkt:
GEOGCS["WGS 84",
    DATUM["WGS_1984",
        SPHEROID["WGS 84",6378137,298.257223563,
            AUTHORITY["EPSG","7030"]],
        AUTHORITY["EPSG","6326"]],
    PRIMEM["Greenwich",0,
        AUTHORITY["EPSG","8901"]],
    UNIT["degree",0.0174532925199433,
        AUTHORITY["EPSG","9122"]],
    AUTHORITY["EPSG","4326"]]

while the PROJ 6.3.0 version is slightly different:

Coordinate Reference System:
  User input: EPSG:4326 
  wkt:
GEOGCRS["WGS 84",
    DATUM["World Geodetic System 1984",
        ELLIPSOID["WGS 84",6378137,298.257223563,
            LENGTHUNIT["metre",1]]],
    PRIMEM["Greenwich",0,
        ANGLEUNIT["degree",0.0174532925199433]],
    CS[ellipsoidal,2],
        AXIS["geodetic latitude (Lat)",north,
            ORDER[1],
            ANGLEUNIT["degree",0.0174532925199433]],
        AXIS["geodetic longitude (Lon)",east,
            ORDER[2],
            ANGLEUNIT["degree",0.0174532925199433]],
    USAGE[
        SCOPE["unknown"],
        AREA["World"],
        BBOX[-90,-180,90,180]],
    ID["EPSG",4326]]

This creates all sort of inconveniences as the st_crs(x) == st_crs(y) is not TRUE check is prone to fail for an object imported from a different PROJ environment.

Is there a way to serve CRS in a way that would be interpreted reliably the same? I am wary of reprojecting everything due to performance reasons.

2 Answers 2

1

Possibly if the CRS from an object has a non-NA $epsg element then the two CRSs must be identical and then you can assign one to the other if a function is complaining that the CRS of x is not the same as the CRS of y.

Otherwise if the string comparison (of the the WKT or the proj4string) returns FALSE and the CRS object comparison returns FALSE then I think you are stuck because the object comparison does a bit of normalisation with the string representation first.

In your example above, both CRS objects should have a $epsg element that is 4326, and so you can safely assign the CRS of one to the other (without reprojecting) in order for them to have identical CRS to satisfy any test in one of the intersection/overlay functions.

1
  • Hard assignment of CRS was the trick. st_crs(my downloaded object) <- 4326 this ensures that the downloaded object has the locally understood wording of WGS84. Risky with objects of unknown origin, but works for me in my use case. Oct 31, 2020 at 8:26
0

After some tweaking I found out that {sf} with PROJ 6.3.0 backend will accept objects created with PROJ 4.9.3 backend without throwing the st_crs(x) == st_crs(y) is not TRUE error (lower version PROJ will throw this error when handling objects originated on a higher version).

So the immediate solution to my problem turned out to be to make sure that all shared objects are created via the lowest commonly used PROJ version. Which is sort of a hassle, but it works.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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