5

I have a PostGIS 2.2 table that stores polygons as geography type. The table contains several polygons that cross the antimeridian. I would like to split these polygons at the antimeridian and convert the table's topology from geography to geometry.

4

In the example below, thepolygon table is called poly_coverage and the output table is fixed_polygons.

Select the polygons that cross the anti-meridian, cast to geom and shift to 0-360, then split them. Select the split polys with any longitude > 180, and shift back 360 deg

Create table beams_fixed_dump4 as Select ST_Translate(foo2.thegeom, -360, 0) as newgeom from
(Select (ST_DUMP(foo.geom)).geom as thegeom from 
(Select ST_Split(shift_lon.geom, line2) as geom 
    from ST_SetSRID(ST_MakeLine(ST_MakePoint(180, -90),ST_MakePoint(180, 90)), 4326) As line2, 
    (Select ST_ShiftLongitude(geog::geometry) as geom from 
        (SELECT ST_Intersects(poly_coverage.geog, line) as inter, poly_coverage.* from poly_coverage, ST_GeogFromText('LINESTRING(180 -90, 180 0, 180 90)') As line) as tmp where inter = TRUE) as shift_lon) as foo) as foo2 where st_xmax(foo2.thegeom) > 180
Union all

Select the polygon that cross the anti-meridian, cast to geom and shift to 0-360, then split them. Select the split polys with all longitude <= 180.

Select newgeom from 
(Select (ST_DUMP(foo.geom)).geom as newgeom from 
(Select ST_Split(shift_lon.geom, line2) as geom from 
    ST_SetSRID(ST_MakeLine(ST_MakePoint(180, -90),ST_MakePoint(180, 90)), 4326) As line2, 
    (Select ST_ShiftLongitude(geog::geometry) as geom from 
        (SELECT ST_Intersects(poly_coverage.geog, line) as inter, poly_coverage.* 
            from poly_coverage, ST_GeogFromText('LINESTRING(180 -90, 180 0, 180 90)') As line) as tmp where inter = TRUE) as shift_lon) as foo)as foo2 where st_xmax(newgeom) <= 180
Union all

Select the polygons that don't cross the anti-meridian and cast to geom.

SELECT foo.geog::geometry from 
(Select ST_Intersects(poly_coverage.geog, line) as inter, poly_coverage.* from 
poly_coverage, ST_GeogFromText('LINESTRING(180 -90, 180 0, 180 90)') As line) as foo where inter = FALSE;
1
  • 4
    This query does not work correctly for polygons that surround one of the poles. Also, I am uncertain if it will process a line correctly with a point exactly at either pole. However, PostGIS seems to correctly convert such polygons from geography to geometry. You can simply exclude these polygons from the first two parts of the query and process them in the third.
    – jgm_GIS
    Dec 23 '16 at 19:45
0

The 180th Meridian was a handy reference for me understanding this problem. I think working with Common Table Expressions (CTE) may give you better performance than all the temporary tables required for SELECT ... FROM (SELECT ... FROM ...). Here is a query I came up with to split multipolygons along the antimeridian/dateline for the table poly_coverage with geometry geom.

WITH t AS (
  SELECT 
   ST_WrapX(
    (
      ST_Dump(ST_Split(ST_ShiftLongitude(geom),
              ST_GeomFromText('LINESTRING(180 -90, 180 90)', 4326)))
    ).geom,
    180, -360) AS geom
  FROM poly_coverage
)
UPDATE poly_coverage
SET geom=ST_Multi(ST_Union(t.geom)))
FROM t;

Working inside-out:

  1. ST_Shift_Longitude maps coordinates to 0-360.
  2. ST_Split cuts along the antimeridian so you now have a GEOMETRYCOLLECTION of two POLYGONs
  3. ST_Dump(...).geom yields the polygons as rows
  4. ST_WrapX maps longitudes back to the range (-180, 180)
  5. ST_Multi(ST_Union(t.geom)) takes the results from the CTE and emits a single MultiPolygon.

This does not address converting your table from geography to geometry. Maybe you could make a new table and update the CTE to insert while casting with geom::geography and/or geog::geometry?

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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