[gdal-dev] ogr2ogr PG->SHP lacks EPSG code

[Even Rouault]

Le samedi 12 novembre 2016 23:27:30, Martin Landa a écrit :
> Hi,
> 
> 2016-11-11 18:56 GMT+01:00 Even Rouault <even.rouault at spatialys.com>:
> > The WKT variant of ESRI Shapefile .prj doesn't have EPSG codes (nor
> > TOWGS84 node by the way). At least that was how it was historically and
> > GDAL "morphs" OGC WKT to ESRI WKT this way (there are other variations
> > in names of projection methods and parameters). Regarding the EPSG code,
> > for some reason, I believe this might perhaps have changed in recent
> > ArcGIS versions but I'm not sure. Perhaps someone could confirm
> 
> you are right, EPSG codes are missing when producing SHP in ArcGIS
> (tested 10.4 version) [1]. The problem is that such Shapefiles are not
> recognized by QGIS as EPSG 5514 (QGIS shows user-defined SRS 1000).
> Any idea? Thanks, Martin

Probably because of the missing TOWGS84 node in the ESRI WKT:

$ gdalsrsinfo ESRI::Export_Output.prj.txt

PROJ.4 : '+proj=krovak +lat_0=49.5 +lon_0=24.83333333333333 
+alpha=30.28813975277778 +k=0.9999 +x_0=0 +y_0=0 +ellps=bessel +units=m 
+no_defs '

versus the srs.db in QGIS resources has

$ ogrinfo resources/srs.db -sql "select parameters from tbl_srs where 
srid=5514" -q

Layer name: SELECT
OGRFeature(SELECT):0
  parameters (String) = +proj=krovak +lat_0=49.5 +lon_0=24.83333333333333 
+alpha=30.28813972222222 +k=0.9999 +x_0=0 +y_0=0 +ellps=bessel 
+towgs84=542.5,89.2,456.9,5.517,2.275,5.516,6.96 +pm=greenwich +units=m 
+no_defs

One of the issue is that qgis srs.db has a towgs84 value different from the 
GDAL one. As QGIS tries to add the missing towgs84 from the EPSG database 
sync'ed frm GDAL from the name of the datum, it tries with 
"+towgs84=589,76,480,0,0,0,0" which is the shift for 
System_Jednotne_Trigonometricke_Site_Katastralni datum, hence the mismatch

(Not sure if the additional +pm=greenwich in qgis srs.db might also be a 
problem in the matching logic)

Looking at git history, I see:

commit 2ffbbded16fa8f4425ccb1be4a04186a43b93754
Author: Radim Blazek <radim.blazek at gmail.com>
Date:   Wed Jun 12 11:12:56 2013 +0200

    Updated towgs84 parameters of JTSK CRS, the towgs84 values were calculated 
for area of Czech and Slovak Republic as it is described here: 
http://grass.fsv.cvut
    
    UPDATE tbl_srs set parameters = '+proj=krovak +lat_0=49.5 
+lon_0=24.83333333333333 +alpha=30.28813972222222 +k=0.9999 +x_0=0 +y_0=0 
+ellps=bessel +towgs84=542.
    
    UPDATE tbl_srs set parameters = '+proj=krovak +lat_0=49.5 
+lon_0=24.83333333333333 +alpha=30.28813972222222 +k=0.9999 +x_0=0 +y_0=0 
+ellps=bessel +towgs84=542.
    
    UPDATE tbl_srs set parameters = '+proj=krovak +lat_0=49.5 +lon_0=42.5 
+alpha=30.28813972222222 +k=0.9999 +x_0=0 +y_0=0 +ellps=bessel 
+towgs84=542.5,89.2,456.9,
    
    UPDATE tbl_srs set parameters = '+proj=krovak +lat_0=49.5 +lon_0=42.5 
+alpha=30.28813972222222 +k=0.9999 +x_0=0 +y_0=0 +ellps=bessel 
+towgs84=542.5,89.2,456.9,
    
    UPDATE tbl_srs set parameters = '+proj=krovak +lat_0=49.5 
+lon_0=24.83333333333333 +alpha=30.28813972222222 +k=0.9999 +x_0=0 +y_0=0 
+ellps=bessel +towgs84=542.

I guess QGIS should also try to match the 'description' column of srs.db with 
the name of the PROJCS gotten from GDAL (with some laundering logic for 
spaces, / to underscore, as done when exporting to ESRI), instead/in addition 
to proj.4 strings

Welcome to SRS hell !


> 
> [1] http://147.32.131.147/landa/Export_Output.prj.txt

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