ST_PointOnSurface handling of invalid.

Prints a loud notice and passes input through MakeValid.

Closes #4103
Closes https://github.com/postgis/postgis/pull/257

git-svn-id: http://svn.osgeo.org/postgis/trunk@16609 b70326c6-7e19-0410-871a-916f4a2858ee

diff --git a/NEWS b/NEWS
index eeb88139d9c69834e3436535d48ceea8f91c1d03..5b61d147249af1274c00ab890cae8633c36c096e 100644 (file)
--- a/NEWS
+++ b/NEWS
@@ -77,6 +77,7 @@ PostGIS 2.5.0
   - #4084: Fixed wrong code-comment regarding front/back of BOX3D (Matthias Bay)
   - #4060, #4094, PostgreSQL JIT support (Raúl Marín, Laurenz Albe)
   - #3960, ST_Centroid now uses lwgeom_centroid (Darafei Praliaskouski)
+  - #4103, ST_PointOnSurface can handle invalid (Darafei Praliaskouski)
 
 PostGIS 2.4.4
 2018/04/08

diff --git a/liblwgeom/liblwgeom.h.in b/liblwgeom/liblwgeom.h.in
index 14101191a969e57a1de56d17bb5d2c41f9c02f08..0291c1b2b568dca84f0406c977d15f961568b69e 100644 (file)
--- a/liblwgeom/liblwgeom.h.in
+++ b/liblwgeom/liblwgeom.h.in
@@ -2275,6 +2275,7 @@ LWGEOM *lwgeom_normalize(const LWGEOM *geom);
 LWGEOM *lwgeom_intersection(const LWGEOM *geom1, const LWGEOM *geom2);
 LWGEOM *lwgeom_difference(const LWGEOM *geom1, const LWGEOM *geom2);
 LWGEOM *lwgeom_symdifference(const LWGEOM* geom1, const LWGEOM* geom2);
+LWGEOM *lwgeom_pointonsurface(const LWGEOM* geom);
 LWGEOM *lwgeom_centroid(const LWGEOM* geom);
 LWGEOM *lwgeom_union(const LWGEOM *geom1, const LWGEOM *geom2);
 LWGEOM *lwgeom_linemerge(const LWGEOM *geom1);

diff --git a/liblwgeom/lwgeom_geos.c b/liblwgeom/lwgeom_geos.c
index 14722e3f20a30260ee77b05396784bc911e2540e..9da485fb6b8be36dea0da9c891bab6a2212b2e67 100644 (file)
--- a/liblwgeom/lwgeom_geos.c
+++ b/liblwgeom/lwgeom_geos.c
@@ -543,7 +543,7 @@ output_geos_as_lwgeom(GEOSGeometry** g, LWGEOM** geom, const int32_t srid, const
        return LW_TRUE;
 }
 
