[FEATURE] New API QgsGeometry::densifyByCount

Densifies a geometry by adding a specified number of vertices to each segment
qgis · Mar 25, 2017 · 5360b79 · 5360b79
1 parent 2dac2d3
commit 5360b79
Show file tree

Hide file tree

Showing 6 changed files with 261 additions and 11 deletions.
diff --git a/python/core/geometry/qgsgeometry.sip b/python/core/geometry/qgsgeometry.sip
@@ -541,17 +541,8 @@ class QgsGeometry
      * @note added in QGIS 3.0
      */
     QgsGeometry extendLine( double startDistance, double endDistance ) const;
-
-    /** Returns a simplified version of this geometry using a specified tolerance value */
-    QgsGeometry simplify( double tolerance )  const;
-
-    /**
-     * Returns the center of mass of a geometry.
-     * @note for line based geometries, the center point of the line is returned,
-     * and for point based geometries, the point itself is returned
-     * @see pointOnSurface()
-     * @see poleOfInaccessibility()
-     */
+    QgsGeometry simplify( double tolerance ) const;
+    QgsGeometry densifyByCount( int extraNodesPerSegment ) const;
     QgsGeometry centroid() const;
 
     /**

diff --git a/src/core/geometry/qgsgeometry.cpp b/src/core/geometry/qgsgeometry.cpp
@@ -1576,6 +1576,13 @@ QgsGeometry QgsGeometry::simplify( double tolerance ) const
   return QgsGeometry( simplifiedGeom );
 }
 
+QgsGeometry QgsGeometry::densifyByCount( int extraNodesPerSegment ) const
+{
+  QgsInternalGeometryEngine engine( *this );
+
+  return engine.densifyByCount( extraNodesPerSegment );
+}
+
 QgsGeometry QgsGeometry::centroid() const
 {
   if ( !d->geometry )

diff --git a/src/core/geometry/qgsgeometry.h b/src/core/geometry/qgsgeometry.h
@@ -611,6 +611,16 @@ class CORE_EXPORT QgsGeometry
     //! Returns a simplified version of this geometry using a specified tolerance value
     QgsGeometry simplify( double tolerance ) const;
 
+    /**
+     * Returns a copy of the geometry which has been densified by adding the specified
+     * number of extra nodes within each segment of the geometry.
+     * If the geometry has z or m values present then these will be linearly interpolated
+     * at the added nodes.
+     * Curved geometry types are automatically segmentized by this routine.
+     * @node added in QGIS 3.0
+     */
+    QgsGeometry densifyByCount( int extraNodesPerSegment ) const;
+
     /**
      * Returns the center of mass of a geometry.
      * @note for line based geometries, the center point of the line is returned,

diff --git a/src/core/geometry/qgsinternalgeometryengine.cpp b/src/core/geometry/qgsinternalgeometryengine.cpp
@@ -476,3 +476,134 @@ QgsGeometry QgsInternalGeometryEngine::orthogonalize( double tolerance, int maxI
     return QgsGeometry( orthogonalizeGeom( mGeometry, maxIterations, tolerance, lowerThreshold, upperThreshold ) );
   }
 }
+
+QgsLineString *doDensifyByCount( QgsLineString *ring, int extraNodesPerSegment )
+{
+  QgsPointSequence out;
+  double multiplier = 1.0 / double( extraNodesPerSegment + 1 );
+
+  int nPoints = ring->numPoints();
+  out.reserve( ( extraNodesPerSegment + 1 ) * nPoints );
+  bool withZ = ring->is3D();
+  bool withM = ring->isMeasure();
+  QgsWkbTypes::Type outType = QgsWkbTypes::Point;
+  if ( ring->is3D() )
+    outType = QgsWkbTypes::addZ( outType );
+  if ( ring->isMeasure() )
+    outType = QgsWkbTypes::addM( outType );
+  double x1 = 0;
+  double x2 = 0;
+  double y1 = 0;
+  double y2 = 0;
+  double z1 = 0;
+  double z2 = 0;
+  double m1 = 0;
+  double m2 = 0;
+  double xOut = 0;
+  double yOut = 0;
+  double zOut = 0;
+  double mOut = 0;
+  for ( int i = 0; i < nPoints - 1; ++i )
+  {
+    x1 = ring->xAt( i );
+    x2 = ring->xAt( i + 1 );
+    y1 = ring->yAt( i );
+    y2 = ring->yAt( i + 1 );
+    if ( withZ )
+    {
+      z1 = ring->zAt( i );
+      z2 = ring->zAt( i + 1 );
+    }
+    if ( withM )
+    {
+      m1 = ring->mAt( i );
+      m2 = ring->mAt( i + 1 );
+    }
+
+    out << QgsPointV2( outType, x1, y1, z1, m1 );
+    for ( int j = 0; j < extraNodesPerSegment; ++j )
+    {
+      double delta = multiplier * ( j + 1 );
+      xOut = x1 + delta * ( x2 - x1 );
+      yOut = y1 + delta * ( y2 - y1 );
+      if ( withZ )
+        zOut = z1 + delta * ( z2 - z1 );
+      if ( withM )
+        mOut = m1 + delta * ( m2 - m1 );
+
+      out << QgsPointV2( outType, xOut, yOut, zOut, mOut );
+    }
+  }
+  out << QgsPointV2( outType, ring->xAt( nPoints - 1 ), ring->yAt( nPoints - 1 ),
+                     withZ ? ring->zAt( nPoints - 1 ) : 0, withM ? ring->mAt( nPoints - 1 ) : 0 );
+
+  QgsLineString *result = new QgsLineString();
+  result->setPoints( out );
+  return result;
+}
+
+QgsAbstractGeometry *densifyGeometryByCount( const QgsAbstractGeometry *geom, int extraNodesPerSegment )
+{
+  std::unique_ptr< QgsAbstractGeometry > segmentizedCopy;
+  if ( QgsWkbTypes::isCurvedType( geom->wkbType() ) )
+  {
+    segmentizedCopy.reset( geom->segmentize() );
+    geom = segmentizedCopy.get();
+  }
+
+  if ( QgsWkbTypes::geometryType( geom->wkbType() ) == QgsWkbTypes::LineGeometry )
+  {
+    return doDensifyByCount( static_cast< QgsLineString * >( geom->clone() ), extraNodesPerSegment );
+  }
+  else
+  {
+    // polygon
+    const QgsPolygonV2 *polygon = static_cast< const QgsPolygonV2 * >( geom );
+    QgsPolygonV2 *result = new QgsPolygonV2();
+
+    result->setExteriorRing( doDensifyByCount( static_cast< QgsLineString * >( polygon->exteriorRing()->clone() ),
+                             extraNodesPerSegment ) );
+    for ( int i = 0; i < polygon->numInteriorRings(); ++i )
+    {
+      result->addInteriorRing( doDensifyByCount( static_cast< QgsLineString * >( polygon->interiorRing( i )->clone() ),
+                               extraNodesPerSegment ) );
+    }
+
+    return result;
+  }
+}
+
+QgsGeometry QgsInternalGeometryEngine::densifyByCount( int extraNodesPerSegment ) const
+{
+  if ( !mGeometry )
+  {
+    return QgsGeometry();
+  }
+
+  if ( QgsWkbTypes::geometryType( mGeometry->wkbType() ) == QgsWkbTypes::PointGeometry )
+  {
+    return QgsGeometry( mGeometry->clone() ); // point geometry, nothing to do
+  }
+
+  if ( const QgsGeometryCollection *gc = dynamic_cast< const QgsGeometryCollection *>( mGeometry ) )
+  {
+    int numGeom = gc->numGeometries();
+    QList< QgsAbstractGeometry * > geometryList;
+    geometryList.reserve( numGeom );
+    for ( int i = 0; i < numGeom; ++i )
+    {
+      geometryList << densifyGeometryByCount( gc->geometryN( i ), extraNodesPerSegment );
+    }
+
+    QgsGeometry first = QgsGeometry( geometryList.takeAt( 0 ) );
+    Q_FOREACH ( QgsAbstractGeometry *g, geometryList )
+    {
+      first.addPart( g );
+    }
+    return first;
+  }
+  else
+  {
+    return QgsGeometry( densifyGeometryByCount( mGeometry, extraNodesPerSegment ) );
+  }
+}
diff --git a/src/core/geometry/qgsinternalgeometryengine.h b/src/core/geometry/qgsinternalgeometryengine.h
@@ -71,6 +71,16 @@ class QgsInternalGeometryEngine
      */
     QgsGeometry orthogonalize( double tolerance = 1.0E-8, int maxIterations = 1000, double angleThreshold = 15.0 ) const;
 
