[postgis-users] polygons of data type geography showing strange behaviour in combination with the date line

[Mike Toews]

Spatial logic with geography types works differently than geometry
types, and you'd need to pin your geometries on a globe with strings
to get the logic correct. It seems from your poly1 and poly2 that
these are interpreted to be rectangular, but unfortunately it is
really difficult to draw a rectangle on a sphere.

For instance, you might think it is 2 degrees tall (because the
coordinates are +/-1 degree N/S), however the geography type treats
all edges as geodesics, and not rhumb lines, so polygons generally
bulge out in each direction. Furthermore, lines that are >180 degrees
long cannot be pinned to a sphere, which is why you are perhaps
struggling with Poly 2.

Your workaround to add more vertices is a good one. This helps
approximate a geometry rectangle as geography rectangle, where the
geodesic lines have shorter distances and minimize warping. You can
densify the geometry points via ST_Segmentize(geometry,
max_segment_length), used here in one of your queries:

postgis=# SELECT ST_DWithin(
postgis(#     'POINT(0 0)'::geography,
postgis(#     ST_Segmentize(
postgis(#         'POLYGON((170 -1,-170 -1,-170 1,170 1,170 -1))'::geometry,
postgis(#         10.0)::geography, 0.0);
 st_dwithin
------------
 t
(1 row)

This example uses a maximum of 10 degrees distance on the geometry
before transferring to geography logic, which looks like:
http://www.gcmap.com/mapui?P=1S+170E-1S+160E-1S+150E-1S+140E-1S+130E-1S+120E-1S+110E-1S+100E-1S+90E-1S+80E-1S+70E-1S+60E-1S+50E-1S+40E-1S+30E-1S+20E-1S+10E-1S+0E-1S+10W-1S+20W-1S+30W-1S+40W-1S+50W-1S+60W-1S+70W-1S+80W-1S+90W-1S+100W-1S+110W-1S+120W-1S+130W-1S+140W-1S+150W-1S+160W-1S+170W-1N+170W-1N+160W-1N+150W-1N+140W-1N+130W-1N+120W-1N+110W-1N+100W-1N+90W-1N+80W-1N+70W-1N+60W-1N+50W-1N+40W-1N+30W-1N+20W-1N+10W-1N+0E-1N+10E-1N+20E-1N+30E-1N+40E-1N+50E-1N+60E-1N+70E-1N+80E-1N+90E-1N+100E-1N+110E-1N+120E-1N+130E-1N+140E-1N+150E-1N+160E-1N+170E-1S+170E&MS=wls&DU=km


On 11 July 2017 at 09:53,  <Christian.Strobl at dlr.de> wrote:
> Dear all,
>
> I define a polygon (Poly 1) which extends over the date line. Point 1
> clearly is inside Poly 1.
>
> geoDB=> SELECT ST_DWITHIN(ST_GeogFromText('POINT(180 0)'),
> ST_GeogFromText('POLYGON((170 -1, 190 -1, 190 1, 170 1, 170 -1))'), 0);
> NOTICE:  Coordinate values were coerced into range [-180 -90, 180 90] for
> GEOGRAPHY
>  st_dwithin
> ------------
>  t
> (1 row)
>
> If I have a closer look at this polygon it shows the representation below.
>
> geoDB=> SELECT ST_ASEWKT(ST_GeogFromText('POLYGON((170 -1, 190 -1, 190 1,
> 170 1, 170 -1))'));
> NOTICE:  Coordinate values were coerced into range [-180 -90, 180 90] for
> GEOGRAPHY
>                         st_asewkt
> ---------------------------------------------------------
>  SRID=4326;POLYGON((170 -1,-170 -1,-170 1,170 1,170 -1))
> (1 row)
>
>
>
>       -180             0               180
>
>         |              |                |
>         |              |                |
> 90 -----|--------------|----------------|----- 90
>         |              |                |
>         |              |                |
>         |              |                |
>         |              | Point2 (0/0)   | Point1 (180/0)
>         |  ____________|______________ ____
>  0 -----|-|____________*______________|_*__|-- 0
>         |              |                |
>         |              | Poly 2         | Poly 1
>         |              |                |
>         |              |                |
>         |              |                |
> -90-----|--------------|----------------|----- -90
>         |              |                |
>         |              |                |
>       -180             0               180
>
>
> If I look at the order of the coordinates it looks like Poly 2.
>
> geoDB=> SELECT ST_DWITHIN(ST_GeogFromText('POINT(0 0)'),
> ST_GeogFromText('POLYGON((170 -1,-170 -1,-170 1,170 1,170 -1))'),0);
>  st_dwithin
> ------------
>  f
> (1 row)
>
> The test with Point 2 shows that this assumption is wrong. Again Point 1 is
> inside.
>
> geoDB=> SELECT ST_DWITHIN(ST_GeogFromText('POINT(180 0)'),
> ST_GeogFromText('POLYGON((170 -1,-170 -1,-170 1,170 1,170 -1))'),0);
>  st_dwithin
> ------------
>  t
> (1 row)
>
> If I use the geometry data type instead of the geography data type Point 2
> is clearly inside Poly 2.
>
> geoDB=> SELECT ST_DWITHIN(ST_GeomFromText('POINT(0 0)'),
> ST_GeomFromText('POLYGON((170 -1,-170 -1,-170 1,170 1,170 -1))'),0);
>  st_dwithin
> ------------
>  t
> (1 row)
>
> And finally here is my question: How do I construct Poly 2 from these
> geographic points?
>
> geoDB=> SELECT ST_ASEWKT(ST_GeogFromText('POLYGON((170 -1,-170 -1,-170 1,170
> 1,170 -1))'));
>                         st_asewkt
> ---------------------------------------------------------
>  SRID=4326;POLYGON((170 -1,-170 -1,-170 1,170 1,170 -1))
>
>
> Best regards and thanks for your help
>
> Christian
>
>
> P.S. A possible workaround is the use of vertices between the corner points,
> e.g.
>
> SELECT ST_ASEWKT(ST_GeogFromText('POLYGON((170 -1, 0 -1, -170 -1, -170 1, 0
> 1, 170 1, 170 -1))'));
>
> geoDB=> SELECT ST_DWITHIN(ST_GeogFromText('POINT(0 0)'),
> ST_GeogFromText('POLYGON((170 -1, 0 -1, -170 -1, -170 1, 0 1, 170 1, 170
> -1))'),0);
>  st_dwithin
> ------------
>  t
> (1 row)
>
> For me this seems like a possible workaround, but not like a answer to the
> question above
>
>
> _______________________________________________
> postgis-users mailing list
> postgis-users at lists.osgeo.org
> https://lists.osgeo.org/mailman/listinfo/postgis-users