<span style="font-weight: bold;"><span style="font-weight: bold;"><span style="font-weight: bold;">we make datum shift with postgis with custom srids, like this... search in the lis for this email and this subject: </span>

</span></span><b>Datum Shift Problems with PostGIS<br><br><br>Regards and contact me for any questions... the following problem was solved and currently we make 3 an 7 parameters datum transformations<br><br><br>REGARDS</b>

<span style="font-weight: bold;"><span style="font-weight: bold;"><span style="font-weight: bold;"> Datum shift Problems with PostGIS<span class="gmail_quote"> Dear List Users, In order to make datum change with proj from command line, i have the following example that works fine.

<br><br>Datum Change from Sad69/utm 19s to wgs84/utm 19s<br> <br># cs2cs +proj=utm +zone=19 +south +ellps=GRS67 <span style="font-weight: bold;">+towgs84=-75,-1,-44</span>  +unit=m +to +proj=utm +zone=19 +south +datum=WGS84 

<br><br>In Order to make the same procces, but using Postgis we have

added a record en the spatial_ref_sys table, with a custom epsg code to

get a personalized transformation, and in en the proj4text field the

following sentence (like the parameters of cs2cs) has been added.

<br><br>Record added to spatial_ref_sys Table

<br><br>srid : <span style="font-weight: bold;">50000</span><br>auth_name : EPSG<br>auth_srid: 50000<br>srtext : Sad69 a wgs84<br>proj4text : +proj=utm +zone=19 +south +ellps=GRS67 <span style="font-weight: bold;">+towgs84=75,1,44

</span>  +unit=m +to +proj=utm +zone=19 +south +datum=WGS84 +no_defs<br>

<br>To run this custom datum transformation we use the transform()

function build in a query inserted as parameter in the pgsql2shp

command.<br><br>pgsql2shp -h localhost -p 5432 -u pgsql -P pgsql -f shapefile_name.shp db_name 

<br>"select transform(setsrid(the_geom<div>,29179),<span style="font-weight: bold;">50000</span>) as the_geom from nombre_tabla"<br><br>(we have our cartography with SRID=-1, reason why we use the setsrid function)

<br>

<br>we get a shapefile correctly transformed to wgs84/utm 19s Datum<br> <br><span style="font-weight: bold;">Here

is the first Question: Why we need to change the sign in the

transformation parameters, +towgs84=-75,-1,-44 in command line and

+towgs84=75,1,44 in the proj4text fiel of the spatial_ref_sys table

?????

<br><br></span>we had many problems with this, we solve it, but we don't know why we need to change it if postgis works with proj.<br><br>Now,

we try to make the inverse process, in order to make a datum shift from

WGS84 utm 19s to sad69 utm 19s, and the proj documentation tells that

the -I parameter is the only that we need to make the inverse datum

shift. With cs2cs command the -I parameter works fine, but the

inclusion of this parameter in the proj4text field in the

spatial_ref_sys table does not make a correct transformation, giving us

around 1000 meters of displacement.

<br><br><span style="font-weight: bold;">Second Question: How can we make this inverse transformation with postgis (wgs84 to psad56) ?</span><br><span><br clear="all"></span></div><br><br><div><span class="gmail_quote">2006/5/15, James G Wilkinson <

<a href="mailto:jgw@alpinegeophysics.com">jgw@alpinegeophysics.com</a>>:</span><blockquote class="gmail_quote" style="border-left: 1px solid rgb(204, 204, 204); margin: 0pt 0pt 0pt 0.8ex; padding-left: 1ex;">Mr. Warmerdam,

<br><br>Many thanks for the follow up.  It addresses many of my concerns.<br><br>On another front, is it possible that anyone out there has developed the<br>particular<br>files to do this datum translation (i.e., geographic NAD83 to LCC

<br>datumless--sphere;<br>LCC NAD83 to geographic datumless--sphere)?  I suspect not, but it never<br>hurts to ask.<br><br>Can anyone point me to some good literature (paper or web-based with<br>web-based<br>preferred) on how to develop three-point and seven-point datum

