0

I have a Python project using GDAL 2.1.3 which creates a CRS like this:

srs = osr.SpatialReference()
srs.SetFromUserInput("EPSG:2226+5714")
print srs.GetLinearUnits()

which displays a value of 0.304800609601 which is correct as the conversion from US Survey Foot to Meter. While upgrading from Python 2.7 and GDAL 2.1.3 to Python 3.10 and GDAL 3.4.3, the same code returns 1.0. So I tried this code:

src_horz = osr.SpatialReference()
src_horz.ImportFromEPSG(2226)
print(f'{src_horz.GetLinearUnits()}, {src_horz.GetLinearUnitsName()}')

src_vert = osr.SpatialReference()
src_vert.ImportFromEPSG(5714)
print(f'{src_vert.GetLinearUnits()}, {src_vert.GetLinearUnitsName()}')

srs = osr.SpatialReference()
srs.SetCompoundCS("", src_horz, src_vert)
print(f'{srs.GetLinearUnits()}, {srs.GetLinearUnitsName()}')

Which gives me the following output:

0.30480060960121924, US survey foot
1.0, Meter
1.0, US survey foot

I believe the last line is wrong since the linear unit for the compound crs should be 0.30480060960121924.

Am I creating the CRS wrong or might I have a version mismatch or possibly a proj db issue?

cary
  • 21
  • 3
  • I think this is a bug in GDAL (or the Python bindings) that was recently fixed with [this commit](https://github.com/OSGeo/gdal/commit/fdd47d8123fa13bc4317f54700235ced788d0646) which was just included in the latest release 3.6.3. Will check that out. – cary Mar 22 '23 at 16:33

0 Answers0