[zonal statistics] add variance statistic

qgis · Jun 29, 2017 · 58f6f93 · 58f6f93
1 parent 8379f9b
commit 58f6f93
Show file tree

Hide file tree

Showing 4 changed files with 29 additions and 5 deletions.
diff --git a/python/analysis/vector/qgszonalstatistics.sip b/python/analysis/vector/qgszonalstatistics.sip
@@ -36,6 +36,7 @@ class QgsZonalStatistics
       Minority,
       Majority,
       Variety,
+      Variance,
       All
     };
     typedef QFlags<QgsZonalStatistics::Statistic> Statistics;

diff --git a/src/analysis/vector/qgszonalstatistics.cpp b/src/analysis/vector/qgszonalstatistics.cpp
@@ -158,6 +158,13 @@ int QgsZonalStatistics::calculateStatistics( QgsFeedback *feedback )
     QgsField varietyField( varietyFieldName, QVariant::Int, QStringLiteral( "int" ) );
     newFieldList.push_back( varietyField );
   }
+  QString varianceFieldName;
+  if ( mStatistics & QgsZonalStatistics::Variance )
+  {
+    varianceFieldName = getUniqueFieldName( mAttributePrefix + "variance", newFieldList );
+    QgsField varianceField( varianceFieldName, QVariant::Double, QStringLiteral( "double precision" ) );
+    newFieldList.push_back( varianceField );
+  }
   vectorProvider->addAttributes( newFieldList );
 
   //index of the new fields
@@ -172,6 +179,7 @@ int QgsZonalStatistics::calculateStatistics( QgsFeedback *feedback )
   int minorityIndex = mStatistics & QgsZonalStatistics::Minority ? vectorProvider->fieldNameIndex( minorityFieldName ) : -1;
   int majorityIndex = mStatistics & QgsZonalStatistics::Majority ? vectorProvider->fieldNameIndex( majorityFieldName ) : -1;
   int varietyIndex = mStatistics & QgsZonalStatistics::Variety ? vectorProvider->fieldNameIndex( varietyFieldName ) : -1;
+  int varianceIndex = mStatistics & QgsZonalStatistics::Variance ? vectorProvider->fieldNameIndex( varianceFieldName ) : -1;
 
   if ( ( mStatistics & QgsZonalStatistics::Count && countIndex == -1 )
        || ( mStatistics & QgsZonalStatistics::Sum && sumIndex == -1 )
@@ -184,6 +192,7 @@ int QgsZonalStatistics::calculateStatistics( QgsFeedback *feedback )
        || ( mStatistics & QgsZonalStatistics::Minority && minorityIndex == -1 )
        || ( mStatistics & QgsZonalStatistics::Majority && majorityIndex == -1 )
        || ( mStatistics & QgsZonalStatistics::Variety && varietyIndex == -1 )
+       || ( mStatistics & QgsZonalStatistics::Variance && varianceIndex == -1 )
      )
   {
     //failed to create a required field
@@ -200,7 +209,8 @@ int QgsZonalStatistics::calculateStatistics( QgsFeedback *feedback )
   QgsFeature f;
 
   bool statsStoreValues = ( mStatistics & QgsZonalStatistics::Median ) ||
-                          ( mStatistics & QgsZonalStatistics::StDev );
+                          ( mStatistics & QgsZonalStatistics::StDev ) ||
+                          ( mStatistics & QgsZonalStatistics::Variance );
   bool statsStoreValueCount = ( mStatistics & QgsZonalStatistics::Minority ) ||
                               ( mStatistics & QgsZonalStatistics::Majority );
 
@@ -288,16 +298,22 @@ int QgsZonalStatistics::calculateStatistics( QgsFeedback *feedback )
         }
         changeAttributeMap.insert( medianIndex, QVariant( medianValue ) );
       }
-      if ( mStatistics & QgsZonalStatistics::StDev )
+      if ( mStatistics & QgsZonalStatistics::StDev || mStatistics & QgsZonalStatistics::Variance )
       {
         double sumSquared = 0;
         for ( int i = 0; i < featureStats.values.count(); ++i )
         {
           double diff = featureStats.values.at( i ) - mean;
           sumSquared += diff * diff;
         }
-        double stdev = qPow( sumSquared / featureStats.values.count(), 0.5 );
-        changeAttributeMap.insert( stdevIndex, QVariant( stdev ) );
+        double variance = sumSquared / featureStats.values.count();
+        if ( mStatistics & QgsZonalStatistics::StDev )
+        {
+          double stdev = qPow( variance, 0.5 );
+          changeAttributeMap.insert( stdevIndex, QVariant( stdev ) );
+        }
+        if ( mStatistics & QgsZonalStatistics::Variance )
+          changeAttributeMap.insert( varianceIndex, QVariant( variance ) );
       }
       if ( mStatistics & QgsZonalStatistics::Min )
         changeAttributeMap.insert( minIndex, QVariant( featureStats.min ) );

