[postgis-tickets] [PostGIS] #4290: bug in Postgis geodetic code impacting distance calculation with geography multipolygon

[PostGIS]

#4290: bug in Postgis geodetic code impacting distance calculation with geography
multipolygon
--------------------------+---------------------------
 Reporter:  dmkaplan2000  |      Owner:  pramsey
     Type:  defect        |     Status:  new
 Priority:  medium        |  Milestone:  PostGIS 2.5.2
Component:  postgis       |    Version:  2.4.x
 Keywords:                |
--------------------------+---------------------------
 I have found a bizarre bug in PostGIS when doing geography-based distance
 calculations using ST_Distance when one of the two polygons is a
 multipolygon. Under certain cases,
 ST_Distance(geography(Polygon,4326),geography(MultiPolygon,4326)) returns
 0 (zero) when the true answer is not zero.

 The table data for a set of polygons that generate the error can be
 downloaded here:

 http://www.davidmkaplan.fr/bad_geography_dist_calcs.csv

 Once you have that data, the code to generate the error is:

 ------------------------------

 CREATE TEMP TABLE tt (gid int PRIMARY KEY, geom geometry);

 \copy tt FROM 'bad_geography_dist_calcs.csv' WITH (FORMAT 'csv',HEADER
 TRUE)

 -- Incorrect dist_geography and dist_geography_tree for t1.gid=1 and
 t2.gid=3
 SELECT t1.gid AS gid1, t2.gid AS gid2,
        ST_Distance(t1.geom,t2.geom) AS dist_lonlat,
        ST_Distance(ST_Transform(t1.geom,26918),
                    ST_Transform(t2.geom,26918)) AS dist_utm,
        ST_Distance(t1.geom::geography,t2.geom::geography) AS
 dist_geography,
        _ST_DistanceTree(t1.geom::geography,t2.geom::geography) AS
 dist_geography_tree,
        _ST_DistanceUncached(t1.geom::geography,t2.geom::geography) AS
 dist_geography_uncached,
 ST_Distance((ST_Dump(t1.geom)).geom::geography,(ST_Dump(t2.geom)).geom::geography)
 AS dist_geography_dump
 FROM tt t1 JOIN tt t2 ON t1.gid<t2.gid
 WHERE t1.gid=1
 ;

 -- No error for postgis 2.4.3, but error for 2.5.0
 SELECT t1.gid AS gid1, t2.gid AS gid2,
        ST_Distance(t1.geom,t2.geom) AS dist_lonlat,
        ST_Distance(ST_Transform(t1.geom,26918),
                    ST_Transform(t2.geom,26918)) AS dist_utm,
        ST_Distance(t1.geom::geography,t2.geom::geography) AS
 dist_geography,
        _ST_DistanceTree(t1.geom::geography,t2.geom::geography) AS
 dist_geography_tree,
        _ST_DistanceUncached(t1.geom::geography,t2.geom::geography) AS
 dist_geography_uncached,
 ST_Distance((ST_Dump(t1.geom)).geom::geography,(ST_Dump(t2.geom)).geom::geography)
 AS dist_geography_dump
 FROM tt t1 JOIN tt t2 ON t1.gid<t2.gid
 WHERE t1.gid=1 AND t2.gid=3
 ;

 ----------------------------

 The specifics of my installation are:

 POSTGIS="2.4.3 r16312" PGSQL="100" GEOS="3.7.0-CAPI-1.11.0 673b9939"
 PROJ="Rel. 4.9.3, 15 August 2016" GDAL="GDAL 2.3.2, released 2018/09/21"
 LIBXML="2.9.4" LIBJSON="0.12.1" LIBPROTOBUF="1.2.1" RASTER

 psql 10.6 (Ubuntu 10.6-0ubuntu0.18.04.1))

 Ubuntu 18.04.1 LTS

 I have also tested this with postgis 2.5.0 r16836 and there the incorrect
 0 is in both queries.

-- 
Ticket URL: <https://trac.osgeo.org/postgis/ticket/4290>
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.