Testing Coordinate Transformation in PostGIS

I spend a lot of time thinking about which projection to use for a certain project,. Sometimes picking the perfect projection bites me when things change (such as the area of interest expands). In my recent case I have been avoiding reprojecting my data because I haven’t done this before in PostGIS and I’ve been convinced that this will be a major undertaking. After a minute’s worth of research I realize I may have been dead wrong. In fact, the process seems so easy and quick that I am leery. So here it is, my attempt at convincing myself that reprojecting in PostGIS is an inconsequential task.

First I will create some test data using a table in UTM 17 N meters (EPSG:26917).

PointsToReproject

Let’s prove this worked:

SELECT 'Table_SRID' as INFO,find_srid('public', 'ProjTest', 'geom')
UNION ALL
SELECT st_AsText(geom), st_SRID(geom) FROM "ProjTest";

TableInfoBefore

Check!

Next we’ll do the transformation to US National Atlas Equal Area

ALTER TABLE "ProjTest"
ALTER COLUMN geom TYPE geometry(Point,2163)
USING ST_Transform(geom,2163);

Then check the results:
TableInfoAfter

Check!

Finally let’s display it in QGIS:

First I’ll leave the project projection in QGIS to UTM 17N, looks as expected.
DisplayAfter26917

Just to double check, I’ll change the project projection to 2163. Surprisingly things look crazy. Notice that the offset is not predictable, meaning that it might not be a datum shift.
DisplayAfter2163

Interestingly, projecting “on the fly” to other coordinate systems does not reproduce this offset. So, perhaps it is the Google Physical web service that I am using that is having problems with reprojecting to US National Atlas Equal Area and not the reproduction of the points. To test that, we will reproject one of the coordinates back to UTM from US National Atlas Equal Area using gdal.


gdaltransform -s_srs EPSG:2163 -t_srs EPSG:26917
1480246.71738073 -990001.583001807

Results:
Original = 250230.930254224, 3854029.95272766
GDAL Transformation = 250230.930251404, 3854029.95272553

Looks like the reprojection was successful, just need to figure out the problems with reprojecting the base map.

In summary, reprojecting data in PostGreSql is trivial. The work can be done in place and scripts are very succinct. I am extremely pleased with this implementation, thanks guys!

Leave a Reply

Your email address will not be published.

− 1 = 1