[FEATURE] Expose GEOS polygonize operation via QgsGeometry

python/core/geometry/qgsgeometry.sip

@@ -777,6 +777,8 @@ class QgsGeometry
      */
     static QgsGeometry unaryUnion( const QList<QgsGeometry>& geometryList );
 
+    static QgsGeometry polygonize( const QList< QgsGeometry>& lines );
+
     /** Converts the geometry to straight line segments, if it is a curved geometry type.
      * @note added in QGIS 2.10
      * @see requiresConversionToStraightSegments

src/core/geometry/qgsgeometry.cpp

@@ -1940,7 +1940,25 @@ bool QgsGeometry::isGeosEqual( const QgsGeometry& g ) const
   return geos.isEqual( *( g.d->geometry ) );
 }
 
-QgsGeometry QgsGeometry::unaryUnion( const QList<QgsGeometry>& geometryList )
+QgsGeometry QgsGeometry::unaryUnion( const QList<QgsGeometry>& geometries )
+{
+  QgsGeos geos( nullptr );
+
+  QList<QgsAbstractGeometry*> geomV2List;
+  QList<QgsGeometry>::const_iterator it = geometries.constBegin();
+  for ( ; it != geometries.constEnd(); ++it )
+  {
+    if ( !(( *it ).isNull() ) )
+    {
+      geomV2List.append(( *it ).geometry() );
+    }
+  }
+
+  QgsAbstractGeometry* geom = geos.combine( geomV2List );
+  return QgsGeometry( geom );
+}
+
+QgsGeometry QgsGeometry::polygonize( const QList<QgsGeometry>& geometryList )
 {
   QgsGeos geos( nullptr );
 
@@ -1954,8 +1972,7 @@ QgsGeometry QgsGeometry::unaryUnion( const QList<QgsGeometry>& geometryList )
     }
   }
 
-  QgsAbstractGeometry* geom = geos.combine( geomV2List );
-  return QgsGeometry( geom );
+  return geos.polygonize( geomV2List );
 }
 
 void QgsGeometry::convertToStraightSegment()

src/core/geometry/qgsgeometry.h

@@ -843,11 +843,21 @@ class CORE_EXPORT QgsGeometry
      **/
     void validateGeometry( QList<Error> &errors );
 
-    /** Compute the unary union on a list of geometries. May be faster than an iterative union on a set of geometries.
-     * @param geometryList a list of QgsGeometry as input
-     * @returns the new computed QgsGeometry, or an empty QgsGeometry
+    /** Compute the unary union on a list of \a geometries. May be faster than an iterative union on a set of geometries.
+     * The returned geometry will be fully noded, i.e. a node will be created at every common intersection of the
+     * input geometries. An empty geometry will be returned in the case of errors.
      */
-    static QgsGeometry unaryUnion( const QList<QgsGeometry>& geometryList );
+    static QgsGeometry unaryUnion( const QList<QgsGeometry>& geometries );
+
+    /**
+     * Creates a GeometryCollection geometry containing possible polygons formed from the constituent
+     * linework of a set of \a geometries. The input geometries must be fully noded (i.e. nodes exist
+     * at every common intersection of the geometries). The easiest way to ensure this is to first
+     * call unaryUnion() on the set of input geometries and then pass the result to polygonize().
+     * An empty geometry will be returned in the case of errors.
+     * @note added in QGIS 3.0
+     */
+    static QgsGeometry polygonize( const QList< QgsGeometry>& geometries );
 
     /** Converts the geometry to straight line segments, if it is a curved geometry type.
      * @note added in QGIS 2.10

src/core/geometry/qgsgeos.cpp

@@ -1941,6 +1941,45 @@ double QgsGeos::lineLocatePoint( const QgsPointV2& point, QString* errorMsg ) co
   return distance;
 }
 
+QgsGeometry QgsGeos::polygonize( const QList<QgsAbstractGeometry*>& geometries, QString* errorMsg )
+{
+  GEOSGeometry** const lineGeosGeometries = new GEOSGeometry*[ geometries.size()];
+  int validLines = 0;
+  Q_FOREACH ( const QgsAbstractGeometry* g, geometries )
+  {
+    GEOSGeometry* l = asGeos( g );
+    if ( l )
+    {
+      lineGeosGeometries[validLines] = l;
+      validLines++;
+    }
+  }
+
+  try
+  {
+    GEOSGeomScopedPtr result( GEOSPolygonize_r( geosinit.ctxt, lineGeosGeometries, validLines ) );
+    for ( int i = 0; i < validLines; ++i )
+    {
+      GEOSGeom_destroy_r( geosinit.ctxt, lineGeosGeometries[i] );
+    }
+    delete[] lineGeosGeometries;
+    return QgsGeometry( fromGeos( result.get() ) );
+  }
+  catch ( GEOSException &e )
+  {
+    if ( errorMsg )
+    {
+      *errorMsg = e.what();
+    }
+    for ( int i = 0; i < validLines; ++i )
+    {
+      GEOSGeom_destroy_r( geosinit.ctxt, lineGeosGeometries[i] );
+    }
+    delete[] lineGeosGeometries;
+    return QgsGeometry();
+  }
+}
+
 
 //! Extract coordinates of linestring's endpoints. Returns false on error.
 static bool _linestringEndpoints( const GEOSGeometry* linestring, double& x1, double& y1, double& x2, double& y2 )

src/core/geometry/qgsgeos.h

@@ -141,6 +141,8 @@ class CORE_EXPORT QgsGeos: public QgsGeometryEngine
      */
     double lineLocatePoint( const QgsPointV2& point, QString* errorMsg = nullptr ) const;
 
