4

I am converting a point from wgs84 into a state plane coord system (epsg: 32051, which is defined in US survey feet) using pyproj, but the transform isn't changing the z values into feet:

point_4326 = (-81.23127933,37.70957554,668.105)
proj1 = pyproj.Proj(init='epsg:4326')
proj2 = pyproj.Proj(init='epsg:32051')
print pyproj.transform(proj1, proj2, *point_44326)

$ (589192.2346128735, 78765.53566868477, 668.105)

I have also tried specifying the exact proj4 string for the new coord system (which includes the USft scaling factor), but get the exact same result:

proj2 = pyproj.Proj('+proj=lcc +lat_1=37.48333333333333 +lat_2=38.88333333333333 +lat_0=37 +lon_0=-81 +x_0=609601.2192024384 +y_0=0 +ellps=clrk66 +datum=NAD27 +to_meter=0.3048006096012192 +no_defs')
print pyproj.transform(proj1, proj2, *point_4326) 

$ (589192.2346128735, 78765.53566868477, 668.105)

Any ideas how I can get pyproj to convert the elevation values to the correct units, or if there are similar tools / libraries which can do this?

4

Look at Proj4 String for NAD83(2011) / Louisiana South (ftUS)

preserve_units=True fix the problem (pyproj silently changes '+units=' parameter)

import pyproj
proj2 = pyproj.Proj(init='epsg:32051',preserve_units=True)
proj1 = pyproj.Proj(init='epsg:4326')
print pyproj.transform(proj1, proj2, *point_4326)
(1933041.523059067, 258416.59493967678, 2191.941154166668)

Control with GDAL/OSR

from osgeo import osr
wgs84 = osr.SpatialReference()
wgs84.ImportFromEPSG(4326)
proj2= osr.SpatialReference()
proj2.ImportFromEPSG(32051)
transformation = osr.CoordinateTransformation(wgs84,proj2)
result = transformation.TransformPoint(*point_4326)
print result
(1933041.5230590631, 258416.59493967678, 2191.941154166668)

And 2192 feets = 668.122 meters

| improve this answer | |

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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