5

I have plain data in a postgresql table which looks like this:

╔════╦═══════╦═════════════╦══════════════╦══════════╦═════════════════╗
║ id ║ index ║   easting   ║   northing   ║ utm_zone ║ utm_zone_letter ║
╠════╬═══════╬═════════════╬══════════════╬══════════╬═════════════════╣
║  1 ║   1   ║ 233410.0000 ║ 1024701.0000 ║       38 ║ N               ║
║  1 ║   2   ║ 213310.0000 ║ 1167201.0000 ║       36 ║ N               ║
║  2 ║   1   ║ 213310.0000 ║ 1167201.0000 ║       31 ║ B               ║
║  2 ║   2   ║ 213310.0000 ║ 1167201.0000 ║       31 ║ B               ║
╚════╩═══════╩═════════════╩══════════════╩══════════╩═════════════════╝

Note: Coordinates are not real, but typed by hand

As you can see these are UTM coordinates but stored in plain db type fields.
I want to convert those last 4 columns to one of Geometry type for PostGis.
To something like this:

╔════╦══════════════════════╗
║ id ║       Geometry       ║
╠════╬══════════════════════╣
║  1 ║ 02347237427342342347 ║
║  2 ║ 27584872345646325863 ║
╚════╩══════════════════════╝

I used query similar to this:


--creates polygon with wrong SRID :(
SELECT z.id, ST_Polygon(z.geom, 3395) geom FROM ( 
    --adds first point of line at the end to close it (requeird by ST_Polygon)
    SELECT l.id id, ST_AddPoint(g.geom, ST_StartPoint(g.geom)) geom FROM (
        -- groups points by ID and creates open line geometry from each group
        SELECT t.id id, ST_MakeLine(t.strPoints) geom FROM (
            --creates POINT(x,y) text from coordinates
            SELECT 
                g.id id, 
                g.index idx, 
                'POINT(' || g.easting || ' ' || g.northing || ')' strPoints
            FROM data g
            ORDER BY g.index
        ) t
        GROUP BY t.id
        HAVING count(t.strPoints) > 2
    ) g

) z
where ST_IsValid(z.geom)

But this does not use utm_zone and utm_zone_letter and points and lines are created without SRID value. Creating polygon from such values and specifing some SRID gives bad results. At same time because points and lines do not have SRID set I cannot use utmzone(geometry) DB function because it requires geometry to have SRID.

I looked at ST_Transform, ST_FromText and other functions but could not find one where I can specify UTM zone number and letter to create Geometry type value.
Creating point without specifying SRID gives me bad results - points fall on penguins in Antarctica which is wrong.

Can you tell me how to convert that data to correct Geometry values with correct SRID specified?

0

3 Answers 3

2

I think you should be able to use this function (from here):

 -- Function: utmzone(geometry)
 -- DROP FUNCTION utmzone(geometry);
 -- Usage: SELECT ST_Transform(the_geom, utmzone(ST_Centroid(the_geom)) )
    FROM sometable;

 CREATE OR REPLACE FUNCTION utmzone(geometry)
   RETURNS integer AS
 $BODY$
 DECLARE
     geomgeog geometry;
     zone int;
     pref int;

 BEGIN
     geomgeog:= ST_Transform($1,4326);

     IF (ST_Y(geomgeog))>0 THEN
        pref:=32600;
     ELSE
        pref:=32700;
     END IF;

     zone:=floor((ST_X(geomgeog)+180)/6)+1;

     RETURN zone+pref;
 END;
 $BODY$ LANGUAGE 'plpgsql' IMMUTABLE
   COST 100;
1
  • Yes I know about that function, but there's one problem: It requires that geometry passed to it must have SRID set. select utmzone(ST_PointFromText('POINT(' || g.easting || ' ' || northing || ')')) from gis g This fails with: ERROR: Input geometry has unknown (0) SRID In my case I guess SRID is defined by utm_zone & utm_zone_letter columns? Then I need to somehow use them to determine correct SRID for each record. Filtering spatial_ref_sys table using those two column values sometimes results in two SRID records. That is not good I guess... need some other way.
    – faskunji
    Commented Jan 12, 2015 at 14:47
1

You could use the information here to generate a SRID - basically add the zone number to 32600 (in the north) and 32700 in the south. This works if you can use WGS84 as the geographic coordinate reference system.

The utm_zone_letter identifies the latitude band. If the data is always in the northern hemisphere, then you could ignore it and use 32600 + utm_zone as the well-known ID. If the data is in both hemispheres, you'll have to check the letters to determine which hemisphere. 'CDEFGHJKLM' are in the southern hemisphere, and 'NPQRSTUVWX' are in the northern hemisphere. A and B are below 80 South and Y and Z are above 84 North.

2
  • Expanded your answer--revert/reword if you want!
    – mkennedy
    Commented Jan 12, 2015 at 19:09
  • Thank you @iant. The link and your previous comment helped me to solve the problem.
    – faskunji
    Commented Mar 19, 2015 at 10:38
1

You'll need to translate the UTM zone to an EPSG code, use ST_FromText to create a geometry object in that projection, and then use ST_Transform to transform your geometry object into whatever SRID you are using for the entire table.

I haven't tried it, but this link gives a description of how to programatically convert from latitude/longitude to UTM zone; I expect it wouldn't be too difficult to reverse the calculation.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.