[postgis-devel] ST_Project inaccuracy? Values out by many miles...

Brent Wood pcreso at pcreso.com
Mon Dec 16 14:00:11 PST 2013






Hi,

I have been looking at using ST_Project() to solve a problem I have (& in theory it is perfect!).

It seems to introduce some inaccuarcies, which can be quite significant in some use cases. The SQL below projects a new point 2km E, then back 2km W. The result should be in the same place we started, but is not, by a larger value than I'm comfortable with. The source point location comes from a real dataset.


select Postgis_full_version();

POSTGIS="2.0.1 r9979" GEOS="3.3.8-CAPI-1.7.8" PROJ="Rel. 4.8.0, 6 March 2012" GDAL="GDAL 1.9.2, released 2012/10/08" LIBXML="2.7.8" TOPOLOGY RASTER



select ST_AsText(point) as point,
       ST_AsText(ST_Project(point, 2000,radians(90))) as proj2k,
       ST_AsText(ST_Project(ST_Project(point, 2000,radians(90)), 2000,radians(270))) as proj2k0
from (select ST_AsText(ST_SetSRID(ST_MakePoint(171.0259489718, -45.2961933568),4326))::geography as point) as foo;

point   | POINT(171.0259489718   -45.2961933568)
proj2k  | POINT(171.051446314948 -45.2961905108324)
proj2k0 | POINT(171.025948973075 -45.296187664865)


The final lat & lon should be the same as the original, but both have changed, the latitude substantially more than the longitude - even though the point latitude should be unchanged (it only moved E & W, not N or S) 

It gets worse with larger distances, eg, at 2000km:

select ST_AsText(point) as point,
       ST_AsText(ST_Project(point, 2000000,radians(90))) as proj2000k,
       ST_AsText(ST_Project(ST_Project(point, 2000000,radians(90)), 2000000,radians(270))) as proj2000k0
from (select ST_AsText(ST_SetSRID(ST_MakePoint(171.0259489718,
 -45.2961933568),4326))::geography as point) as foo;

point      | POINT(171.0259489718 -45.2961933568)
proj2000k  | POINT(-164.264848338683 -42.5385273590618)
proj2000k0 | POINT(172.015779821441 -40.0225960863261)

it causes over 5 degree change in the latitude.

Is this a bug I should be submitting or have I missed something?


Thanks

  Brent Wood
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/postgis-devel/attachments/20131216/f4ae9bff/attachment.html>


More information about the postgis-devel mailing list