-/* Output encoder and sanity checker for GEOS wrappers */
+/* Clean up and return NULL */
 inline static LWGEOM*
 geos_clean_and_fail(GEOSGeometry* g1, GEOSGeometry* g2, GEOSGeometry* g3, const char* funcname)
 {
@@ -554,7 +554,7 @@ geos_clean_and_fail(GEOSGeometry* g1, GEOSGeometry* g2, GEOSGeometry* g3, const
        return NULL;
 }
 
-/* Output encoder and sanity checker for GEOS wrappers */
+/* Clean up */
 inline static void
 geos_clean(GEOSGeometry* g1, GEOSGeometry* g2, GEOSGeometry* g3)
 {
@@ -829,6 +829,55 @@ lwgeom_centroid(const LWGEOM* geom)
        return result;
 }
 
+LWGEOM *
+lwgeom_pointonsurface(const LWGEOM *geom)
+{
+       LWGEOM *result;
+       int32_t srid = get_result_srid(geom, NULL, __func__);
+       uint8_t is3d = FLAGS_GET_Z(geom->flags);
+       GEOSGeometry *g1, *g3;
+
+       if (srid == SRID_INVALID) return NULL;
+
+       if (lwgeom_is_empty(geom))
+       {
+               LWPOINT *lwp = lwpoint_construct_empty(srid, is3d, lwgeom_has_m(geom));
+               return lwpoint_as_lwgeom(lwp);
+       }
+
+       initGEOS(lwnotice, lwgeom_geos_error);
+
+       if (!input_lwgeom_to_geos(&g1, geom, __func__)) return NULL;
+
+       g3 = GEOSPointOnSurface(g1);
+
+       if (!g3)
+       {
+               GEOSGeometry *g1v;
+               lwnotice("%s: GEOS Error: %s", __func__, lwgeom_geos_errmsg);
+
+               if (!GEOSisValid(g1))
+               {
+                       lwnotice(
+                           "Your geometry dataset is not valid per OGC Specification. "
+                           "Please fix it with manual review of entries that are not ST_IsValid(geom). "
+                           "Retrying GEOS operation with ST_MakeValid of your input.");
+                       g1v = LWGEOM_GEOS_makeValid(g1);
+                       g3 = GEOSPointOnSurface(g1v);
+                       geos_clean(g1v, NULL, NULL);
+               }
+       }
+
+       if (!g3) return geos_clean_and_fail(g1, NULL, NULL, __func__);
+
+       if (!output_geos_as_lwgeom(&g3, &result, srid, is3d, __func__))
+               return geos_clean_and_fail(g1, NULL, g3, __func__);
+
+       geos_clean(g1, NULL, g3);
+
+       return result;
+}
+
 LWGEOM*
 lwgeom_union(const LWGEOM* geom1, const LWGEOM* geom2)
 {

diff --git a/postgis/lwgeom_geos.c b/postgis/lwgeom_geos.c
index cf915affa9ddcde9e695ca048f4271c918f238c6..7ad59262ece0649f91f4e79b56f4a8bb7e3718d4 100644 (file)
--- a/postgis/lwgeom_geos.c
+++ b/postgis/lwgeom_geos.c
@@ -1250,55 +1250,20 @@ Datum geos_difference(PG_FUNCTION_ARGS)
 PG_FUNCTION_INFO_V1(pointonsurface);
 Datum pointonsurface(PG_FUNCTION_ARGS)
 {
-       GSERIALIZED *geom;
-       GEOSGeometry *g1, *g3;
-       GSERIALIZED *result;
+       GSERIALIZED *geom, *result;
+       LWGEOM *lwgeom, *lwresult;
 
        geom = PG_GETARG_GSERIALIZED_P(0);
 
-       /* Empty.PointOnSurface == Point Empty */
-       if ( gserialized_is_empty(geom) )
-       {
-               LWPOINT *lwp = lwpoint_construct_empty(
-                                  gserialized_get_srid(geom),
-                                  gserialized_has_z(geom),
-                                  gserialized_has_m(geom));
-               result = geometry_serialize(lwpoint_as_lwgeom(lwp));
-               lwpoint_free(lwp);
-               PG_RETURN_POINTER(result);
-       }
-
-       initGEOS(lwpgnotice, lwgeom_geos_error);
-
-       g1 = POSTGIS2GEOS(geom);
-
-       if (!g1) HANDLE_GEOS_ERROR(__func__);
-
-       g3 = GEOSPointOnSurface(g1);
-
-       if (!g3)
-       {
-               GEOSGeom_destroy(g1);
-               HANDLE_GEOS_ERROR("GEOSPointOnSurface");
-       }
-
-       GEOSSetSRID(g3, gserialized_get_srid(geom));
-
-       result = GEOS2POSTGIS(g3, gserialized_has_z(geom));
-
-       if (!result)
-       {
-               GEOSGeom_destroy(g1);
-               GEOSGeom_destroy(g3);
-               elog(ERROR,"GEOS pointonsurface() threw an error (result postgis geometry formation)!");
-               PG_RETURN_NULL(); /* never get here */
-       }
-
-       GEOSGeom_destroy(g1);
-       GEOSGeom_destroy(g3);
-
+       lwgeom = lwgeom_from_gserialized(geom);
+       lwresult = lwgeom_pointonsurface(lwgeom);
+       lwgeom_free(lwgeom);
        PG_FREE_IF_COPY(geom, 0);
 
+       if (!lwresult) PG_RETURN_NULL();
+
+       result = geometry_serialize(lwresult);
+       lwgeom_free(lwresult);
        PG_RETURN_POINTER(result);
 }
 

diff --git a/regress/tickets.sql b/regress/tickets.sql
index 9560f3866e8bc7df711bcae8c655e20417d954cf..9e8a7bb98a812bbee542b69b30c7070d3915ec0b 100644 (file)
--- a/regress/tickets.sql
+++ b/regress/tickets.sql
@@ -1088,5 +1088,8 @@ SELECT '#4055b', ST_SRID(unnest(ST_ClusterIntersecting(ARRAY['SRID=4326;POINT (3
 --#4089
 select '#4089', st_astext(st_geomfromtwkb(st_AsTWKB(st_GeometryFromText('LINESTRING Z(1 1 1, 3 3 1)'), 1, 0, 0, false, true)));
 
+select '#4103', ST_Intersects(ST_PointOnSurface(geom), geom)
+from (select '0103000020110F0000010000000A000000000000C41E644741000000EEA2A75A41000000F420644741000000629EA75A410000007A2D644741000000E49FA75A41000000C02E644741000000409DA75A41000000286A64474100000064A4A75A410000007867644741000000FAA9A75A41000000E82B644741000000D2A2A75A41000000222D64474100000046A0A75A41000000242B6447410000006CA4A75A41000000C41E644741000000EEA2A75A41'::geometry geom) z;
+
 -- Clean up
 DELETE FROM spatial_ref_sys;

diff --git a/regress/tickets_expected b/regress/tickets_expected
index 26873ddfdb289d9086eebd349c35ace27d7da4d7..e727304ea34f27732bbaec84e8952d00e958d68a 100644 (file)
--- a/regress/tickets_expected
+++ b/regress/tickets_expected
@@ -349,3 +349,8 @@ NOTICE:  Self-intersection
 #4055a|4326
 #4055b|4326
 #4089|LINESTRING Z (1 1 1,3 3 1)
+NOTICE:  lwgeom_pointonsurface: GEOS Error: TopologyException: Input geom 1 is invalid: Self-intersection
+NOTICE:  Self-intersection
+NOTICE:  Your geometry dataset is not valid per OGC Specification. Please fix it with manual review of entries that are not ST_IsValid(geom). Retrying GEOS operation with ST_MakeValid of your input.
+NOTICE:  Self-intersection
+#4103|t