+    static QgsGeometry polygonize( const QList<QgsAbstractGeometry*>& geometries, QString* errorMsg = nullptr );
+
     /** Create a geometry from a GEOSGeometry
      * @param geos GEOSGeometry. Ownership is NOT transferred.
      */

tests/src/python/test_qgsgeometry.py

@@ -3724,6 +3724,37 @@ class TestQgsGeometry(unittest.TestCase):
         self.assertTrue(compareWkt(result, exp, 0.00001),
                         "orthogonalize: mismatch Expected:\n{}\nGot:\n{}\n".format(exp, result))
 
+    def testPolygonize(self):
+        o = QgsGeometry.polygonize([])
+        self.assertFalse(o)
+        empty = QgsGeometry()
+        o = QgsGeometry.polygonize([empty])
+        self.assertFalse(o)
+        line = QgsGeometry.fromWkt('LineString()')
+        o = QgsGeometry.polygonize([line])
+        self.assertFalse(o)
+
+        l1 = QgsGeometry.fromWkt("LINESTRING (100 180, 20 20, 160 20, 100 180)")
+        l2 = QgsGeometry.fromWkt("LINESTRING (100 180, 80 60, 120 60, 100 180)")
+        o = QgsGeometry.polygonize([l1, l2])
+        exp = "GeometryCollection(POLYGON ((100 180, 160 20, 20 20, 100 180), (100 180, 80 60, 120 60, 100 180)),POLYGON ((100 180, 120 60, 80 60, 100 180)))"
+        result = o.exportToWkt()
+        self.assertTrue(compareWkt(result, exp, 0.00001),
+                        "polygonize: mismatch Expected:\n{}\nGot:\n{}\n".format(exp, result))
+
+        lines = [QgsGeometry.fromWkt('LineString(0 0, 1 1)'),
+                 QgsGeometry.fromWkt('LineString(0 0, 0 1)'),
+                 QgsGeometry.fromWkt('LineString(0 1, 1 1)'),
+                 QgsGeometry.fromWkt('LineString(1 1, 1 0)'),
+                 QgsGeometry.fromWkt('LineString(1 0, 0 0)'),
+                 QgsGeometry.fromWkt('LineString(5 5, 6 6)'),
+                 QgsGeometry.fromWkt('Point(0, 0)')]
+        o = QgsGeometry.polygonize(lines)
+        exp = "GeometryCollection (Polygon ((0 0, 1 1, 1 0, 0 0)),Polygon ((1 1, 0 0, 0 1, 1 1)))"
+        result = o.exportToWkt()
+        self.assertTrue(compareWkt(result, exp, 0.00001),
+                        "polygonize: mismatch Expected:\n{}\nGot:\n{}\n".format(exp, result))
+
 
 if __name__ == '__main__':
     unittest.main()