QGIS API Documentation 3.40.0-Bratislava (b56115d8743)
Loading...
Searching...
No Matches
qgsalgorithmrasterlayeruniquevalues.cpp
Go to the documentation of this file.
1/***************************************************************************
2 qgsalgorithmrasterlayeruniquevalues.cpp
3 ---------------------
4 begin : April 2017
5 copyright : (C) 2017 by Mathieu Pellerin
6 email : nirvn dot asia at gmail dot com
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
19#include "qgsstringutils.h"
20#include "qgsunittypes.h"
21
22#include <QTextStream>
23
25
26QString QgsRasterLayerUniqueValuesReportAlgorithm::name() const
27{
28 return QStringLiteral( "rasterlayeruniquevaluesreport" );
29}
30
31QString QgsRasterLayerUniqueValuesReportAlgorithm::displayName() const
32{
33 return QObject::tr( "Raster layer unique values report" );
34}
35
36QStringList QgsRasterLayerUniqueValuesReportAlgorithm::tags() const
37{
38 return QObject::tr( "count,area,statistics" ).split( ',' );
39}
40
41QString QgsRasterLayerUniqueValuesReportAlgorithm::group() const
42{
43 return QObject::tr( "Raster analysis" );
44}
45
46QString QgsRasterLayerUniqueValuesReportAlgorithm::groupId() const
47{
48 return QStringLiteral( "rasteranalysis" );
49}
50
51void QgsRasterLayerUniqueValuesReportAlgorithm::initAlgorithm( const QVariantMap & )
52{
53 addParameter( new QgsProcessingParameterRasterLayer( QStringLiteral( "INPUT" ),
54 QObject::tr( "Input layer" ) ) );
55 addParameter( new QgsProcessingParameterBand( QStringLiteral( "BAND" ),
56 QObject::tr( "Band number" ), 1, QStringLiteral( "INPUT" ) ) );
57 addParameter( new QgsProcessingParameterFileDestination( QStringLiteral( "OUTPUT_HTML_FILE" ),
58 QObject::tr( "Unique values report" ), QObject::tr( "HTML files (*.html)" ), QVariant(), true ) );
59 addParameter( new QgsProcessingParameterFeatureSink( QStringLiteral( "OUTPUT_TABLE" ),
60 QObject::tr( "Unique values table" ), Qgis::ProcessingSourceType::Vector, QVariant(), true, false ) );
61
62 addOutput( new QgsProcessingOutputString( QStringLiteral( "EXTENT" ), QObject::tr( "Extent" ) ) );
63 addOutput( new QgsProcessingOutputString( QStringLiteral( "CRS_AUTHID" ), QObject::tr( "CRS authority identifier" ) ) );
64 addOutput( new QgsProcessingOutputNumber( QStringLiteral( "WIDTH_IN_PIXELS" ), QObject::tr( "Width in pixels" ) ) );
65 addOutput( new QgsProcessingOutputNumber( QStringLiteral( "HEIGHT_IN_PIXELS" ), QObject::tr( "Height in pixels" ) ) );
66 addOutput( new QgsProcessingOutputNumber( QStringLiteral( "TOTAL_PIXEL_COUNT" ), QObject::tr( "Total pixel count" ) ) );
67 addOutput( new QgsProcessingOutputNumber( QStringLiteral( "NODATA_PIXEL_COUNT" ), QObject::tr( "NoData pixel count" ) ) );
68}
69
70QString QgsRasterLayerUniqueValuesReportAlgorithm::shortHelpString() const
71{
72 return QObject::tr( "This algorithm returns the count and area of each unique value in a given raster layer. "
73 "The area calculation is done in the area unit of the layer's CRS." );
74}
75
76QgsRasterLayerUniqueValuesReportAlgorithm *QgsRasterLayerUniqueValuesReportAlgorithm::createInstance() const
77{
78 return new QgsRasterLayerUniqueValuesReportAlgorithm();
79}
80
81bool QgsRasterLayerUniqueValuesReportAlgorithm::prepareAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback * )
82{
83 QgsRasterLayer *layer = parameterAsRasterLayer( parameters, QStringLiteral( "INPUT" ), context );
84 const int band = parameterAsInt( parameters, QStringLiteral( "BAND" ), context );
85
86 if ( !layer )
87 throw QgsProcessingException( invalidRasterError( parameters, QStringLiteral( "INPUT" ) ) );
88
89 mBand = parameterAsInt( parameters, QStringLiteral( "BAND" ), context );
90 if ( mBand < 1 || mBand > layer->bandCount() )
91 throw QgsProcessingException( QObject::tr( "Invalid band number for BAND (%1): Valid values for input raster are 1 to %2" ).arg( mBand )
92 .arg( layer->bandCount() ) );
93
94 mInterface.reset( layer->dataProvider()->clone() );
95 mHasNoDataValue = layer->dataProvider()->sourceHasNoDataValue( band );
96 mLayerWidth = layer->width();
97 mLayerHeight = layer->height();
98 mExtent = layer->extent();
99 mCrs = layer->crs();
100 mRasterUnitsPerPixelX = layer->rasterUnitsPerPixelX();
101 mRasterUnitsPerPixelY = layer->rasterUnitsPerPixelY();
102 mSource = layer->source();
103
104 return true;
105}
106
107QVariantMap QgsRasterLayerUniqueValuesReportAlgorithm::processAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
108{
109 const QString outputFile = parameterAsFileOutput( parameters, QStringLiteral( "OUTPUT_HTML_FILE" ), context );
110
111 QString areaUnit = QgsUnitTypes::toAbbreviatedString( QgsUnitTypes::distanceToAreaUnit( mCrs.mapUnits() ) );
112
113 QString tableDest;
114 std::unique_ptr< QgsFeatureSink > sink;
115 if ( parameters.contains( QStringLiteral( "OUTPUT_TABLE" ) ) && parameters.value( QStringLiteral( "OUTPUT_TABLE" ) ).isValid() )
116 {
117 QgsFields outFields;
118 outFields.append( QgsField( QStringLiteral( "value" ), QMetaType::Type::Double, QString(), 20, 8 ) );
119 outFields.append( QgsField( QStringLiteral( "count" ), QMetaType::Type::LongLong, QString(), 20 ) );
120 outFields.append( QgsField( areaUnit.replace( QStringLiteral( "²" ), QStringLiteral( "2" ) ), QMetaType::Type::Double, QString(), 20, 8 ) );
121 sink.reset( parameterAsSink( parameters, QStringLiteral( "OUTPUT_TABLE" ), context, tableDest, outFields, Qgis::WkbType::NoGeometry, QgsCoordinateReferenceSystem() ) );
122 if ( !sink )
123 throw QgsProcessingException( invalidSinkError( parameters, QStringLiteral( "OUTPUT_TABLE" ) ) );
124 }
125
126 QHash< double, qgssize > uniqueValues;
127 qgssize noDataCount = 0;
128
129 const qgssize layerSize = static_cast< qgssize >( mLayerWidth ) * static_cast< qgssize >( mLayerHeight );
132 const int nbBlocksWidth = std::ceil( 1.0 * mLayerWidth / maxWidth );
133 const int nbBlocksHeight = std::ceil( 1.0 * mLayerHeight / maxHeight );
134 const int nbBlocks = nbBlocksWidth * nbBlocksHeight;
135
136 QgsRasterIterator iter( mInterface.get() );
137 iter.startRasterRead( mBand, mLayerWidth, mLayerHeight, mExtent );
138
139 int iterLeft = 0;
140 int iterTop = 0;
141 int iterCols = 0;
142 int iterRows = 0;
143 bool isNoData = false;
144 std::unique_ptr< QgsRasterBlock > rasterBlock;
145 while ( iter.readNextRasterPart( mBand, iterCols, iterRows, rasterBlock, iterLeft, iterTop ) )
146 {
147 feedback->setProgress( 100 * ( ( iterTop / maxHeight * nbBlocksWidth ) + iterLeft / maxWidth ) / nbBlocks );
148 for ( int row = 0; row < iterRows; row++ )
149 {
150 if ( feedback->isCanceled() )
151 break;
152 for ( int column = 0; column < iterCols; column++ )
153 {
154 const double value = rasterBlock->valueAndNoData( row, column, isNoData );
155 if ( mHasNoDataValue && isNoData )
156 {
157 noDataCount++;
158 }
159 else
160 {
161 uniqueValues[ value ]++;
162 }
163 }
164 }
165 if ( feedback->isCanceled() )
166 break;
167 }
168
169 QMap< double, qgssize > sortedUniqueValues;
170 for ( auto it = uniqueValues.constBegin(); it != uniqueValues.constEnd(); ++it )
171 {
172 sortedUniqueValues.insert( it.key(), it.value() );
173 }
174
175 QVariantMap outputs;
176 outputs.insert( QStringLiteral( "EXTENT" ), mExtent.toString() );
177 outputs.insert( QStringLiteral( "CRS_AUTHID" ), mCrs.authid() );
178 outputs.insert( QStringLiteral( "WIDTH_IN_PIXELS" ), mLayerWidth );
179 outputs.insert( QStringLiteral( "HEIGHT_IN_PIXELS" ), mLayerHeight );
180 outputs.insert( QStringLiteral( "TOTAL_PIXEL_COUNT" ), layerSize );
181 outputs.insert( QStringLiteral( "NODATA_PIXEL_COUNT" ), noDataCount );
182
183 const double pixelArea = mRasterUnitsPerPixelX * mRasterUnitsPerPixelY;
184
185 if ( !outputFile.isEmpty() )
186 {
187 QFile file( outputFile );
188 if ( file.open( QIODevice::WriteOnly | QIODevice::Truncate ) )
189 {
190 const QString encodedAreaUnit = QgsStringUtils::ampersandEncode( areaUnit );
191
192 QTextStream out( &file );
193#if QT_VERSION < QT_VERSION_CHECK(6, 0, 0)
194 out.setCodec( "UTF-8" );
195#endif
196 out << QStringLiteral( "<html><head><meta http-equiv=\"Content-Type\" content=\"text/html;charset=utf-8\"/></head><body>\n" );
197 out << QStringLiteral( "<p>%1: %2 (%3 %4)</p>\n" ).arg( QObject::tr( "Analyzed file" ), mSource, QObject::tr( "band" ) ).arg( mBand );
198 out << QObject::tr( "<p>%1: %2</p>\n" ).arg( QObject::tr( "Extent" ), mExtent.toString() );
199 out << QObject::tr( "<p>%1: %2</p>\n" ).arg( QObject::tr( "Projection" ), mCrs.userFriendlyIdentifier() );
200 out << QObject::tr( "<p>%1: %2 (%3 %4)</p>\n" ).arg( QObject::tr( "Width in pixels" ) ).arg( mLayerWidth ).arg( QObject::tr( "units per pixel" ) ).arg( mRasterUnitsPerPixelX );
201 out << QObject::tr( "<p>%1: %2 (%3 %4)</p>\n" ).arg( QObject::tr( "Height in pixels" ) ).arg( mLayerHeight ).arg( QObject::tr( "units per pixel" ) ).arg( mRasterUnitsPerPixelY );
202 out << QObject::tr( "<p>%1: %2</p>\n" ).arg( QObject::tr( "Total pixel count" ) ).arg( layerSize );
203 if ( mHasNoDataValue )
204 out << QObject::tr( "<p>%1: %2</p>\n" ).arg( QObject::tr( "NoData pixel count" ) ).arg( noDataCount );
205 out << QStringLiteral( "<table><tr><td>%1</td><td>%2</td><td>%3 (%4)</td></tr>\n" ).arg( QObject::tr( "Value" ), QObject::tr( "Pixel count" ), QObject::tr( "Area" ), encodedAreaUnit );
206
207 for ( auto it = sortedUniqueValues.constBegin(); it != sortedUniqueValues.constEnd(); ++it )
208 {
209 const double area = it.value() * pixelArea;
210 out << QStringLiteral( "<tr><td>%1</td><td>%2</td><td>%3</td></tr>\n" ).arg( it.key() ).arg( it.value() ).arg( QString::number( area, 'g', 16 ) );
211 }
212 out << QStringLiteral( "</table>\n</body></html>" );
213 outputs.insert( QStringLiteral( "OUTPUT_HTML_FILE" ), outputFile );
214 }
215 }
216
217 if ( sink )
218 {
219 for ( auto it = sortedUniqueValues.constBegin(); it != sortedUniqueValues.constEnd(); ++it )
220 {
221 QgsFeature f;
222 const double area = it.value() * pixelArea;
223 f.setAttributes( QgsAttributes() << it.key() << it.value() << area );
224 if ( !sink->addFeature( f, QgsFeatureSink::FastInsert ) )
225 throw QgsProcessingException( writeFeatureError( sink.get(), parameters, QStringLiteral( "OUTPUT_TABLE" ) ) );
226 }
227 outputs.insert( QStringLiteral( "OUTPUT_TABLE" ), tableDest );
228 }
229
230 return outputs;
231}
232
233
@ Vector
Tables (i.e. vector layers with or without geometry). When used for a sink this indicates the sink ha...
@ NoGeometry
No geometry.
A vector of attributes.
This class represents a coordinate reference system (CRS).
@ FastInsert
Use faster inserts, at the cost of updating the passed features to reflect changes made at the provid...
The feature class encapsulates a single feature including its unique ID, geometry and a list of field...
Definition qgsfeature.h:58
void setAttributes(const QgsAttributes &attrs)
Sets the feature's attributes.
bool isCanceled() const
Tells whether the operation has been canceled already.
Definition qgsfeedback.h:53
void setProgress(double progress)
Sets the current progress for the feedback object.
Definition qgsfeedback.h:61
Encapsulate a field in an attribute table or data source.
Definition qgsfield.h:53
Container of fields for a vector layer.
Definition qgsfields.h:46
bool append(const QgsField &field, Qgis::FieldOrigin origin=Qgis::FieldOrigin::Provider, int originIndex=-1)
Appends a field.
Definition qgsfields.cpp:69
virtual QgsRectangle extent() const
Returns the extent of the layer.
QString source() const
Returns the source for the layer.
QgsCoordinateReferenceSystem crs
Definition qgsmaplayer.h:83
Contains information about the context in which a processing algorithm is executed.
Custom exception class for processing related exceptions.
Base class for providing feedback from a processing algorithm.
A numeric output for processing algorithms.
A string output for processing algorithms.
A raster band parameter for Processing algorithms.
A feature sink output for processing algorithms.
A generic file based destination parameter, for specifying the destination path for a file (non-map l...
A raster layer parameter for processing algorithms.
QgsRasterDataProvider * clone() const override=0
Clone itself, create deep copy.
virtual bool sourceHasNoDataValue(int bandNo) const
Returns true if source band has no data value.
Iterator for sequentially processing raster cells.
static const int DEFAULT_MAXIMUM_TILE_WIDTH
Default maximum tile width.
static const int DEFAULT_MAXIMUM_TILE_HEIGHT
Default maximum tile height.
Represents a raster layer.
int height() const
Returns the height of the (unclipped) raster.
int bandCount() const
Returns the number of bands in this layer.
double rasterUnitsPerPixelX() const
Returns the number of raster units per each raster pixel in X axis.
QgsRasterDataProvider * dataProvider() override
Returns the source data provider.
double rasterUnitsPerPixelY() const
Returns the number of raster units per each raster pixel in Y axis.
int width() const
Returns the width of the (unclipped) raster.
static QString ampersandEncode(const QString &string)
Makes a raw string safe for inclusion as a HTML/XML string literal.
static Q_INVOKABLE QString toAbbreviatedString(Qgis::DistanceUnit unit)
Returns a translated abbreviation representing a distance unit.
static Q_INVOKABLE Qgis::AreaUnit distanceToAreaUnit(Qgis::DistanceUnit distanceUnit)
Converts a distance unit to its corresponding area unit, e.g., meters to square meters.
unsigned long long qgssize
Qgssize is used instead of size_t, because size_t is stdlib type, unknown by SIP, and it would be har...
Definition qgis.h:6465