Improve numerical stability of QgsAbstractGeometryV2::centroid

qgis · May 28, 2016 · 857d544e5564b29477cce26638bfb4bcc29bd230 · 857d544
1 parent 149ab30
commit 857d544e5564b29477cce26638bfb4bcc29bd230
Showing with 6 additions and 1 deletion.

+6 −1 src/core/geometry/qgsabstractgeometryv2.cpp
diff --git a/src/core/geometry/qgsabstractgeometryv2.cpp b/src/core/geometry/qgsabstractgeometryv2.cpp
@@ -174,6 +174,7 @@ QgsPointV2 QgsAbstractGeometryV2::centroid() const
   double A = 0.;
   double Cx = 0.;
   double Cy = 0.;
+  QgsPointV2 v0 = vertexAt( QgsVertexId( 0, 0, 0 ) );
   int i = 0, j = 1;
   if ( vertexAt( QgsVertexId( 0, 0, 0 ) ) != vertexAt( QgsVertexId( 0, 0, n - 1 ) ) )
   {
@@ -184,6 +185,10 @@ QgsPointV2 QgsAbstractGeometryV2::centroid() const
   {
     QgsPointV2 vi = vertexAt( QgsVertexId( 0, 0, i ) );
     QgsPointV2 vj = vertexAt( QgsVertexId( 0, 0, j ) );
+    vi.rx() -= v0.x();
+    vi.ry() -= v0.y();
+    vj.rx() -= v0.x();
+    vj.ry() -= v0.y();
     double d = vi.x() * vj.y() - vj.x() * vi.y();
     A += d;
     Cx += ( vi.x() + vj.x() ) * d;
@@ -203,7 +208,7 @@ QgsPointV2 QgsAbstractGeometryV2::centroid() const
   }
   else
   {
-    return QgsPointV2( Cx / ( 3. * A ), Cy / ( 3. * A ) );
+    return QgsPointV2( v0.x() + Cx / ( 3. * A ), v0.y() + Cy / ( 3. * A ) );
   }
 }