Opened 10 years ago

Closed 10 years ago

Last modified 10 years ago

#2638 closed defect (fixed)

geography ST_Intersects bugginess with Polygon/multilinestring M

Reported by: robe Owned by: pramsey
Priority: medium Milestone: PostGIS 2.0.5
Component: postgis Version: 2.1.x
Keywords: Cc:

Description

SELECT ST_Intersects(geog1, geog2) As g1_in_g2, ST_Intersects(geog2,geog1) As g2_in_g1, _ST_DistanceUncached(geog1,geog2) As dist_unc
 , _ST_DistanceTree(geog1,geog2) As dtree_g1_g2 
 , _ST_DistanceTree(geog2,geog1) As dtree_g2_g1
FROM ST_GeogFromText('MULTILINESTRING M ((-71.0821 50.3036 1,50 -71 1),(-56.0821 50.3036 1,50 -56 1),(-41.0821 50.3036 1,50 -41 1),(-26.0821 50.3036 1,50 -26 1),(-11.0821 50.3036 1,50 -11 1),(3.9179 50.3036 1,50 4 1),(18.9179 50.3036 1,50 19 1),(33.9179 50.3036 1,50 34 1),(48.9179 50.3036 1,50 49 1),(-71.0821 50.3036 2,50 -71 2),(-56.0821 50.3036 2,50 -56 2),(-41.0821 50.3036 2,50 -41 2),(-26.0821 50.3036 2,50 -26 2),(-11.0821 50.3036 2,50 -11 2),(3.9179 50.3036 2,50 4 2),(18.9179 50.3036 2,50 19 2),(33.9179 50.3036 2,50 34 2),(48.9179 50.3036 2,50 49 2))') As geog1,
ST_GeogFromText('MULTIPOLYGON M (((0 0 2,10 0 1,10 10 -2,0 10 -5,0 0 -5),(5 5 6,7 5 6,7 7 6,5 7 10,5 5 -2)))') As geog2;

POSTGIS="2.1.2dev r12220" GEOS="3.4.2-CAPI-1.8.2 r0" PROJ="Rel. 4.8.0, 6 March 2012" GDAL="GDAL 1.10.0, released 2013/04/24" LIBXML="2.7.8" LIBJSON="UNKNOWN" RASTER PostgreSQL 9.2.4, compiled by Visual C++ build 1600, 32-bit

 g1_in_g2 | g2_in_g1 | dist_unc | dtree_g1_g2 | dtree_g2_g1
----------+----------+----------+-------------+-------------
 t        | f        |        0 |           0 |           0
(1 row)

In PostGIS 2.0 POSTGIS="2.0.4" GEOS="3.4.2-CAPI-1.8.2 r0" PROJ="Rel. 4.8.0, 6 March 2012" GDAL="GDAL 1.10.0, released 2013/04/24" LIBXML="2.7.8" LIBJSON="UNKNOWN" (core procs from "2.0.4" need upgrade) RASTER (raster procs from "2.0.4" need upgrade) PostgreSQL 9.2.4, compiled by Visual C++ build 1600, 32-bit

SELECT ST_Intersects(geog1, geog2) As g1_in_g2, ST_Intersects(geog2,geog1) As g2_in_g1
 , ST_Distance(geog1,geog2) As dist_g1_g2, ST_Distance(geog2,geog1) As dist_g2_g1
FROM ST_GeogFromText('MULTILINESTRING M ((-71.0821 50.3036 1,50 -71 1),(-56.0821 50.3036 1,50 -56 1),(-41.0821 50.3036 1,50 -41 1),(-26.0821 50.3036 1,50 -26 1),(-11.0821 50.3036 1,50 -11 1),(3.9179 50.3036 1,50 4 1),(18.9179 50.3036 1,50 19 1),(33.9179 50.3036 1,50 34 1),(48.9179 50.3036 1,50 49 1),(-71.0821 50.3036 2,50 -71 2),(-56.0821 50.3036 2,50 -56 2),(-41.0821 50.3036 2,50 -41 2),(-26.0821 50.3036 2,50 -26 2),(-11.0821 50.3036 2,50 -11 2),(3.9179 50.3036 2,50 4 2),(18.9179 50.3036 2,50 19 2),(33.9179 50.3036 2,50 34 2),(48.9179 50.3036 2,50 49 2))') As geog1,
ST_GeogFromText('MULTIPOLYGON M (((0 0 2,10 0 1,10 10 -2,0 10 -5,0 0 -5),(5 5 6,7 5 6,7 7 6,5 7 10,5 5 -2)))') As geog2;

 g1_in_g2 | g2_in_g1 |    dist_g1_g2    |    dist_g2_g1
----------+----------+------------------+------------------
 t        | f        | 184534.350787176 | 184534.350787176

