Rasterio: WKT representation of CRS is broken

[drnextgis]

In a previous versions of rasterio WKT representation of EPSG:3857 looks good:

>>> import rasterio                                                                                                                                                                                             
>>> rasterio.__version__                                                                                                                                                                                        
'1.0.13'
>>> from rasterio.crs import CRS                                                                                                                                                                                
>>> crs = CRS(init="epsg:3857")                                                                                                                                                                                 
>>> crs.wkt 
'PROJCS["WGS 84 / Pseudo-Mercator",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"]],PROJECTION["Mercator_1SP"],PARAMETER["central_meridian",0],PARAMETER["scale_factor",1],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["X",EAST],AXIS["Y",NORTH],EXTENSION["PROJ4","+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +wktext +no_defs"],AUTHORITY["EPSG","3857"]]'

In 1.0.17 it seems like something was broken:

>>> import rasterio
>>> rasterio.__version__ 
'1.0.17'
>>> from rasterio.crs import CRS
>>> crs = CRS(init="epsg:3857")
>>> crs.wkt
'PROJCS["unnamed",GEOGCS["unnamed ellipse",DATUM["unknown",SPHEROID["unnamed",6378137,0],EXTENSION["PROJ4_GRIDS","@null"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433]],PROJECTION["Mercator_2SP"],PARAMETER["standard_parallel_1",0],PARAMETER["central_meridian",0],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["Meter",1],AUTHORITY["EPSG","3857"]]'

[sgillies]

Yes. 1.0.13 with GDAL 2.4.0 gives

>>> CRS(init='epsg:3857').wkt
'PROJCS["WGS 84 / Pseudo-Mercator",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"]],PROJECTION["Mercator_1SP"],PARAMETER["central_meridian",0],PARAMETER["scale_factor",1],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["X",EAST],AXIS["Y",NORTH],EXTENSION["PROJ4","+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +wktext +no_defs"],AUTHORITY["EPSG","3857"]]'

so this is a change in Rasterio.

I've always been baffled by the OGR SRS API. I think I had better ask Even for some help.

[sgillies]

Here's how the regression came about: in 1.0.13, when the CRS constructor gets init='epsg:3857' we picked off the 3857 and then called OSRImportFromEPSG(3857), which does the right thing. I'm not sure how we determined that was the correct pattern, it might have been accidental. In 1.0.16, we do not pick the epsg code out of the keyword argument value, instead we concatenate all the keyword args and call OSRImportFromPROJ4("+init=epsg:3857"), which does the wrong (or unexpected) thing.

I asked on gdal-dev why there is a difference and hope to get some insight tomorrow.

[sgillies]

I got the clarification I was looking for: use of the older PROJ4 +init=epsg:3857 syntax within GDAL is deprecated. We'll go back to always using OSRImportFromEPSG when possible and this regression will be fixed.

[sgillies]

Done.

[drnextgis]

@sgillies thank you!