deprecate exportToProj4

[edzer]

See here.

I guess this means that when specifying a CRS with an EPSG code, we then pick up the WKT and show (parts of) that on print, and ignore the proj4string that exportToProj4 would give altogether, unless this is explicitly asked for e.g. in the st_as_text.crs method.
#1146
r-spatial/discuss#28

Pretty drastic; discussion welcome.

[mdsumner]

By deprecate do you mean to ignore the proj4string by not creating it at all? Or, maintain the creation but start concealing it for a future hard-removal?

[rsbivand]

For sp in rgdal, I'm looking at keeping the existing CRS object, but piggy-backing a comment with the WKT2_2018 string. The PROJ string in the CRS object may be displayed, but the WKT2 string will be used, not the PROJ string. I'm also looking at instantiating CRS from a PROJ string as before, but checking will generate a WKT2 comment, and as importFromProj4() doesn't drop the tags that we need, the WKT2 should represent the input PROJ string OK. More to follow - we can test things before seeing what to do with sf.

[mdsumner]

Cool thanks

[edzer]

@rsbivand did you btw try to add a slot, rather than a comment, to CRS objects in sp? Does this break any dependencies?

[edzer]

@mdsumner not sure yet - crs objects will need to become very clear about what is being pushed down to PROJ.

[rsbivand]

Wrt. S4 slots, my understanding is that the main specific danger is that serialised objects with class definitions that differ from the available package version will be invalid. Since we do not know that users' workflows do not involve saving and loading serialised objects, we should avoid breaking changes involving S4 objects. Since S3 objects need to test for NULL components and attributes anyway, S3 class changes are only breaking if badly written (that is fail to check for NULL returns on access). The comment() trick also depends on the same need to check for NULL on return, but avoids breaking S4 changes. I'll see if I have time to show this, I'm looking at the Green book and the 2000 Venables & Ripley S Programming book.

The serialisation issue is real:

library(methods)
setClass("CRS", slots = c(projargs = "character"))
setValidity("CRS", function(object) {
		if (length(object@projargs) != 1)
			return("projargs must be of length 1")
	})
(res <- new("CRS", projargs="+proj=longlat +ellps=WGS84"))
validObject(res)
tf <- tempfile(fileext="rds")
saveRDS(res, file=tf)
rm(res)
removeClass("CRS")
setClass("CRS", slots = c(projargs = "character", WKT2 = "character"))
setValidity("CRS", function(object) {
		if (length(object@projargs) != 1 || !nzchar(object@projargs))
			return("projargs must be of length 1")
		if (length(object@WKT2) != 1 || !nzchar(object@WKT2))
			return("WKT2 must be of length 1")
	})
res_ser <- readRDS(tf)
try(validObject(res_ser))

with:

> res_ser <- readRDS(tf)
Warning: Not a validObject(): no slot of name "WKT2" for this object of class "CRS"
> try(validObject(res_ser))
Error in validObject(res_ser) : 
  invalid class “CRS” object: slots in class definition but not in object: "WKT2"

[rsbivand]

Another point: rgdal now can work with old CRS (for GDAL < 3 & PROJ < 6), and with piggy-backed WKT2 CRS (sp version unchanged, GDAL >=3 & PROJ >= 6). If we change the CRS definition in sp, we'll need to set NA on a new WKT2 slot for GDAL < 3 & PROJ < 6, and make rgdal depend on the sp with the new CRS definition.

[rsbivand]

This is the WIP-vignette from R-Forge/rgdal 1.5-1 rev 888. Comments welcome!
PROJ6_GDAL3.zip

RFC now on http://rgdal.r-forge.r-project.org/articles/PROJ6_GDAL3.html, no need to unzip.

[rsbivand]

Note that *.rda/*.rds objects can contain "CRS" etc. objects with PROJ strings with leading spaces - such spaces have been stripped off for many years, but some legacy data() objects still exist. They are caught in rgdal::showSRID() with stop("string starts with space"). r-spatial/discuss#28

[rsbivand]

Video of talk on this problem on: https://www.youtube.com/playlist?list=PLXUoTpMa_9s10NVk4dBQljNOaOXAOhcE0, third in playlist; presentation at https://rsbivand.github.io/ECS530_h19/ECS530_III.html.

[edzer]

See also here: https://lists.osgeo.org/pipermail/gdal-dev/2019-December/051344.html

From reading that, my feeling is that we're going to have to not even create proj4strings when they were not used to initialize the CRS (e.g. read data from file, from EPSG, or from WKT).

[rsbivand]

I'll shift from GDAL 3.0.2 to current master to see. Looks as though Frank Warmerdam's fix via the raster OVERRIDE_PROJ_DATUM_WITH_TOWGS84 setting on read has reached EOL. Means that raster will be as severely affected as vector - I'm unsure whether we should use the proposed tricks to keep the +towgs84= tag, or just jump to WKT2 and choice of pipeline?

[rsbivand]

The master make failed for me, I'll continue to try to get it built.

Update, I forgot to make distclean, built now.

[rsbivand]

I guess for the +towgs84 changes we need to condition on GDAL > 3.0.2 ? Master is 3.1.0.

[rsbivand]

@edzer Practical question: when I declare an OGRSpatialReference directly:

OGRSpatialReference hSRS = (OGRSpatialReference) NULL;
hSRS.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);

compiles without trouble.

But when declared as a pointer, what I thought was equivalent hSRS-> leads to grief. Mostly in ogrP4S in rgdal/src/ogr_proj.cpp and elsewhere in rgdal/src/OGR_write.cpp (local, not committed to SVN). I guess I need a break and a coffee.

[edzer]

You can't -> on a NULL pointer: you need to allocate an object first. See e.g. here.

[rsbivand]

The problem was on read; I'd forgotten to check whether import... returned NULL, so the object pointer became NULL when a shapefile had no *.prj file. On testing more carefully, now OK. Committed as rev. 902 with lots of exposed Rprintf() to show GetAxisMappingStrategy() before and after setting to OAMS_TRADITIONAL_GIS_ORDER. I'll comment out this noise when I can see the outcomes better.