There seems to be a bug in both 2.1.2 and 2.0.4

For 2.0.4 I was able to get the ST_Intersects to flip flop (f, f) and then when I ran it again one became a true as you see above

Anyway clearly the answer is flawed however you look at it.

I still need to check if the issue shows with plain 2D for this. I didn't see it in the regression output but it may have just been dum luck

Change History (10)

comment:1 by robe, 10 years ago

Just tested by collapsing to 2D.

SELECT ST_Intersects(geog1, geog2) As g1_in_g2, ST_Intersects(geog2,geog1) As g2_in_g1
 , ST_Distance(geog1,geog2) As dist_g1_g2, ST_Distance(geog2,geog1) As dist_g2_g1
FROM ST_GeogFromText('MULTILINESTRING((-71.0821 50.3036,50 -71),(-56.0821 50.3036,50 -56),(-41.0821 50.3036,50 -41),(-26.0821 50.3036,50 -26),(-11.0821 50.3036,50 -11),(3.9179 50.3036,50 4),(18.9179 50.3036,50 19),(33.9179 50.3036,50 34),(48.9179 50.3036,50 49),(-71.0821 50.3036,50 -71),(-56.0821 50.3036,50 -56),(-41.0821 50.3036,50 -41),(-26.0821 50.3036,50 -26),(-11.0821 50.3036,50 -11),(3.9179 50.3036,50 4),(18.9179 50.3036,50 19),(33.9179 50.3036,50 34),(48.9179 50.3036,50 49))') As geog1,
ST_GeogFromText('MULTIPOLYGON(((0 0,10 0,10 10,0 10,0 0),(5 5,7 5,7 7,5 7,5 5)))') As geog2;

The answers agree and intersect true always.

 g1_in_g2 | g2_in_g1 | dist_g1_g2 | dist_g2_g1
----------+----------+------------+------------
 t        | t        |          0 |          0
(1 row)

So somehow it seems the M is getting in the way when it really shouldn't be considered at all and is causing some sort of instability.

comment:2 by robe, 10 years ago

definitely some instability here — run on postgis 2.2

PostgreSQL 9.2.4, compiled by Visual C++ build 1600, 32-bit POSTGIS="2.2.0dev r12228" GEOS="3.4.2-CAPI-1.8.2 r0" PROJ="Rel. 4.8.0, 6 March 2012" GDAL="GDAL 1.10.0, released 2013/04/24" LIBXML="2.7.8" LIBJSON="UNKNOWN" RASTER
SELECT ST_Intersects(geog1, geog2) As g1_in_g2, ST_Intersects(geog2,geog1) As g2_in_g1, _ST_DistanceUncached(geog1,geog2) As dist_unc
 , _ST_DistanceTree(geog1,geog2) As dtree_g1_g2 
 , _ST_DistanceTree(geog2,geog1) As dtree_g2_g1
FROM ST_GeogFromText('MULTILINESTRING M ((-71.0821 50.3036 1,50 -71 1),(-56.0821 50.3036 1,50 -56 1),(-41.0821 50.3036 1,50 -41 1),(-26.0821 50.3036 1,50 -26 1),(-11.0821 50.3036 1,50 -11 1),(3.9179 50.3036 1,50 4 1),(18.9179 50.3036 1,50 19 1),(33.9179 50.3036 1,50 34 1),(48.9179 50.3036 1,50 49 1),(-71.0821 50.3036 2,50 -71 2),(-56.0821 50.3036 2,50 -56 2),(-41.0821 50.3036 2,50 -41 2),(-26.0821 50.3036 2,50 -26 2),(-11.0821 50.3036 2,50 -11 2),(3.9179 50.3036 2,50 4 2),(18.9179 50.3036 2,50 19 2),(33.9179 50.3036 2,50 34 2),(48.9179 50.3036 2,50 49 2))') As geog1,
ST_GeogFromText('MULTIPOLYGON M (((0 0 2,10 0 1,10 10 -2,0 10 -5,0 0 -5),(5 5 6,7 5 6,7 7 6,5 7 10,5 5 -2)))') As geog2;

I got the same answer as 2.1.2dev the first run and on second run got this magical answer

 g1_in_g2 | g2_in_g1 | dist_unc             | dtree_g1_g2 | dtree_g2_g1
----------+----------+----------------------+-------------+-------------
 t        | f        | 184534.350787176;0;0 |           0 |           0

which suggests the uncached distance may be unstable as well for this M thingy

comment:3 by robe, 10 years ago

2.1.0 also exhibits same flip flopping PostgreSQL 9.2.4, compiled by Visual C++ build 1600, 32-bit POSTGIS="2.1.0 r11822" GEOS="3.4.2-CAPI-1.8.2 r0" PROJ="Rel. 4.8.0, 6 March 2012" GDAL="GDAL 1.10.0, released 2013/04/24" LIBXML="2.7.8" LIBJSON="UNKNOWN" RASTER

