Wednesday, April 22, 2009

Using PostGIS to convert SVY21 easting/northing to latitude/longitude coordinates

Besides commercial software like GeoMedia, Global Mapper, FME, ArcGis, etc. there are also free open source tools that can be used to convert SVY21 easting/northing coordinates to geographical latitudes/longitude coordinates. A great open source software is PostGIS and in this posting I will show how to use it without any programming or customization to convert SVY21 coordinates in CSV (comma separated values) files to geographical coordinates.

Note: Prior to doing the conversion, the SVY21 projection parameters must be added to the PostGIS spatial reference system table as shown in the previous posting here

The source CSV file is shown in the screen shot below. The file has only three columns - a running number, the easting coordinate and the northing coordinate. 

  1. Open up the pgAdmin GUI by selecting Start | All Programs | PostgreSQL 8.3 | pgAdmin III.

    The pgAdmin III application appears.

  2. Connect to your PostgreSQL server and select your database e.g. svy21.

  3. Select Tools | Query Tool.

    The Query dialog box appears.

  4. Enter the following SQL command to create a database table to store the records from the CSV file.

    CREATE TABLE svy21coords   ( id SERIAL NOT NULL,   easting DOUBLE PRECISION,   northing DOUBLE PRECISION,
    CONSTRAINT id_pkey PRIMARY KEY (id)  ); 
  5. On the keyboard, press F5 to execute the command.

    The table svy21coords is created.

  6. Select Edit | Clear Window. Enter the following SQL command to load in the SVY21 records from the CSV file. Note: replace C:/temp/svy21coords.csv with the name of your file.

    COPY svy21coords (id, easting, northing)  FROM 'c:/temp/svy21coords.csv'   WITH CSV HEADER ; 
  7. Press F5 to execute the command.

    The records are loaded to the table.

  8. Select Edit | Clear Window. Enter the following SQL command to create a geometry field that is associated with the SVY21 coordinate system.

    SELECT AddGeometryColumn 
    ( 'public', 'svy21coords', 'the_geom', 93414, 'POINT', 2 ); 
    Note: 93414 is the SRID of the SVY21 coordinate system in the spatial_ref_sys table.

  9. Press F5 to execute the command.

    The geometry field is appended to the table.

  10. Select Edit | Clear Window. Enter the following SQL command.

    SELECT 'UPDATE svy21coords SET the_geom=GeomFromText(''POINT('||easting||' '||northing||')'',93414) WHERE id='||id||';'  FROM svy21coords; 
  11. Press F5 to execute the query.

    The query results are displayed in the bottom pane of the Query window.


  12. Select all the query results by clicking the first row and shift-clicking the last row. Press CTRL-C to copy the selected results to the Windows Clipboard.

  13. Click the top pane. Select Edit | Clear Window. Select Edit | Paste.

    The results are pasted into the top pane of the query window. Note the double quotation marks on every line. 
    "UPDATE svy21coords SET the_geom=GeomFromText('POINT(39815 35200)',93414) WHERE id=1;"
    "UPDATE svy21coords SET the_geom=GeomFromText('POINT(39800 35200)',93414) WHERE id=2;"
    "UPDATE svy21coords SET the_geom=GeomFromText('POINT(43200 35900)',93414) WHERE id=3;"


  14. Select Edit | Find and Replace.

    The Find and Replace dialog box appears.


  15. In the Find what field, enter ". Click Replace All. Close the Find and Replace dialog box.

    All the double quote marks are removed.
    UPDATE svy21coords1 SET the_geom=GeomFromText('POINT(39815 35200)',93414) WHERE id=1;
    UPDATE svy21coords1 SET the_geom=GeomFromText('POINT(39800 35200)',93414) WHERE id=2;
    UPDATE svy21coords1 SET the_geom=GeomFromText('POINT(43200 35900)',93414) WHERE id=3;


  16. Press F5 to execute the SQL command.

    The geometry field the_geom is updated with the SVY21 point coordinates.

  17. Select Edit | Clear Window. Enter the following SQL command to transform the SVY21 geometries to geographical latitude/longitude and write the records to a temporary table CSVOUTPUT.

    CREATE TABLE CSVOUTPUT  AS SELECT id,easting,northing, x(ST_Transform(the_geom,4326)),y(ST_Transform(the_geom,4326))  from svy21coords  Note: 4326 is the SRID of the geographical lat/lng coordinate system. 
  18. Press F5.

    The SVY21 points are transformed into latitude/longitude and written to the output table.

  19. Select Edit | Clear Window. Enter the following SQL command to save the converted coordinates to an output CSV file. Note: replace C:/temp/out/svy21coords.csv with your output file name.

    COPY CSVOUTPUT TO 'C:/TEMP/OUT/svy21coords.csv' WITH CSV HEADER
  20. Press F5.

    The output file is created as shown in the Excel screen shot below.

No comments: