6

I have some points in EPSG:2284 and I'm trying to convert them to EPSG:4326 using pyproj, however I get weird results. When I go to the spatial reference website I see that there is an endpoint that does the conversion. When I use that instead I get the correct result.

Here is some code illustrating the different results:

from pyproj import Proj, transform
import requests

if __name__ == '__main__':

    # points I want to convert from EPSG:2284 -> EPSG:4326
    targets = [
        (3669486.63386,11293013.10592,709.08004),
        (3669559.13811,11292972.72831,711.60055),
        (3669639.51308,11292851.48264,712.22258),
        (3669800.62304,11292949.38118,714.75766),
    ]

    for northing, easting, up in targets:
        print "Converting", easting, northing

        # Transform using pyproj (gives wrong answer)
        WGS84 = Proj(init='EPSG:4326')
        inp = Proj(init='EPSG:2284')
        x, y = transform(inp, WGS84, easting, northing)
        print x, y, "(with proj4)"

        # Transform using spatialreference.org URL (gives correct answer)
        url = "http://spatialreference.org/projection/?json=%7B%22type%22%3A%22Feature%22%2C%20%22geometry%22%3A%7B%22type%22%3A%22Point%22%2C%20%22coordinates%22%3A%5B{lon}%2C%20{lat}%5D%7D%2C%22properties%22%3A%7B%7D%7D&inref=EPSG:{inp}&outref=epsg:{outp}"
        url = url.format(lon=easting, lat=northing, inp='2284', outp='4326')
        x, y = requests.get(url).json()['coordinates']
        print x, y, '(using spatialreference.org URL)'

Does anyone know why they are giving different result and what is wrong with the Pyproj code?

8

Simply use (Converting elevations into correct units with pyproj?, Proj4 String for NAD83(2011) / Louisiana South (ftUS), ...)

preserve_units=True (as you say, pyproj assumes that your coordinates are in meters, therefore)

for northing, easting, up in targets:
    print "Converting", easting, northing
    # Transform using pyproj (gives wrong answer)
    WGS84 = Proj(init='EPSG:4326')
    inp = Proj(init='EPSG:2284',preserve_units=True)
    ...

Results

Converting 11293013.1059 3669486.63386
-79.1537547735 37.3989918968 (with proj4)
-79.153755 37.398992 (using spatialreference.org URL)
Converting 11292972.7283 3669559.13811
-79.1538955003 37.3991902547 (with proj4)
-79.153896 37.39919 (using spatialreference.org URL)
Converting 11292851.4826 3669639.51308
-79.1543148015 37.3994086895 (with proj4)
-79.154315 37.399409 (using spatialreference.org URL)
Converting 11292949.3812 3669800.62304
-79.1539816309 37.3998530254 (with proj4)
-79.153982 37.399853 (using spatialreference.org URL)

Supplementary control with GDAL/OSR

from osgeo import osr
wgs84 = osr.SpatialReference()
wgs84.ImportFromEPSG(4326)
inp= osr.SpatialReference()
inp.ImportFromEPSG(2284)
transformation = osr.CoordinateTransformation(inp,wgs84)
for northing, easting, up in targets:
    print transformation.TransformPoint(easting,northing)

(-79.15375477347949, 37.3989918968496, 0.0)
(-79.15389550031416, 37.39919025468435, 0.0)
(-79.15431480154861, 37.39940868945351, 0.0)
(-79.15398163085317, 37.39985302535012, 0.0)
1

Pyproj expects meters not feet.

We're looking for long answers that provide some explanation and context. Don't just give a one-line answer; explain why your answer is right, ideally with citations. Answers that don't include explanations may be removed.

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.