Fix failing cran checks

[sckott]

Please correct before 2020-04-09 to safely retain your package on CRAN.
This apparently needs proj7 to reproduce.

failing on two platforms

https://www.r-project.org/nosvn/R.check/r-patched-linux-x86_64/scrubr-00check.html
https://www.r-project.org/nosvn/R.check/r-devel-linux-x86_64-debian-gcc/scrubr-00check.html

[rsbivand]

@sckott I see the same CRAN errors with PROJ 7 on spatial revdeps, also with github master dplyr. The problem I see is (@edzer) that in st_geos_binop(), sf 0.9-0:

Browse[4]> st_crs(x) == st_crs(y)
[1] FALSE
Browse[4]> st_crs(x)
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]]
Browse[4]> st_crs(y)
Coordinate Reference System:
  User input: +proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs 
  wkt:
GEOGCRS["unknown",
    DATUM["World Geodetic System 1984",
        ELLIPSOID["WGS 84",6378137,298.257223563,
            LENGTHUNIT["metre",1]],
        ID["EPSG",6326]],
    PRIMEM["Greenwich",0,
        ANGLEUNIT["degree",0.0174532925199433],
        ID["EPSG",8901]],
    CS[ellipsoidal,2],
        AXIS["longitude",east,
            ORDER[1],
            ANGLEUNIT["degree",0.0174532925199433,
                ID["EPSG",9122]]],
        AXIS["latitude",north,
            ORDER[2],
            ANGLEUNIT["degree",0.0174532925199433,
                ID["EPSG",9122]]]]

where x is N,E, and y is E,N, so axis-swapping in the wild? It seems to happen here z <- sf::st_set_crs(z, 4326).

[edzer]

That's right, however actual axis swapping will not occur until the user sets st_axis_order(TRUE).

[rsbivand]

So the crs object expressses N/E when the data are actually still E/N? That would explain why the bbox values are unchanged, I'd wondered. How should z <- sf::st_set_crs(z, 4326) on an NA crs be written to use GIS/visualization order rather than authority order? Use sf::st_axis_order(authority_compliant = FALSE) after conditioning on sf version before st_set_crs()? The subsequent error is caused by the two crs being axis-order different.

sf::st_axis_order() is FALSE on load, so adding the line makes no difference and the errors persist. z <- sf::st_set_crs(z, "+init=epsg:4326") does work, though, but I don't know why (line 264 in scrubr/R/coord-funs.R).

[edzer]

sf uses this for == comparison of crs objects, I assumed the default would ignore axis order difference.

[rsbivand]

That is what is used in st_geom_binop(), and while GDAL only returns one test value, PROJ (now) returns three, with the third claiming to be equivalent and axis-orientation neutral (rgdal::.compare_CRS() However IIRC it still saw a difference because one has a name and a scope and the other does not - I don't know whether the name makes a difference, but the scope might, as the internal representation has more nodes. Even would know, but I hesitate to ping him. Shall I investigate further - I need to be inside running scrubr?

[edzer]

That's OK, I only assumed the reprex to be sth like

library(sf)
st_crs(4326) == st_crs("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
# [1] FALSE

[rsbivand]

I looked anyway, and had bot check properly before, the PROJ test does discriminate correctly, but cannot protect against the data having a different axis order from the CRS definition:

> st_a <- st_crs(4326)
> st_b <- st_crs("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
> st_a == st_b
[1] FALSE
> b <- sp::CRS(SRS_string=st_b$wkt)
> a <- sp::CRS(SRS_string=st_a$wkt)
> rgdal::compare_CRS(a, b)
                      strict                   equivalent 
                       FALSE                        FALSE 
equivalent_except_axis_order 
                        TRUE 

So the PROJ function used in rgdal::compare_CRS() does see that they are equivalent but not identical except axis order. I don't think that helps, though, the fix in the scrubr function is to use a Proj4 string, also as I found above "+init=epsg:4326". I feel I'm becoming deranged ...

[edzer]

I think this fixes the first one:

library(sf)
# Linking to GEOS 3.8.0, GDAL 3.0.4, PROJ 7.0.0
st_a <- st_crs(4326)
st_b <- st_crs("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
st_a == st_b
# [1] TRUE
st_axis_order(TRUE)
st_a == st_b
# [1] TRUE

I'm just suprised that the second one now is also TRUE...

[rsbivand]

For me:

library(sf)
# Linking to GEOS 3.8.1, GDAL 3.1.0dev-bc1aa43f1d, PROJ 7.0.0
st_a <- st_crs(4326)
st_b <- st_crs("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
st_a == st_b
# [1] FALSE
st_axis_order(TRUE)
st_a == st_b
# [1] TRUE

With github/master GDAL, and PROJ I think 7.0.0 released. sf is 0.9-1 github/master from last Saturday I think, not folllowing your commit above. Curious what revdeps throw up ...

[edzer]

Could you pls correct the output above (double outputs)?

[edzer]

We now see (with GDAL >= 3):

library(sf)
# Linking to GEOS 3.8.0, GDAL 3.0.4, PROJ 7.0.0
st_axis_order(FALSE)
st_a <- st_crs(4326)
st_b <- st_crs("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
st_a == st_b
# [1] TRUE
st_axis_order(TRUE)
st_a == st_b
# [1] TRUE
st_axis_order(FALSE)
st_a == st_b
# [1] TRUE

and, eh, so be it?

[rsbivand]

OK, but maybe this needs to be an sf issue, so that it doesn't get neglected?

[sckott]

thanks for this. so is the solution then to use sf::st_set_crs(z, "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs") instead of sf::st_set_crs(z, 4326)?