Carlos
Carlos

Reputation: 4637

Different calculations between ST_Transform and CoordTransform that involves presition

I am working with GeoDjango and Postgis, and doing transformation on a point from 4326 to 24877, when I have checked the result, it is not the same

My values are:

Latitude -16.42238172128686
Longitude -71.47541752550751

Initial position

Using long lat as x, y coord

In [63]: pnt = Point(-71.47541752550751, -16.42238172128686, srid=4326)
In [64]: pnt.transform(CoordTransform(SpatialReference(4326), SpatialReference(24877)))
In [65]: pnt.x, pnt.y, pnt.srid
Out[65]: (1521142.247877425, 8160547.8770980425, 24877)

Inverse

In [66]: pnt = Point( 1521142.247877425,8160547.8770980425, srid=24877 )
In [67]: pnt.transform(CoordTransform(SpatialReference(24877), SpatialReference(4326)))
In [68]: pnt.x, pnt.y, pnt.srid
Out[68]: (-71.47541753138003, -16.42238165040434, 4326)

Position after calculus

SQL on Postgres

SELECT ST_AsText(ST_Transform(ST_SetSRID(ST_MakePoint(-71.47541752550751, -16.42238172128686), 4326),24877))

POINT(1520903.86571082 8160169.90886929)

Taking data from postgis results

In [69]: pgs = Point(1520903.86571082, 8160169.90886929, srid=24877)
In [70]: pgs.transform(CoordTransform(SpatialReference(24877), SpatialReference(4326)))
In [71]: pgs.x, pnt.y, pgs.srid
Out[71]: (-71.47745375612124, -16.42238165040434, 4326)

Image on sql results

As you could see there is a difference between geodjango and postgres, on at least 300 meters, and that means different places

After reading and asking on some forums and chats, I have found that this is related to presition, I was asking more about it on internet, as a newbie on this topic, this is getting some heard to fin.

I am looking for information that give me enough knowledge for continue working on it, for example:

  1. The right way on do that
  2. The reason on this difference
  3. how many meters are the difference between systems or the right way on calculate that
  4. how to validate this positions
  5. what is the best way on work this information

Upvotes: 2

Views: 807

Answers (2)

Paul Norman
Paul Norman

Reputation: 201

If you look at http://spatialreference.org/ref/epsg/24877/ you can see that EPSG 24877 is only valid between -84.0000, -10.4000, -78.0000, 0.0000, and the point you're trying to transform is outside that (-71.475 isn't within -84,-78 and -16.422 isn't within -10.4,0)

A secondary issue might be if a transform from WGS84 to PSAD56. There's a grid shift involved there. Unless you have some data which is already in PSAD56 I'd suggest WGS 84 / UTM zone 19S (EPSG:32719) which has bounds of -72.0000, -80.0000, -66.0000, 0.0000.

If you absolutely need PSAD56 you might have more luck with PSAD56 / UTM zone 19S (EPSG:24879) but there still might be grid shift issues.

Upvotes: 3

Urth
Urth

Reputation: 21

csotelo: you have lat/long switched around, that point doesn't seem to exist in 24877

>>> p = Point(-16.42238172128686,-71.47541752550751, srid=4326)
>>> p.transform(24877)
>>> print p
POINT (2388367.8123248801566660 916017.4938517492264509)

# SELECT ST_AsText(ST_Transform(ST_SetSRID(ST_MakePoint(-16.42238172128686,-71.47541752550751), 4326),24877));
                st_astext                 
------------------------------------------
 POINT(2388367.81232488 916017.493851751)

Point(-71.47541752550751, -16.42238172128686, srid=4326) is outside the region defined by 24877 and therefore can only be approximated. You'll find that that point is actually near antarctica. https://maps.google.nl/maps?q=-71.47541752550751,+-16.42238172128686&hl=en&ll=-71.475417,-16.422382&spn=60.209188,270.527344&sll=53.157645,5.636521&sspn=1.462449,4.22699&t=m&z=3

Upvotes: 2

Related Questions