[gdal-dev] OGRCreateCoordinateTransformation problem on GDAL 3.0.3 and PROJ 6.3.0

Even Rouault even.rouault at spatialys.com
Thu Jan 23 08:58:21 PST 2020


> Thanks! Yes, with SRS->importFromUserInput("EPSG:3857"); etc this works
> fine; I'm OK with the warning, I just don't see how the subsequent error
> relates to syntax that is deprecated.

I'm looking at this, but haven't yet an explanation. This is very subtle,
and seems to be linked specifically to EPSG:3857

I've modified your reproducer as the following:

int main() {

   OGRSpatialReference *aSRS = new OGRSpatialReference;

   aSRS->importFromWkt("GEOGCRS[\"WGS 84\",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]]]]");

   OGRSpatialReference *bSRS = new OGRSpatialReference;
#ifdef WORKING
   OGRSpatialReference srs;
   srs.importFromProj4("+init=epsg:3857");
   char* wktb = NULL;
   srs.exportToPrettyWkt(&wktb);
#else
   bSRS->importFromProj4("+init=epsg:3857");
   const char* wktb = "PROJCRS[\"WGS 84 / Pseudo-Mercator\",BASEGEOGCRS[\"WGS 84\",DATUM[\"World Geodetic System 1984\",ELLIPSOID[\"WGS 84\",6378137,298.257223563,LENGTHUNIT[\"metre\",1]]],PRIMEM[\"Greenwich\",0,ANGLEUNIT[\"degree\",0.0174532925199433]],ID[\"EPSG\",4326]],CONVERSION[\"unnamed\",METHOD[\"Popular Visualisation Pseudo Mercator\",ID[\"EPSG\",1024]],PARAMETER[\"Latitude of natural origin\",0,ANGLEUNIT[\"degree\",0.0174532925199433],ID[\"EPSG\",8801]],PARAMETER[\"Longitude of natural origin\",0,ANGLEUNIT[\"degree\",0.0174532925199433],ID[\"EPSG\",8802]],PARAMETER[\"False easting\",0,LENGTHUNIT[\"metre\",1],ID[\"EPSG\",8806]],PARAMETER[\"False northing\",0,LENGTHUNIT[\"metre\",1],ID[\"EPSG\",8807]]],CS[Cartesian,2],AXIS[\"(E)\",east,ORDER[1],LENGTHUNIT[\"metre\",1,ID[\"EPSG\",9001]]],AXIS[\"(N)\",north,ORDER[2],LENGTHUNIT[\"metre\",1,ID[\"EPSG\",9001]]]]";
#endif
   bSRS->importFromWkt((const char *) wktb);

   OGRCoordinateTransformation *ct =
OGRCreateCoordinateTransformation(aSRS, bSRS);
   if (ct == NULL) {
        printf("ct NULL\n");
        exit(1);
   }
   exit(0);
}

So the error is linked to having importFromProj4() and then importFromWkt()
This doesn't make sense as the later should cancel the former, so there's some
side effect of importFromProj4() that has later consequences.
If you define WORKING, which will do the importFromProj4() + exportToPrettyWkt() 
in a temporary object, and import the resulting WKT in the final bSRS, then
it works. Doesn't make sense either but could be a workaround

But Sean's advice is definitely the way forward . Use importFromEPSG(XXXX) or
SetFromUserInput("EPSG:XXXX")
And SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER) if you don't
want EPSG compliant axis definitions.

Even

-- 
Spatialys - Geospatial professional services
http://www.spatialys.com


More information about the gdal-dev mailing list