diff --git a/src/analysis/vector/qgszonalstatistics.h b/src/analysis/vector/qgszonalstatistics.h
@@ -55,7 +55,8 @@ class ANALYSIS_EXPORT QgsZonalStatistics
       Minority = 256, //!< Minority of pixel values
       Majority = 512, //!< Majority of pixel values
       Variety = 1024, //!< Variety (count of distinct) pixel values
-      All = Count | Sum | Mean | Median | StDev | Max | Min | Range | Minority | Majority | Variety
+      Variance = 2048, //!< Variance of pixel values
+      All = Count | Sum | Mean | Median | StDev | Max | Min | Range | Minority | Majority | Variety | Variance
     };
     Q_DECLARE_FLAGS( Statistics, Statistic )
 

diff --git a/tests/src/analysis/testqgszonalstatistics.cpp b/tests/src/analysis/testqgszonalstatistics.cpp
@@ -94,6 +94,7 @@ void TestQgsZonalStatistics::testStatistics()
   QCOMPARE( f.attribute( "minority" ).toDouble(), 0.0 );
   QCOMPARE( f.attribute( "majority" ).toDouble(), 1.0 );
   QCOMPARE( f.attribute( "variety" ).toDouble(), 2.0 );
+  QCOMPARE( f.attribute( "variance" ).toDouble(), 0.222222222222222 );
 
   request.setFilterFid( 1 );
   fetched = mVectorLayer->getFeatures( request ).nextFeature( f );
@@ -109,6 +110,7 @@ void TestQgsZonalStatistics::testStatistics()
   QCOMPARE( f.attribute( "minority" ).toDouble(), 0.0 );
   QCOMPARE( f.attribute( "majority" ).toDouble(), 1.0 );
   QCOMPARE( f.attribute( "variety" ).toDouble(), 2.0 );
+  QCOMPARE( f.attribute( "variance" ).toDouble(), 0.24691358024691 );
 
   request.setFilterFid( 2 );
   fetched = mVectorLayer->getFeatures( request ).nextFeature( f );
@@ -124,6 +126,7 @@ void TestQgsZonalStatistics::testStatistics()
   QCOMPARE( f.attribute( "minority" ).toDouble(), 0.0 );
   QCOMPARE( f.attribute( "majority" ).toDouble(), 1.0 );
   QCOMPARE( f.attribute( "variety" ).toDouble(), 2.0 );
+  QCOMPARE( f.attribute( "variance" ).toDouble(), 0.13888888888889 );
 
   // same with long prefix to ensure that field name truncation handled correctly
   QgsZonalStatistics zsl( mVectorLayer, mRasterLayer, QStringLiteral( "myqgis2_" ), 1, QgsZonalStatistics::All );
@@ -143,6 +146,7 @@ void TestQgsZonalStatistics::testStatistics()
   QCOMPARE( f.attribute( "myqgis2__2" ).toDouble(), 0.0 );
   QCOMPARE( f.attribute( "myqgis2__3" ).toDouble(), 1.0 );
   QCOMPARE( f.attribute( "myqgis2_va" ).toDouble(), 2.0 );
+  QCOMPARE( f.attribute( "myqgis2__4" ).toDouble(), 0.222222222222222 );
 
   request.setFilterFid( 1 );
   fetched = mVectorLayer->getFeatures( request ).nextFeature( f );
@@ -158,6 +162,7 @@ void TestQgsZonalStatistics::testStatistics()
   QCOMPARE( f.attribute( "myqgis2__2" ).toDouble(), 0.0 );
   QCOMPARE( f.attribute( "myqgis2__3" ).toDouble(), 1.0 );
   QCOMPARE( f.attribute( "myqgis2_va" ).toDouble(), 2.0 );
+  QCOMPARE( f.attribute( "myqgis2__4" ).toDouble(), 0.24691358024691 );
 
   request.setFilterFid( 2 );
   fetched = mVectorLayer->getFeatures( request ).nextFeature( f );
@@ -173,6 +178,7 @@ void TestQgsZonalStatistics::testStatistics()
   QCOMPARE( f.attribute( "myqgis2__2" ).toDouble(), 0.0 );
   QCOMPARE( f.attribute( "myqgis2__3" ).toDouble(), 1.0 );
   QCOMPARE( f.attribute( "myqgis2_va" ).toDouble(), 2.0 );
+  QCOMPARE( f.attribute( "myqgis2__4" ).toDouble(), 0.13888888888889 );
 }
 
 QGSTEST_MAIN( TestQgsZonalStatistics )