Add OGRGeometry::SetPrecision() / OGR_G_SetPrecision(), as wrapper of…

… GEOSGeom_setPrecision_r()

autotest/ogr/ogr_geos.py

@@ -624,3 +624,13 @@ def test_ogr_geos_prepared_geom():
     # Test workaround for https://github.com/libgeos/geos/pull/423
     assert not pg.Intersects(ogr.CreateGeometryFromWkt("POINT EMPTY"))
     assert not pg.Contains(ogr.CreateGeometryFromWkt("POINT EMPTY"))
+
+
+###############################################################################
+
+
+def test_ogr_geos_set_precision():
+
+    g = ogr.CreateGeometryFromWkt("LINESTRING (1 1,9 9)")
+    g = g.SetPrecision(10)
+    assert g.ExportToWkt() == "LINESTRING (0 0,10 10)"

ogr/ogr_api.h

@@ -295,6 +295,17 @@ OGRGeometryH CPL_DLL OGR_G_Normalize(OGRGeometryH) CPL_WARN_UNUSED_RESULT;
 int CPL_DLL OGR_G_IsSimple(OGRGeometryH);
 int CPL_DLL OGR_G_IsRing(OGRGeometryH);
 
+/** This option causes OGR_G_SetPrecision()
+  * to not attempt at preserving the topology */
+#define OGR_GEOS_PREC_NO_TOPO (1 << 0)
+
+/** This option causes OGR_G_SetPrecision()
+  * to retain collapsed elements */
+#define OGR_GEOS_PREC_KEEP_COLLAPSED (1 << 1)
+
+OGRGeometryH CPL_DLL OGR_G_SetPrecision(OGRGeometryH, double dfGridSize,
+                                        int nFlags) CPL_WARN_UNUSED_RESULT;
+
 OGRGeometryH CPL_DLL OGR_G_Polygonize(OGRGeometryH) CPL_WARN_UNUSED_RESULT;
 
 /*! @cond Doxygen_Suppress */

ogr/ogr_geometry.h

@@ -594,6 +594,8 @@ class CPL_DLL OGRGeometry
 
     virtual double Distance3D(const OGRGeometry *poOtherGeom) const;
 
+    OGRGeometry *SetPrecision(double dfGridSize, int nFlags) const;
+
     //! @cond Doxygen_Suppress
     // backward compatibility to non-standard method names.
     OGRBoolean Intersect(OGRGeometry *) const

ogr/ogrgeometry.cpp

@@ -6144,6 +6144,88 @@ OGRGeometryH OGR_G_SimplifyPreserveTopology(OGRGeometryH hThis,
         OGRGeometry::FromHandle(hThis)->SimplifyPreserveTopology(dTolerance));
 }
 
+/************************************************************************/
+/*                           SetPrecision()                             */
+/************************************************************************/
+
+/** Set the geometry's precision, rounding all its coordinates to the precision
+ * grid.
+ *
+ * This function is the same as the C function OGR_G_SetPrecision().
+ *
+ * This function is built on the GEOS library, check it for the definition
+ * of the geometry operation.
+ * If OGR is built without the GEOS library, this function will always fail,
+ * issuing a CPLE_NotSupported error.
+ *
+ * @param dfGridSize size of the precision grid, or 0 for FLOATING
+ *                 precision.
+ * @param nFlags The bitwise OR of zero, one or several of OGR_GEOS_PREC_NO_TOPO
+ *               and OGR_GEOS_PREC_KEEP_COLLAPSED
+ *
+ * @return a new geometry or NULL if an error occurs.
+ *
+ * @since GDAL 3.9
+ */
+
+OGRGeometry *OGRGeometry::SetPrecision(UNUSED_IF_NO_GEOS double dfGridSize,
+                                       UNUSED_IF_NO_GEOS int nFlags) const
+{
+#ifndef HAVE_GEOS
+    CPLError(CE_Failure, CPLE_NotSupported, "GEOS support not enabled.");
+    return nullptr;
+
+#else
+    OGRGeometry *poOGRProduct = nullptr;
+
+    GEOSContextHandle_t hGEOSCtxt = createGEOSContext();
+    GEOSGeom hThisGeosGeom = exportToGEOS(hGEOSCtxt);
+    if (hThisGeosGeom != nullptr)
+    {
+        GEOSGeom hGeosProduct = GEOSGeom_setPrecision_r(
+            hGEOSCtxt, hThisGeosGeom, dfGridSize, nFlags);
+        GEOSGeom_destroy_r(hGEOSCtxt, hThisGeosGeom);
+        poOGRProduct =
+            BuildGeometryFromGEOS(hGEOSCtxt, hGeosProduct, this, nullptr);
+    }
+    freeGEOSContext(hGEOSCtxt);
+    return poOGRProduct;
+
+#endif  // HAVE_GEOS
+}
+
+/************************************************************************/
+/*                         OGR_G_SetPrecision()                         */
+/************************************************************************/
+
+/** Set the geometry's precision, rounding all its coordinates to the precision
+ * grid.
+ *
+ * This function is the same as the C++ method OGRGeometry::Simplify().
+ *
+ * This function is built on the GEOS library, check it for the definition
+ * of the geometry operation.
+ * If OGR is built without the GEOS library, this function will always fail,
+ * issuing a CPLE_NotSupported error.
+ *
+ * @param hThis the geometry.
+ * @param dfGridSize size of the precision grid, or 0 for FLOATING
+ *                 precision.
+ * @param nFlags The bitwise OR of zero, one or several of OGR_GEOS_PREC_NO_TOPO
+ *               and OGR_GEOS_PREC_KEEP_COLLAPSED
+ *
+ * @return a new geometry or NULL if an error occurs.
+ *
+ * @since GDAL 3.9
+ */
+OGRGeometryH OGR_G_SetPrecision(OGRGeometryH hThis, double dfGridSize,
+                                int nFlags)
+{
+    VALIDATE_POINTER1(hThis, "OGR_G_SetPrecision", nullptr);
+    return OGRGeometry::ToHandle(
+        OGRGeometry::FromHandle(hThis)->SetPrecision(dfGridSize, nFlags));
+}
+
 /************************************************************************/
 /*                         DelaunayTriangulation()                      */
 /************************************************************************/

swig/include/ogr.i

@@ -532,6 +532,9 @@ typedef void retGetPoints;
 
 %constant char *OLMD_FID64             = "OLMD_FID64";
 
+%constant int GEOS_PREC_NO_TOPO = 1;
+%constant int GEOS_PREC_KEEP_COLLAPSED = 2;
+
 #else
 typedef int OGRErr;
 
@@ -583,6 +586,9 @@ typedef int OGRErr;
 
 #define OLMD_FID64             "OLMD_FID64"
 
+#define GEOS_PREC_NO_TOPO          1
+#define GEOS_PREC_KEEP_COLLAPSED   2
+
 #endif
 
 #if defined(SWIGCSHARP) || defined(SWIGJAVA) || defined(SWIGPYTHON)
@@ -3654,6 +3660,11 @@ public:
     return (OGRGeometryShadow*) OGR_G_MakeValidEx(self, options);
   }
 
+  %newobject SetPrecision;
+  OGRGeometryShadow* SetPrecision(double gridSize, int flags = 0) {
+    return (OGRGeometryShadow*) OGR_G_SetPrecision(self, gridSize, flags);
+  }
+
   %newobject Normalize;
   OGRGeometryShadow* Normalize() {
     return (OGRGeometryShadow*) OGR_G_Normalize(self);