[gdal-dev] Differences between OSRImportFromEPSG and OSRImportFromProj4

classic Classic list List threaded Threaded
3 messages Options
Reply | Threaded
Open this post in threaded view
|

[gdal-dev] Differences between OSRImportFromEPSG and OSRImportFromProj4

Sean Gillies-3
Hi all,

I observe a difference between the WKT exports of OSRImportFromEPSG(osr, 3857) and OSRImportFromProj4(osr, "+init=epsg:3857"). The former gives a named coordinate system, the latter gives an unnamed coordinate system.

Case 1:

OSRImportFromEPSG(osr, 3857)
OSRExportToWkt(osr)

'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"]]'

Case 2:

OSRImportFromProj4(osr, "+init=epsg:3857")
OSRExportToWkt(osr)

'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"]]'

Am I wrong to understand the semantics of "+init=epsg:3857" as "open the epsg database, get entry 3857, and use that"? Is this different from specifying the entry directly?

Thanks in advance,

--
Sean Gillies

_______________________________________________
gdal-dev mailing list
[hidden email]
https://lists.osgeo.org/mailman/listinfo/gdal-dev
Reply | Threaded
Open this post in threaded view
|

Re: Differences between OSRImportFromEPSG and OSRImportFromProj4

Even Rouault-2
Hi Sean,

>
> I observe a difference between the WKT exports of OSRImportFromEPSG(osr,
> 3857) and OSRImportFromProj4(osr, "+init=epsg:3857"). The former gives a
> named coordinate system, the latter gives an unnamed coordinate system.
>
>
> Am I wrong to understand the semantics of "+init=epsg:3857" as "open the
> epsg database, get entry 3857, and use that"? Is this different from
> specifying the entry directly?

That's an case where my recent work will help reducing the confusion
(hopefully)
In current released versions, OSRImportFromEPSG() read in the GDAL data/*.csv
files a fully-fledged description of the geodetic object.
OSRImportromProj4("+init=epsg:XXXX") goes to PROJ which opens the $prefix/
share/proj/epsg file, identifies the PROJ string in it that corresponds to the
code, returns it to GDAL, and GDAL reconstructs a OGRSpatialReference by
analysing the PROJ4 string, which has no names.

With the RFC73 related work, now there's a single database. The
+init=epsg:XXXX syntax is still supported as a legacy mechanism by PROJ, but
it no longer opens a 'epsg' file but rather goes to the proj.db database. Note
that for backward compatibility purposes, if you use the +init=epsg:XXXX
syntax, PROJ will return CRS objects with possibly a non-EPSG compliant axis
order (for example it will return longitude-latitude for geographic CRS,
whereas EPSG mandates latitude-longitude). GDAL (with the commit I just pushed
since the +init=epsg:XXXX case was actually broken in master) will now emit a
warning if you use that legacy syntax.

$ gdalsrsinfo EPSG:4326 -o WKT2

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]]

vs

$ gdalsrsinfo +init=epsg:4326 -o WKT2
Warning 1: +init=epsg:XXXX syntax is deprecated. It might return "
a CRS with a non-EPSG compliant axis order.

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]]],
    USAGE[
        SCOPE["unknown"],
        AREA["World"],
        BBOX[-90,-180,90,180]]]

In the later, the axis roder is long/lat and the ID["EPSG",4326] node has been
stripped off since the object is not conformant with EPSG.

Even

--
Spatialys - Geospatial professional services
http://www.spatialys.com
_______________________________________________
gdal-dev mailing list
[hidden email]
https://lists.osgeo.org/mailman/listinfo/gdal-dev
Reply | Threaded
Open this post in threaded view
|

Re: Differences between OSRImportFromEPSG and OSRImportFromProj4

Sean Gillies-3
Hi Even,

On Thu, Feb 7, 2019 at 2:13 AM Even Rouault <[hidden email]> wrote:
Hi Sean,

>
> I observe a difference between the WKT exports of OSRImportFromEPSG(osr,
> 3857) and OSRImportFromProj4(osr, "+init=epsg:3857"). The former gives a
> named coordinate system, the latter gives an unnamed coordinate system.
>
>
> Am I wrong to understand the semantics of "+init=epsg:3857" as "open the
> epsg database, get entry 3857, and use that"? Is this different from
> specifying the entry directly?

That's an case where my recent work will help reducing the confusion
(hopefully)
In current released versions, OSRImportFromEPSG() read in the GDAL data/*.csv
files a fully-fledged description of the geodetic object.
OSRImportromProj4("+init=epsg:XXXX") goes to PROJ which opens the $prefix/
share/proj/epsg file, identifies the PROJ string in it that corresponds to the
code, returns it to GDAL, and GDAL reconstructs a OGRSpatialReference by
analysing the PROJ4 string, which has no names.

With the RFC73 related work, now there's a single database. The
+init=epsg:XXXX syntax is still supported as a legacy mechanism by PROJ, but
it no longer opens a 'epsg' file but rather goes to the proj.db database. Note
that for backward compatibility purposes, if you use the +init=epsg:XXXX
syntax, PROJ will return CRS objects with possibly a non-EPSG compliant axis
order (for example it will return longitude-latitude for geographic CRS,
whereas EPSG mandates latitude-longitude). GDAL (with the commit I just pushed
since the +init=epsg:XXXX case was actually broken in master) will now emit a
warning if you use that legacy syntax.

$ gdalsrsinfo EPSG:4326 -o WKT2

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]]

vs

$ gdalsrsinfo +init=epsg:4326 -o WKT2
Warning 1: +init=epsg:XXXX syntax is deprecated. It might return "
a CRS with a non-EPSG compliant axis order.

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]]],
    USAGE[
        SCOPE["unknown"],
        AREA["World"],
        BBOX[-90,-180,90,180]]]

In the later, the axis roder is long/lat and the ID["EPSG",4326] node has been
stripped off since the object is not conformant with EPSG.

Thank you for the explanation and the fix. Thanks, too, for the details in https://trac.osgeo.org/gdal/wiki/rfc73_proj6_wkt2_srsbarn#OGRCoordinateTransformationchanges and https://trac.osgeo.org/gdal/wiki/rfc73_proj6_wkt2_srsbarn#Axisorderissues which make a lot of sense to me and answer all of the follow up questions that came to my mind.

--
Sean Gillies

_______________________________________________
gdal-dev mailing list
[hidden email]
https://lists.osgeo.org/mailman/listinfo/gdal-dev