ST_Intersects incorrect for MultiPoints

[dbaston]

Same root cause as #3433; the problem occurs when prepared geometry use kicks in. Example case:

WITH
a AS (SELECT ST_Collect(ST_MakePoint(random(), random())) as geom from generate_series(1, 10)),
b AS (SELECT (ST_Dump(geom)).geom FROM a)
SELECT count(*) FROM b INNER JOIN a ON ST_Intersects(a.geom, b.geom);

This should return 10 rows but only returns one. Seems to return the expected result for GEOS < 3.5, perhaps because PreparedPoint implementation isn't actually invoked.

[dbaston]

The issue is resolved with the proposed fix to ​https://trac.osgeo.org/geos/ticket/764, though I wonder if this should be handled at the PostGIS level too.

[pramsey]

I'm assigning to you, something for Paris perhaps. I agree we should probably have an #ifdef'ed fix in PostGIS we can tear out when we finally get a fixed version of GEOS as our baseline of support.

[dbaston]

Since all access to PreparedGeometry is through a cache, I just modified the cache code to avoid generating a prepared geometry around a Point or a MultiPoint. This feels like an odd place to make the fix, but it at least avoids littering the code with extra conditions every time we try to use a GEOSPrepared* function. I guess this could be #ifdef'd out allow prepared geometry creation for certain GEOS versions, but I'm not sure which released versions both (a) actually invoke PreparedPoint::intersects through the C API and (b) provide correct results for both MultiPoints and Points.

​https://github.com/postgis/postgis/pull/87

[dbaston]

After GEOS 3.6 is released, consider re-enabling PreparedIntersects for Points/MultiPoints if there is a performance benefit.

[dbaston]

Interim (or final?) fix applied to trunk at r14686, 2.2 at r14687

[Algunenano]

After GEOS 3.6 is released, consider re-enabling PreparedIntersects? for Points/MultiPoints? if there is a performance benefit.

I've just tested it:

Current trunk:

template_postgis=# WITH
a AS (SELECT ST_Collect(ST_MakePoint(random(), random())) as geom from generate_series(1, 1000)),
b AS (SELECT (ST_Dump(geom)).geom FROM a)
SELECT count(*) FROM b INNER JOIN a ON ST_Intersects(a.geom, b.geom);
 count 
-------
  1000
(1 row)

Time: 2366.630 ms (00:02.367)
template_postgis=# WITH
a AS (SELECT ST_Collect(ST_MakePoint(random(), random())) as geom from generate_series(1, 1000)),
b AS (SELECT (ST_Dump(geom)).geom FROM a)
SELECT count(*) FROM b INNER JOIN a ON ST_Intersects(a.geom, b.geom);
 count 
-------
  1000
(1 row)

Time: 2355.091 ms (00:02.355)
template_postgis=# WITH
a AS (SELECT ST_Collect(ST_MakePoint(random(), random())) as geom from generate_series(1, 1000)),
b AS (SELECT (ST_Dump(geom)).geom FROM a)
SELECT count(*) FROM b INNER JOIN a ON ST_Intersects(a.geom, b.geom);
 count 
-------
  1000
(1 row)

Time: 2344.796 ms (00:02.345)

Removing the patch since GEOS 3.6 is now the minimum requirement:

template_postgis=# WITH
a AS (SELECT ST_Collect(ST_MakePoint(random(), random())) as geom from generate_series(1, 1000)),
b AS (SELECT (ST_Dump(geom)).geom FROM a)
SELECT count(*) FROM b INNER JOIN a ON ST_Intersects(a.geom, b.geom);
 count 
-------
  1000
(1 row)

Time: 89.477 ms
template_postgis=# WITH
a AS (SELECT ST_Collect(ST_MakePoint(random(), random())) as geom from generate_series(1, 1000)),
b AS (SELECT (ST_Dump(geom)).geom FROM a)
SELECT count(*) FROM b INNER JOIN a ON ST_Intersects(a.geom, b.geom);
 count 
-------
  1000
(1 row)

Time: 81.209 ms
template_postgis=# WITH
a AS (SELECT ST_Collect(ST_MakePoint(random(), random())) as geom from generate_series(1, 1000)),
b AS (SELECT (ST_Dump(geom)).geom FROM a)
SELECT count(*) FROM b INNER JOIN a ON ST_Intersects(a.geom, b.geom);
 count 
-------
  1000
(1 row)

Time: 78.410 ms
template_postgis=# WITH
a AS (SELECT ST_Collect(ST_MakePoint(random(), random())) as geom from generate_series(1, 1000)),
b AS (SELECT (ST_Dump(geom)).geom FROM a)
SELECT count(*) FROM b INNER JOIN a ON ST_Intersects(a.geom, b.geom);
 count 
-------
  1000
(1 row)

Time: 82.055 ms

It's a nice 30x that we have there. I'll remove it now

[algunenano]

In 17554:

Use GEOSPrepare also for Points and Multipoints

This was removed waiting for a fix that was introduced
in GEOS 3.6, which is now the minimum required version.

Basic tests show this makes ST_Intersects with cached
points 30x faster.

Closes #3437

[Algunenano]

I've only applied it to trunk. If required it could be backported to previous stable releases depending on the GEOS release.