[postgis-tickets] [PostGIS] #5437: ST_Intersects polygon geography empty results

[PostGIS]

#5437: ST_Intersects polygon geography empty results
-----------------------------------+---------------------------
  Reporter:  royjacksonalertmedia  |      Owner:  pramsey
      Type:  defect                |     Status:  new
  Priority:  medium                |  Milestone:  PostGIS 3.3.4
 Component:  postgis               |    Version:  3.3.x
Resolution:                        |   Keywords:
-----------------------------------+---------------------------
Description changed by royjacksonalertmedia:

Old description:

> I have a geography dataset containing various (valid) unioned polygons,
> including the United States (incl Alaska and Hawaii).
>
> I use a materialized view with a spatial index to st_union the shapes for
> a given id, which include one or many states, countries, circles.
>
> I published this view as a feature service with arcgis enterprise 11.1,
> and noticed a feature tile was returning no data.  Specifically, a large
> chunk of alaska is missing.
>
> I enabled logging in the database to figure out the query being used:
>
> select objectid,shape from maps_db.public.threat_polygon_shapes_layer
> where (id='faee6fc1d8f2f583577f85c54da3900b') AND
> (public.ST_Intersects(shape,public.ST_setSRID(public.ST_GeogFromWKB('\x01030000000100000005000000386094f6ff7f66c000415442d9a05040386094f6ff7f66c0047fa5b145435540fcffffffff7f56c0047fa5b145435540fcffffffff7f56c000415442d9a05040386094f6ff7f66c000415442d9a05040'),'4326'))
> = 't') ORDER BY objectid ASC
>
> The threat polygon here (alaska + US) is st_valid. The srids are the
> same.  The extent polygon is touching the northern end of the globe.
>
> I can view the two shapes together in pgadmin, and they overlap.  I can
> also cast as geometry and union them together.
>
> If I cast the shapes as geometry in the ST_Intersects query, the query
> returns expected results.
>
> select objectid,shape from maps_db.public.threat_polygon_shapes_layer
> where (id='faee6fc1d8f2f583577f85c54da3900b') AND
> (public.ST_Intersects(shape::geometry,public.ST_setSRID(public.ST_GeogFromWKB('\x01030000000100000005000000386094f6ff7f66c000415442d9a05040386094f6ff7f66c0047fa5b145435540fcffffffff7f56c0047fa5b145435540fcffffffff7f56c000415442d9a05040386094f6ff7f66c000415442d9a05040'),'4326')::geometry)
> = 't') ORDER BY objectid ASC
>
> The binary representation of the US + Alaska polygon is found in an
> attached file.
>
> PostGIS version details:
>
> POSTGIS="3.3.2 0" [EXTENSION] PGSQL="140" GEOS="3.10.3-CAPI-1.16.1"
> PROJ="9.1.0" LIBXML="2.9.9" LIBJSON="0.12.99" LIBPROTOBUF="1.3.0"
> WAGYU="0.5.0 (Internal)"
>
> Same results with:
> POSTGIS="3.3.3 2355e8e" [EXTENSION] PGSQL="140" GEOS="3.9.0-CAPI-1.16.2"
> PROJ="7.2.1" LIBXML="2.9.10" LIBJSON="0.15" LIBPROTOBUF="1.3.3"
> WAGYU="0.5.0 (Internal)"

New description:

 I have a geography dataset containing various (valid) unioned polygons,
 including the United States (incl Alaska and Hawaii).

 I use a materialized view with a spatial index to st_union the shapes for
 a given id, which include one or many states, countries, circles.

 I published this view as a feature service with arcgis enterprise 11.1,
 and noticed a feature tile was returning no data.  Specifically, a large
 chunk of alaska is missing.

 I enabled logging in the database to figure out the query being used:

 select objectid,shape from maps_db.public.threat_polygon_shapes_layer
 where (id='faee6fc1d8f2f583577f85c54da3900b') AND
 (public.ST_Intersects(shape,public.ST_setSRID(public.ST_GeogFromWKB('\x01030000000100000005000000386094f6ff7f66c000415442d9a05040386094f6ff7f66c0047fa5b145435540fcffffffff7f56c0047fa5b145435540fcffffffff7f56c000415442d9a05040386094f6ff7f66c000415442d9a05040'),'4326'))
 = 't') ORDER BY objectid ASC

 The threat polygon here (alaska + US) is st_valid. The srids are the same.
 The extent polygon is touching the northern end of the globe.

 I can view the two shapes together in pgadmin, and they overlap.  I can
 also cast as geometry and union them together.

 If I cast the shapes as geometry in the ST_Intersects query, the query
 returns expected results.

 select objectid,shape from maps_db.public.threat_polygon_shapes_layer
 where (id='faee6fc1d8f2f583577f85c54da3900b') AND
 (public.ST_Intersects(shape::geometry,public.ST_setSRID(public.ST_GeogFromWKB('\x01030000000100000005000000386094f6ff7f66c000415442d9a05040386094f6ff7f66c0047fa5b145435540fcffffffff7f56c0047fa5b145435540fcffffffff7f56c000415442d9a05040386094f6ff7f66c000415442d9a05040'),'4326')::geometry)
 = 't') ORDER BY objectid ASC

 The binary representation of the US + Alaska polygon is found in an
 attached file.

 PostgreSQL 14.7 on aarch64-unknown-linux-gnu, compiled by aarch64-unknown-
 linux-gnu-gcc (GCC) 9.5.0, 64-bit

 PostGIS version details:

 POSTGIS="3.3.2 0" [EXTENSION] PGSQL="140" GEOS="3.10.3-CAPI-1.16.1"
 PROJ="9.1.0" LIBXML="2.9.9" LIBJSON="0.12.99" LIBPROTOBUF="1.3.0"
 WAGYU="0.5.0 (Internal)"

 Same results with:
 POSTGIS="3.3.3 2355e8e" [EXTENSION] PGSQL="140" GEOS="3.9.0-CAPI-1.16.2"
 PROJ="7.2.1" LIBXML="2.9.10" LIBJSON="0.15" LIBPROTOBUF="1.3.3"
 WAGYU="0.5.0 (Internal)"

--
-- 
Ticket URL: <https://trac.osgeo.org/postgis/ticket/5437#comment:1>
PostGIS <http://trac.osgeo.org/postgis/>
The PostGIS Trac is used for bug, enhancement & task tracking, a user and developer wiki, and a view into the subversion code repository of PostGIS project.