Can somebody explain why the result of osr.SpatialReference() differ between Python 2 and 3?

sr = osr.SpatialReference()
sr.SetFromUserInput('EPSG:28992')
sr.ExportToProj4()

In Python 2.7 this gives me

Out[15]: '+proj=sterea +lat_0=52.15616055555555 +lon_0=5.38763888888889 +k=0.9999079 +x_0=155000 +y_0=463000 +ellps=bessel +towgs84=565.417,50.3319,465.552,-0.398957,0.343988,-1.8774,4.0725 +units=m +no_defs '

In Python 3 though, I get a different answer:

Out[14]: '+proj=sterea +lat_0=52.15616055555555 +lon_0=5.38763888888889 +k=0.9999079 +x_0=155000 +y_0=463000 +ellps=bessel +towgs84=565.2369,50.0087,465.658,-0.406857,0.350733,-1.87035,4.0812 +units=m +no_defs'
up vote 7 down vote accepted

The projections differ in the datum shift values for the conversion from Amersfoort datum to WGS84.

The first one is tfm code 4833, and the second is tfm code 15934.

GDAL 2.1.0 uses tfm code 4833, while GDAL 2.2.0 and later uses tfm code 15934. The change was done in https://trac.osgeo.org/gdal/changeset/37081

According to the remarks, 4833 is the latest definition. See also https://trac.osgeo.org/gdal/ticket/2487 for the difference.

If you want higher precision, you could use a grid shift file for the transformation, see http://osgeo-org.1560.x6.nabble.com/Amersfoort-RD-New-td5293753.html

  • so his system uses different gdal versions? – nickves Nov 16 at 11:12
  • Unless they compiled the Python bindings themselves, that seems likely. – bugmenot123 Nov 16 at 12:29
  • Excellent answer. – LarsVegas Nov 16 at 18:32
  • 1
    @bugmenot123 it is not necessary to compile the Python bindings yourself, you can just exchange the csv files to get the transformation synchronous in both pythons. – AndreJ Nov 19 at 10:36

Your Answer

By clicking "Post Your Answer", you acknowledge that you have read our updated terms of service, privacy policy and cookie policy, and that your continued use of the website is subject to these policies.

Not the answer you're looking for? Browse other questions tagged or ask your own question.