+    /**
+     * Densifies the geometry by adding the specified number of extra nodes within each
+     * segment of the geometry.
+     * If the geometry has z or m values present then these will be linearly interpolated
+     * at the added nodes.
+     * Curved geometry types are automatically segmentized by this routine.
+     * @node added in QGIS 3.0
+     */
+    QgsGeometry densifyByCount( int extraNodesPerSegment ) const;
+
   private:
     const QgsAbstractGeometry *mGeometry = nullptr;
 };

diff --git a/tests/src/python/test_qgsgeometry.py b/tests/src/python/test_qgsgeometry.py
@@ -3894,6 +3894,107 @@ def testVoronoi(self):
         self.assertTrue(compareWkt(result, exp, 0.00001),
                         "delaunay: mismatch Expected:\n{}\nGot:\n{}\n".format(exp, result))
 
+    def testDensifyByCount(self):
+
+        empty = QgsGeometry()
+        o = empty.densifyByCount(4)
+        self.assertFalse(o)
+
+        # point
+        input = QgsGeometry.fromWkt("PointZ( 1 2 3 )")
+        o = input.densifyByCount(100)
+        exp = "PointZ( 1 2 3 )"
+        result = o.exportToWkt()
+        self.assertTrue(compareWkt(result, exp, 0.00001),
+                        "densify by count: mismatch Expected:\n{}\nGot:\n{}\n".format(exp, result))
+        input = QgsGeometry.fromWkt(
+            "MULTIPOINT ((155 271), (150 360), (260 360), (271 265), (280 260), (270 370), (154 354), (150 260))")
+        o = input.densifyByCount(100)
+        exp = "MULTIPOINT ((155 271), (150 360), (260 360), (271 265), (280 260), (270 370), (154 354), (150 260))"
+        result = o.exportToWkt()
+        self.assertTrue(compareWkt(result, exp, 0.00001),
+                        "densify by count: mismatch Expected:\n{}\nGot:\n{}\n".format(exp, result))
+
+        # line
+        input = QgsGeometry.fromWkt("LineString( 0 0, 10 0, 10 10 )")
+        o = input.densifyByCount(0)
+        exp = "LineString( 0 0, 10 0, 10 10 )"
+        result = o.exportToWkt()
+        self.assertTrue(compareWkt(result, exp, 0.00001),
+                        "densify by count: mismatch Expected:\n{}\nGot:\n{}\n".format(exp, result))
+        o = input.densifyByCount(1)
+        exp = "LineString( 0 0, 5 0, 10 0, 10 5, 10 10 )"
+        result = o.exportToWkt()
+        self.assertTrue(compareWkt(result, exp, 0.00001),
+                        "densify by count: mismatch Expected:\n{}\nGot:\n{}\n".format(exp, result))
+        o = input.densifyByCount(3)
+        exp = "LineString( 0 0, 2.5 0, 5 0, 7.5 0, 10 0, 10 2.5, 10 5, 10 7.5, 10 10 )"
+        result = o.exportToWkt()
+        self.assertTrue(compareWkt(result, exp, 0.00001),
+                        "densify by count: mismatch Expected:\n{}\nGot:\n{}\n".format(exp, result))
+        input = QgsGeometry.fromWkt("LineStringZ( 0 0 1, 10 0 2, 10 10 0)")
+        o = input.densifyByCount(1)
+        exp = "LineStringZ( 0 0 1, 5 0 1.5, 10 0 2, 10 5 1, 10 10 0 )"
+        result = o.exportToWkt()
+        self.assertTrue(compareWkt(result, exp, 0.00001),
+                        "densify by count: mismatch Expected:\n{}\nGot:\n{}\n".format(exp, result))
+        input = QgsGeometry.fromWkt("LineStringM( 0 0 0, 10 0 2, 10 10 0)")
+        o = input.densifyByCount(1)
+        exp = "LineStringM( 0 0 0, 5 0 1, 10 0 2, 10 5 1, 10 10 0 )"
+        result = o.exportToWkt()
+        self.assertTrue(compareWkt(result, exp, 0.00001),
+                        "densify by count: mismatch Expected:\n{}\nGot:\n{}\n".format(exp, result))
+        input = QgsGeometry.fromWkt("LineStringZM( 0 0 1 10, 10 0 2 8, 10 10 0 4)")
+        o = input.densifyByCount(1)
+        exp = "LineStringZM( 0 0 1 10, 5 0 1.5 9, 10 0 2 8, 10 5 1 6, 10 10 0 4 )"
+        result = o.exportToWkt()
+        self.assertTrue(compareWkt(result, exp, 0.00001),
+                        "densify by count: mismatch Expected:\n{}\nGot:\n{}\n".format(exp, result))
+
+        # polygon
+        input = QgsGeometry.fromWkt("Polygon(( 0 0, 10 0, 10 10, 0 0 ))")
+        o = input.densifyByCount(0)
+        exp = "Polygon(( 0 0, 10 0, 10 10, 0 0 ))"
+        result = o.exportToWkt()
+        self.assertTrue(compareWkt(result, exp, 0.00001),
+                        "densify by count: mismatch Expected:\n{}\nGot:\n{}\n".format(exp, result))
+
+        input = QgsGeometry.fromWkt("PolygonZ(( 0 0 1, 10 0 2, 10 10 0, 0 0 1 ))")
+        o = input.densifyByCount(1)
+        exp = "PolygonZ(( 0 0 1, 5 0 1.5, 10 0 2, 10 5 1, 10 10 0, 5 5 0.5, 0 0 1 ))"
+        result = o.exportToWkt()
+        self.assertTrue(compareWkt(result, exp, 0.00001),
+                        "densify by count: mismatch Expected:\n{}\nGot:\n{}\n".format(exp, result))
+
+        input = QgsGeometry.fromWkt("PolygonZM(( 0 0 1 4, 10 0 2 6, 10 10 0 8, 0 0 1 4 ))")
+        o = input.densifyByCount(1)
+        exp = "PolygonZM(( 0 0 1 4, 5 0 1.5 5, 10 0 2 6, 10 5 1 7, 10 10 0 8, 5 5 0.5 6, 0 0 1 4 ))"
+        result = o.exportToWkt()
+        self.assertTrue(compareWkt(result, exp, 0.00001),
+                        "densify by count: mismatch Expected:\n{}\nGot:\n{}\n".format(exp, result))
+        # (not strictly valid, but shouldn't matter!
+        input = QgsGeometry.fromWkt("PolygonZM(( 0 0 1 4, 10 0 2 6, 10 10 0 8, 0 0 1 4 ), ( 0 0 1 4, 10 0 2 6, 10 10 0 8, 0 0 1 4 ) )")
+        o = input.densifyByCount(1)
+        exp = "PolygonZM(( 0 0 1 4, 5 0 1.5 5, 10 0 2 6, 10 5 1 7, 10 10 0 8, 5 5 0.5 6, 0 0 1 4 ),( 0 0 1 4, 5 0 1.5 5, 10 0 2 6, 10 5 1 7, 10 10 0 8, 5 5 0.5 6, 0 0 1 4 ))"
+        result = o.exportToWkt()
+        self.assertTrue(compareWkt(result, exp, 0.00001),
+                        "densify by count: mismatch Expected:\n{}\nGot:\n{}\n".format(exp, result))
+
+        # multi line
+        input = QgsGeometry.fromWkt("MultiLineString(( 0 0, 5 0, 10 0, 10 5, 10 10), (20 0, 25 0, 30 0, 30 5, 30 10 ) )")
+        o = input.densifyByCount(1)
+        exp = "MultiLineString(( 0 0, 2.5 0, 5 0, 7.5 0, 10 0, 10 2.5, 10 5, 10 7.5, 10 10 ),( 20 0, 22.5 0, 25 0, 27.5 0, 30 0, 30 2.5, 30 5, 30 7.5, 30 10 ))"
+        result = o.exportToWkt()
+        self.assertTrue(compareWkt(result, exp, 0.00001),
+                        "densify by count: mismatch Expected:\n{}\nGot:\n{}\n".format(exp, result))
+
+        # multipolygon
+        input = QgsGeometry.fromWkt("MultiPolygonZ(((  0 0 1, 10 0 2, 10 10 0, 0 0 1)),((  0 0 1, 10 0 2, 10 10 0, 0 0 1 )))")
+        o = input.densifyByCount(1)
+        exp = "MultiPolygonZ((( 0 0 1, 5 0 1.5, 10 0 2, 10 5 1, 10 10 0, 5 5 0.5, 0 0 1 )),(( 0 0 1, 5 0 1.5, 10 0 2, 10 5 1, 10 10 0, 5 5 0.5, 0 0 1 )))"
+        result = o.exportToWkt()
+        self.assertTrue(compareWkt(result, exp, 0.00001),
+                        "densify by count: mismatch Expected:\n{}\nGot:\n{}\n".format(exp, result))
 
 if __name__ == '__main__':
     unittest.main()