libgeos / geos

Geometry Engine, Open Source
https://libgeos.org
GNU Lesser General Public License v2.1
1.19k stars 358 forks source link

Within function doesn't consider POINT A is within a GEOMETRYCOLLECTION containing A. #982

Open cuteDen-ECNU opened 11 months ago

cuteDen-ECNU commented 11 months ago

Consider the following statement:

SELECT ST_Within(a1, a2), ST_Within(a1, a3)
FROM ST_GeomFromText('POINT(0 0)') As a1
, ST_GeomFromText('GEOMETRYCOLLECTION(POINT(0 0), LINESTRING(0 0, 1 0))') As a2
, ST_GeomFromText('GEOMETRYCOLLECTION(LINESTRING(0 0, 1 0), POINT(0 0))') As a3;
--expected{t, t}, actual{f, t}

According to the definition of ST_Within and interior of POINT, geometry a1 is within a2

ST_Within(A, B) ⇔ (A ⋂ B = A) ∧ (Int(A) ⋂ Int(B) ≠ ∅) For POINTs, the interior is the point itself.

If the POINT(0 0) is behind LINESTRING(0 0, 1 0) (such as a3), Postgis gives the correct answer for ST_Within(a1, a3). When the POINT(0 0) is in front of LINESTRING(0 0, 1 0) that is shown in a2, Postgis doesn't consider the 'Within' relationship.

ST_With(a1,a2) and ST_Within (a1, a2) give different results while the only difference between a2 and a3 is the order of the linestring and the point. It makes me believe it is a functional issue.

The version is the lastest: POSTGIS="3.5.0dev 3.4.0rc1-705-g5c3ec8392" [EXTENSION] PGSQL="170" GEOS="3.13.0dev-CAPI-1.18.0" PROJ="8.2.1 NETWORK_ENABLED=OFF URL_ENDPOINT=https://cdn.proj.org USER_WRITABLE_DIRECTORY=/tmp/proj DATABASE_PATH=/usr/share/proj/proj.db" LIBXML="2.9.13"

pramsey commented 11 months ago

Reproduces for me under GEOS test.

diff --git a/tests/unit/geom/GeometryCollectionTest.cpp b/tests/unit/geom/GeometryCollectionTest.cpp
index 41477bb2f..e42e7f358 100644
--- a/tests/unit/geom/GeometryCollectionTest.cpp
+++ b/tests/unit/geom/GeometryCollectionTest.cpp
@@ -181,4 +181,27 @@ void object::test<9>
     ensure(gc->hasDimension(geos::geom::Dimension::A));
 }

+
+// https://github.com/libgeos/geos/issues/982
+template<>
+template<>
+void object::test<10>
+()
+{
+    auto g1 = readWKT("POINT(0 0)");
+    auto g2 = readWKT("GEOMETRYCOLLECTION(POINT(0 0), LINESTRING(0 0, 1 0))");
+    auto g3 = readWKT("GEOMETRYCOLLECTION(LINESTRING(0 0, 1 0), POINT(0 0))");
+
+    // std::cout << "g1->within(g2.get()) = "   << g1->within(g2.get()) << std::endl;
+    // std::cout << "g1->within(g3.get()) = "   << g1->within(g3.get()) << std::endl;
+    // std::cout << "g2->contains(g1.get()) = " << g2->contains(g1.get()) << std::endl;
+    // std::cout << "g3->contains(g1.get()) = " << g3->contains(g1.get()) << std::endl;
+
+    ensure(g1->within(g2.get()));
+    ensure(g1->within(g3.get()));
+    ensure(g2->contains(g1.get()));
+    ensure(g3->contains(g1.get()));
+}
+
pramsey commented 11 months ago

The control runs all the way down into RelateOp, and the DE9IM matrix returned looks bad to me.

select st_relate('POINT(0 0)', 'GEOMETRYCOLLECTION(POINT(0 0), LINESTRING(0 0, 1 0))');

F0F
FFF
102

Says there's no interior/interior interaction. Meanwhile, just looking at points.

select st_relate('POINT(0 0)', 'POINT(0 0)');

0FF
FFF
FF2

Same in JTS @dr-jts ?

dr-jts commented 11 months ago

Same in JTS?

JTS doesn't support predicates on GeometryCollections - that is a GEOS extension. And it's tricky to do - as evidenced by this bug.

dr-jts commented 11 months ago

This is one of those tricky cases which arise when working with mixed-type GeometryCollections. It's also running into the well-known quirk of the OGC semantics for contain and within.

Under the OGC definition of within, a geometry is within another only if their Interiors have at least one point in common. In this case, POINT (0 0) lies on the Boundary of the LINESTRING(0 0, 0 1) element. And it lies in the Interior of the POINT (0 0) element of the collection. So the relate logic has to make a choice about how to evaluate the topological location of the shared point (0 0).

The OGC SF doesn't really define the semantics for the topology of GeometryCollections. A reasonable approach is to treat the GC as if it is the Union of its elements. This gives a well-founded basis for evaluating the topology. In this case, that determines that the LineString "absorbs" the Point, so the location (0 0) is in the Boundary of the GC. This gives:

Within(POINT(0 0), GEOMETRYCOLLECTION(POINT(0 0), LINESTRING(0 0, 1 0))) == FALSE

However, at the moment the GEOS logic uses a "last-one-wins" approach to computing the boundary points. This is what is causing the asymmetry of the Within result. This needs to be fixed to use "Boundary-priority" semantics. (This is relative easy to do for the Point/Line case. The case where an element falls inside the interior of a polygon is potentially a lot more computationally expensive, however).