<br>translations that I<br>can then plug into PROJ4?  Anything that I may develop I will be happy<br>to share<br>with the community.<br><br>Regards,<br><br>Jim Wilkinson<br><br>Frank Warmerdam wrote:<br><br>> James G Wilkinson wrote:

<br>><br>>> A situation has occurred that I am having difficulty getting my hands<br>>> around.<br>>> I have a single point with the following attributes:<br>>><br>>> lon = -115.10101<br>>> lat = 

36.15526<br>>> PROJ4 entry in SPATIAL_REF_SYS<br>>> +proj=longlat +ellps=GRS80 +datum=NAD83 +no_defs<br>>><br>>> In ARC/Info, when I project this point to the following:<br>>><br>>> projection: LAMBERT

<br>>> datum: NONE<br>>> units: meters<br>>> spheroid: sphere<br>>> 1st parallel: 33.0<br>>> 2nd parallel: 45.0<br>>> central meridian: -118.0<br>>> latitude of origin: 37.0<br>>> false easting: 

0.0<br>>> false northing: 0.0<br>>><br>>> I get back a point at  y-coord = -89372 m.<br>>><br>>> The above projection is entered into my SPATIAL_REF_SYS table<br>>> with the following PROJ4 parameters (SRID=888889):

<br>>> +proj=lcc +lat_1=33.0 +lat_2=45.0 +lat_0=37.0 +lon_0=-118.0 +x_0=0<br>>> +y_0=0 +a=6370997.0 +b=6370997.0 +units=m<br>>><br>>> Now with PostgreSQL v.7.4.8, PostGIS v.0.9.1, and PROJ4 v.4.4.7

,<br>>> this same point using SELECT Y ( TRANSFORM (the_geom, 888889) ) FROM<br>>> spatial.pt_test<br>>> returns -89074 (which is close enough for my work.....for now).<br>>><br>>> When I migrated to PostgreSQL 

v.8.1.3, PostGIS v.1.1.2, and PROJ4<br>>> v.4.4.9,<br>>> my SELECT returns -109641 (!! and this is what baffles me !!).<br>>><br>>> Any help that folks can lend will be greatly appreciated.<br>>

> > Jim, > > I believe the problem is that PROJ.4 is attempting to correct for the > difference in ellipsoid when it should really just treat your NAD83 > lat/long location as a value on the defined sphere (without conversion).

<br>><br>> What it is doing now is taking the location in NAD83, converting it to<br>> a geocentric x/y/z relative to the center of the earth.  Then converting<br>> that location onto the sphere.  Because of the substantial difference

<br>> between GRS80 ellipsoid and a sphere this produces a large apparent<br>> shift in lat/long between the two earth models.<br>><br>> My suggestion is that if you want to treat your lat/long coordinate as<br>

> being on the sphere, then use the input coordinate system:<br>>   +proj=latlong +a=6370997.0 +b=6370997.0<br>><br>> I will add that the PROJ.4 behavior, while seemingly "proper", has been<br>> no end of problem and I will likely alter it in some future version.  It

<br>> seems the common practice in the GIS industry is to do what ESRI does.<br>> If not provided with datum shift information to just ignore the<br>> difference<br>> in ellipsoids.<br>><br>> Best regards,

<br><br>_______________________________________________<br>postgis-users mailing list<br><a href="mailto:postgis-users@postgis.refractions.net">postgis-users@postgis.refractions.net</a><br><a href="http://postgis.refractions.net/mailman/listinfo/postgis-users">

http://postgis.refractions.net/mailman/listinfo/postgis-users</a><br></blockquote></div><br><br clear="all"><br>-- <br>René F. Viáncos S.<br>Director de Geomática y TIC<br>Unidad de Ambiente, Región Territorio y Espacio, UNARTE

<br>Facultad de Ciencias<br>Universidad de Chile<br>Tel (56-2) 632 62 09<br>Cel (56 9) 933 72 66<br><a href="mailto:rviancos@uchile.cl">rviancos@uchile.cl</a><br><a href="mailto:rviancos@gmail.com">rviancos@gmail.com</a>