— first run gave me

 g1_in_g2 | g2_in_g1 | dist_unc | dtree_g1_g2 | dtree_g2_g1
----------+----------+----------+-------------+-------------
 t        | f        |        0 |           0 |           0
(1 row)

second run gave me

 g1_in_g2 | g2_in_g1 | dist_unc         | dtree_g1_g2 | dtree_g2_g1
----------+----------+------------------+-------------+-------------
 t        | f        | 184534.350787176 |           0 |           0
(1 row)

comment:4 by pramsey, 10 years ago

Distance *uncached* is flipping on you? This just keeps getting harder and harder to understand.

comment:5 by robe, 10 years ago

only using M as far as I can tell (with 2D and 2DZ haven't witnessed the issue). Why it even cares about M seems pretty puzzling to me.

comment:6 by robe, 10 years ago

I just realized something — I wasn't testing dist_unc_g2_g1. That one seems to be always consistently the larger number though the g1_g2 changes back and forth between 0 and non-zero.

SELECT ST_Intersects(geog1, geog2) As g1_in_g2, ST_Intersects(geog2,geog1) As g2_in_g1
 , _ST_DistanceUncached(geog1,geog2) As dist_unc_g1_g2
 , _ST_DistanceUncached(geog1,geog2) As dist_unc_g2_g1
 , _ST_DistanceTree(geog1,geog2) As dtree_g1_g2 
 , _ST_DistanceTree(geog2,geog1) As dtree_g2_g1
FROM ST_GeogFromText('MULTILINESTRING M ((-71.0821 50.3036 1,50 -71 1),(-56.0821 50.3036 1,50 -56 1),(-41.0821 50.3036 1,50 -41 1),(-26.0821 50.3036 1,50 -26 1),(-11.0821 50.3036 1,50 -11 1),(3.9179 50.3036 1,50 4 1),(18.9179 50.3036 1,50 19 1),(33.9179 50.3036 1,50 34 1),(48.9179 50.3036 1,50 49 1),(-71.0821 50.3036 2,50 -71 2),(-56.0821 50.3036 2,50 -56 2),(-41.0821 50.3036 2,50 -41 2),(-26.0821 50.3036 2,50 -26 2),(-11.0821 50.3036 2,50 -11 2),(3.9179 50.3036 2,50 4 2),(18.9179 50.3036 2,50 19 2),(33.9179 50.3036 2,50 34 2),(48.9179 50.3036 2,50 49 2))') As geog1,
ST_GeogFromText('MULTIPOLYGON M (((0 0 2,10 0 1,10 10 -2,0 10 -5,0 0 -5),(5 5 6,7 5 6,7 7 6,5 7 10,5 5 -2)))') As geog2;

g1_in_g2 | g2_in_g1 | dist_unc_g1_g2 |  dist_unc_g2_g1  | dtree_g1_g2 | dtree_g_g1
---------+----------+----------------+------------------+-------------+------------
t        | f        |              0 | 184534.350787176 |           0 |

Is it possible you have some assumption in your 2.1 code that ST_Distance(g1,g2) = ST_Distance(g2,g1) (so you are arbitrarily choosing one to check )

That would explain the flipflopping. Note the 2.0.4 is wrong but it is always consistently wrong (for distance where as ST_Intersects it seems to be treating g1 xg2 different from g2 x g1 (which is BIZARRE when the distance is always the higher) and ST_Intersects is just short circuit for

SELECT $1 && $2 AND _ST_Distance($1, $2, 0.0, false) < 0.00001

so how on earth can that ever be true when the distance it returns is never 0. Might be having same behavior but just seeing it less often in 2.0.4 because of some other counter balancing.

comment:7 by robe, 10 years ago

I just realized that ST_Intersects uses the sphere distance rather than the default spheroid distance and we already know from #2634 , #2636 that using sphere distance is unstable. So this may be a different manifestation of the same bug looked at a different angle. You said 2.1 now reads from spatial_ref_sys instead of a hard-coded sphere. I wonder if that is somehow adding to the perceived instability.

anyrate it sure seems to be tripped up by the M coordinate when it shouldn't care.

comment:8 by pramsey, 10 years ago

I'm seeing the same results as you, so there's some hope. I may home some plane time tomorrow to pick at this scab.

comment:9 by pramsey, 10 years ago

Resolution: fixed
Status: newclosed

Got a free hour on the plane home, found it, and yes, it's longstanding, a small error in handling gbox overlap tests. Fixed in 2.0 r12298

In 2.1, r12296

in trunk r12297

comment:10 by robe, 10 years ago

You should just stay up in the air you seem to be much more productive :)

Note: See TracTickets for help on using tickets.