Branch data Line data Source code
1 : : /***************************************************************************
2 : : qgszonalstatistics.h - description
3 : : ----------------------------
4 : : begin : August 29th, 2009
5 : : copyright : (C) 2009 by Marco Hugentobler
6 : : email : marco at hugis dot net
7 : : ***************************************************************************/
8 : :
9 : : /***************************************************************************
10 : : * *
11 : : * This program is free software; you can redistribute it and/or modify *
12 : : * it under the terms of the GNU General Public License as published by *
13 : : * the Free Software Foundation; either version 2 of the License, or *
14 : : * (at your option) any later version. *
15 : : * *
16 : : ***************************************************************************/
17 : :
18 : : #ifndef QGSZONALSTATISTICS_H
19 : : #define QGSZONALSTATISTICS_H
20 : :
21 : : #include <QString>
22 : : #include <QMap>
23 : :
24 : : #include <limits>
25 : : #include <cfloat>
26 : :
27 : : #include "qgis_analysis.h"
28 : : #include "qgsfeedback.h"
29 : : #include "qgscoordinatereferencesystem.h"
30 : : #include "qgsfields.h"
31 : :
32 : : class QgsGeometry;
33 : : class QgsVectorLayer;
34 : : class QgsRasterLayer;
35 : : class QgsRasterInterface;
36 : : class QgsRasterDataProvider;
37 : : class QgsRectangle;
38 : : class QgsField;
39 : : class QgsFeatureSink;
40 : : class QgsFeatureSource;
41 : :
42 : : /**
43 : : * \ingroup analysis
44 : : * \brief A class that calculates raster statistics (count, sum, mean) for a polygon or multipolygon layer and appends the results as attributes.
45 : : */
46 : 0 : class ANALYSIS_EXPORT QgsZonalStatistics
47 : : {
48 : : public:
49 : :
50 : : //! Enumeration of flags that specify statistics to be calculated
51 : : enum Statistic
52 : : {
53 : : Count = 1, //!< Pixel count
54 : : Sum = 2, //!< Sum of pixel values
55 : : Mean = 4, //!< Mean of pixel values
56 : : Median = 8, //!< Median of pixel values
57 : : StDev = 16, //!< Standard deviation of pixel values
58 : : Min = 32, //!< Min of pixel values
59 : : Max = 64, //!< Max of pixel values
60 : : Range = 128, //!< Range of pixel values (max - min)
61 : : Minority = 256, //!< Minority of pixel values
62 : : Majority = 512, //!< Majority of pixel values
63 : : Variety = 1024, //!< Variety (count of distinct) pixel values
64 : : Variance = 2048, //!< Variance of pixel values
65 : : All = Count | Sum | Mean | Median | StDev | Max | Min | Range | Minority | Majority | Variety | Variance
66 : : };
67 : : Q_DECLARE_FLAGS( Statistics, Statistic )
68 : :
69 : : //! Error codes for calculation
70 : : enum Result
71 : : {
72 : : Success = 0, //!< Success
73 : : LayerTypeWrong = 1, //!< Layer is not a polygon layer
74 : : LayerInvalid, //!< Layer is invalid
75 : : RasterInvalid, //!< Raster layer is invalid
76 : : RasterBandInvalid, //!< The raster band does not exist on the raster layer
77 : : FailedToCreateField = 8, //!< Output fields could not be created
78 : : Canceled = 9 //!< Algorithm was canceled
79 : : };
80 : :
81 : : /**
82 : : * Convenience constructor for QgsZonalStatistics, using an input raster layer.
83 : : *
84 : : * The raster layer must exist for the lifetime of the zonal statistics calculation.
85 : : *
86 : : * \warning Constructing QgsZonalStatistics using this method is not thread safe, and
87 : : * the constructor which accepts a QgsRasterInterface should be used instead.
88 : : */
89 : : QgsZonalStatistics( QgsVectorLayer *polygonLayer,
90 : : QgsRasterLayer *rasterLayer,
91 : : const QString &attributePrefix = QString(),
92 : : int rasterBand = 1,
93 : : QgsZonalStatistics::Statistics stats = QgsZonalStatistics::Statistics( QgsZonalStatistics::Count | QgsZonalStatistics::Sum | QgsZonalStatistics::Mean ) );
94 : :
95 : : /**
96 : : * Constructor for QgsZonalStatistics, using a QgsRasterInterface.
97 : : *
98 : : * The \a polygonLayer gives the vector layer containing the (multi)polygon features corresponding to the
99 : : * different zones. This layer will be modified, adding extra attributes for each of the zonal statistics
100 : : * calculated.
101 : : *
102 : : * Pixel values for each zone are taken from the raster \a rasterInterface. The constructor must also
103 : : * be given various properties relating to the input raster, such as the raster CRS (\a rasterCrs),
104 : : * and the size (X and Y) in map units for each raster pixel. The source raster band is specified
105 : : * via \a rasterBand, where a value of 1 corresponds to the first band.
106 : : *
107 : : * If the CRS of the \a polygonLayer and \a rasterCrs differ, the calculation will automatically
108 : : * reproject the zones to ensure valid results are calculated.
109 : : *
110 : : * The \a attributePrefix argument specifies an optional prefix to use when creating the
111 : : * new fields for each calculated statistic.
112 : : *
113 : : * Finally, the calculated statistics can be set via the \a stats argument. A new field will be
114 : : * added to \a polygonLayer for each statistic calculated.
115 : : *
116 : : * \warning The raster interface must exist for the lifetime of the zonal statistics calculation. For thread
117 : : * safe use, always use a cloned raster interface.
118 : : *
119 : : * \since QGIS 3.2
120 : : */
121 : : QgsZonalStatistics( QgsVectorLayer *polygonLayer,
122 : : QgsRasterInterface *rasterInterface,
123 : : const QgsCoordinateReferenceSystem &rasterCrs,
124 : : double rasterUnitsPerPixelX,
125 : : double rasterUnitsPerPixelY,
126 : : const QString &attributePrefix = QString(),
127 : : int rasterBand = 1,
128 : : QgsZonalStatistics::Statistics stats = QgsZonalStatistics::Statistics( QgsZonalStatistics::Count | QgsZonalStatistics::Sum | QgsZonalStatistics::Mean ) );
129 : :
130 : :
131 : : /**
132 : : * Runs the calculation.
133 : : */
134 : : QgsZonalStatistics::Result calculateStatistics( QgsFeedback *feedback );
135 : :
136 : : /**
137 : : * Returns the friendly display name for a \a statistic.
138 : : * \see shortName()
139 : : * \since QGIS 3.12
140 : : */
141 : : static QString displayName( QgsZonalStatistics::Statistic statistic );
142 : :
143 : : /**
144 : : * Returns a short, friendly display name for a \a statistic, suitable for use in a field name.
145 : : * \see displayName()
146 : : * \since QGIS 3.12
147 : : */
148 : : static QString shortName( QgsZonalStatistics::Statistic statistic );
149 : :
150 : : /**
151 : : * Calculates the specified \a statistics for the pixels of \a rasterBand
152 : : * in \a rasterInterface (a raster layer dataProvider() ) within polygon \a geometry.
153 : : *
154 : : * Returns a map of statistic to result value.
155 : : *
156 : : * \since QGIS 3.16
157 : : */
158 : : static QMap<QgsZonalStatistics::Statistic, QVariant> calculateStatistics( QgsRasterInterface *rasterInterface, const QgsGeometry &geometry, double cellSizeX, double cellSizeY, int rasterBand, QgsZonalStatistics::Statistics statistics );
159 : :
160 : : private:
161 : : QgsZonalStatistics() = default;
162 : :
163 : 0 : class FeatureStats
164 : : {
165 : : public:
166 : 0 : FeatureStats( bool storeValues = false, bool storeValueCounts = false )
167 : 0 : : mStoreValues( storeValues )
168 : 0 : , mStoreValueCounts( storeValueCounts )
169 : : {
170 : 0 : }
171 : :
172 : 0 : void reset()
173 : : {
174 : 0 : sum = 0;
175 : 0 : count = 0;
176 : 0 : max = std::numeric_limits<double>::lowest();
177 : 0 : min = std::numeric_limits<double>::max();
178 : 0 : valueCount.clear();
179 : 0 : values.clear();
180 : 0 : }
181 : :
182 : 0 : void addValue( double value, double weight = 1.0 )
183 : : {
184 : 0 : if ( weight < 1.0 )
185 : : {
186 : 0 : sum += value * weight;
187 : 0 : count += weight;
188 : 0 : }
189 : : else
190 : : {
191 : 0 : sum += value;
192 : 0 : ++count;
193 : : }
194 : 0 : min = std::min( min, value );
195 : 0 : max = std::max( max, value );
196 : 0 : if ( mStoreValueCounts )
197 : 0 : valueCount.insert( value, valueCount.value( value, 0 ) + 1 );
198 : 0 : if ( mStoreValues )
199 : 0 : values.append( value );
200 : 0 : }
201 : 0 : double sum = 0.0;
202 : 0 : double count = 0.0;
203 : 0 : double max = std::numeric_limits<double>::lowest();
204 : 0 : double min = std::numeric_limits<double>::max();
205 : : QMap< double, int > valueCount;
206 : : QList< double > values;
207 : :
208 : : private:
209 : : bool mStoreValues = false;
210 : : bool mStoreValueCounts = false;
211 : : };
212 : :
213 : : QString getUniqueFieldName( const QString &fieldName, const QList<QgsField> &newFields );
214 : :
215 : : QgsRasterInterface *mRasterInterface = nullptr;
216 : : QgsCoordinateReferenceSystem mRasterCrs;
217 : :
218 : : double mCellSizeX = 0;
219 : : double mCellSizeY = 0;
220 : :
221 : : //! Raster band to calculate statistics
222 : : int mRasterBand = 0;
223 : : QgsVectorLayer *mPolygonLayer = nullptr;
224 : : QString mAttributePrefix;
225 : : Statistics mStatistics = QgsZonalStatistics::All;
226 : : };
227 : :
228 : : Q_DECLARE_OPERATORS_FOR_FLAGS( QgsZonalStatistics::Statistics )
229 : :
230 : : // clazy:excludeall=qstring-allocations
231 : :
232 : : #endif // QGSZONALSTATISTICS_H
|
