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
18 #include "qgsalgorithmrasterlayeruniquevalues.h"
19 #include "qgsstringutils.h"
20 #include <QTextStream>
21
22 ///@cond PRIVATE
23
name() const24 QString QgsRasterLayerUniqueValuesReportAlgorithm::name() const
25 {
26 return QStringLiteral( "rasterlayeruniquevaluesreport" );
27 }
28
displayName() const29 QString QgsRasterLayerUniqueValuesReportAlgorithm::displayName() const
30 {
31 return QObject::tr( "Raster layer unique values report" );
32 }
33
tags() const34 QStringList QgsRasterLayerUniqueValuesReportAlgorithm::tags() const
35 {
36 return QObject::tr( "count,area,statistics" ).split( ',' );
37 }
38
group() const39 QString QgsRasterLayerUniqueValuesReportAlgorithm::group() const
40 {
41 return QObject::tr( "Raster analysis" );
42 }
43
groupId() const44 QString QgsRasterLayerUniqueValuesReportAlgorithm::groupId() const
45 {
46 return QStringLiteral( "rasteranalysis" );
47 }
48
initAlgorithm(const QVariantMap &)49 void QgsRasterLayerUniqueValuesReportAlgorithm::initAlgorithm( const QVariantMap & )
50 {
51 addParameter( new QgsProcessingParameterRasterLayer( QStringLiteral( "INPUT" ),
52 QObject::tr( "Input layer" ) ) );
53 addParameter( new QgsProcessingParameterBand( QStringLiteral( "BAND" ),
54 QObject::tr( "Band number" ), 1, QStringLiteral( "INPUT" ) ) );
55 addParameter( new QgsProcessingParameterFileDestination( QStringLiteral( "OUTPUT_HTML_FILE" ),
56 QObject::tr( "Unique values report" ), QObject::tr( "HTML files (*.html)" ), QVariant(), true ) );
57 addParameter( new QgsProcessingParameterFeatureSink( QStringLiteral( "OUTPUT_TABLE" ),
58 QObject::tr( "Unique values table" ), QgsProcessing::TypeVector, QVariant(), true, false ) );
59
60 addOutput( new QgsProcessingOutputString( QStringLiteral( "EXTENT" ), QObject::tr( "Extent" ) ) );
61 addOutput( new QgsProcessingOutputString( QStringLiteral( "CRS_AUTHID" ), QObject::tr( "CRS authority identifier" ) ) );
62 addOutput( new QgsProcessingOutputNumber( QStringLiteral( "WIDTH_IN_PIXELS" ), QObject::tr( "Width in pixels" ) ) );
63 addOutput( new QgsProcessingOutputNumber( QStringLiteral( "HEIGHT_IN_PIXELS" ), QObject::tr( "Height in pixels" ) ) );
64 addOutput( new QgsProcessingOutputNumber( QStringLiteral( "TOTAL_PIXEL_COUNT" ), QObject::tr( "Total pixel count" ) ) );
65 addOutput( new QgsProcessingOutputNumber( QStringLiteral( "NODATA_PIXEL_COUNT" ), QObject::tr( "NODATA pixel count" ) ) );
66 }
67
shortHelpString() const68 QString QgsRasterLayerUniqueValuesReportAlgorithm::shortHelpString() const
69 {
70 return QObject::tr( "This algorithm returns the count and area of each unique value in a given raster layer." );
71 }
72
createInstance() const73 QgsRasterLayerUniqueValuesReportAlgorithm *QgsRasterLayerUniqueValuesReportAlgorithm::createInstance() const
74 {
75 return new QgsRasterLayerUniqueValuesReportAlgorithm();
76 }
77
prepareAlgorithm(const QVariantMap & parameters,QgsProcessingContext & context,QgsProcessingFeedback *)78 bool QgsRasterLayerUniqueValuesReportAlgorithm::prepareAlgorithm( const QVariantMap ¶meters, QgsProcessingContext &context, QgsProcessingFeedback * )
79 {
80 QgsRasterLayer *layer = parameterAsRasterLayer( parameters, QStringLiteral( "INPUT" ), context );
81 const int band = parameterAsInt( parameters, QStringLiteral( "BAND" ), context );
82
83 if ( !layer )
84 throw QgsProcessingException( invalidRasterError( parameters, QStringLiteral( "INPUT" ) ) );
85
86 mBand = parameterAsInt( parameters, QStringLiteral( "BAND" ), context );
87 if ( mBand < 1 || mBand > layer->bandCount() )
88 throw QgsProcessingException( QObject::tr( "Invalid band number for BAND (%1): Valid values for input raster are 1 to %2" ).arg( mBand )
89 .arg( layer->bandCount() ) );
90
91 mInterface.reset( layer->dataProvider()->clone() );
92 mHasNoDataValue = layer->dataProvider()->sourceHasNoDataValue( band );
93 mLayerWidth = layer->width();
94 mLayerHeight = layer->height();
95 mExtent = layer->extent();
96 mCrs = layer->crs();
97 mRasterUnitsPerPixelX = layer->rasterUnitsPerPixelX();
98 mRasterUnitsPerPixelY = layer->rasterUnitsPerPixelY();
99 mSource = layer->source();
100
101 return true;
102 }
103
processAlgorithm(const QVariantMap & parameters,QgsProcessingContext & context,QgsProcessingFeedback * feedback)104 QVariantMap QgsRasterLayerUniqueValuesReportAlgorithm::processAlgorithm( const QVariantMap ¶meters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
105 {
106 const QString outputFile = parameterAsFileOutput( parameters, QStringLiteral( "OUTPUT_HTML_FILE" ), context );
107
108 QString areaUnit = QgsUnitTypes::toAbbreviatedString( QgsUnitTypes::distanceToAreaUnit( mCrs.mapUnits() ) );
109
110 QString tableDest;
111 std::unique_ptr< QgsFeatureSink > sink;
112 if ( parameters.contains( QStringLiteral( "OUTPUT_TABLE" ) ) && parameters.value( QStringLiteral( "OUTPUT_TABLE" ) ).isValid() )
113 {
114 QgsFields outFields;
115 outFields.append( QgsField( QStringLiteral( "value" ), QVariant::Double, QString(), 20, 8 ) );
116 outFields.append( QgsField( QStringLiteral( "count" ), QVariant::Int, QString(), 20 ) );
117 outFields.append( QgsField( areaUnit.replace( QStringLiteral( "²" ), QStringLiteral( "2" ) ), QVariant::Double, QString(), 20, 8 ) );
118 sink.reset( parameterAsSink( parameters, QStringLiteral( "OUTPUT_TABLE" ), context, tableDest, outFields, QgsWkbTypes::NoGeometry, QgsCoordinateReferenceSystem() ) );
119 if ( !sink )
120 throw QgsProcessingException( invalidSinkError( parameters, QStringLiteral( "OUTPUT_TABLE" ) ) );
121 }
122
123 QHash< double, qgssize > uniqueValues;
124 qgssize noDataCount = 0;
125
126 const qgssize layerSize = static_cast< qgssize >( mLayerWidth ) * static_cast< qgssize >( mLayerHeight );
127 const int maxWidth = QgsRasterIterator::DEFAULT_MAXIMUM_TILE_WIDTH;
128 const int maxHeight = QgsRasterIterator::DEFAULT_MAXIMUM_TILE_HEIGHT;
129 const int nbBlocksWidth = std::ceil( 1.0 * mLayerWidth / maxWidth );
130 const int nbBlocksHeight = std::ceil( 1.0 * mLayerHeight / maxHeight );
131 const int nbBlocks = nbBlocksWidth * nbBlocksHeight;
132
133 QgsRasterIterator iter( mInterface.get() );
134 iter.startRasterRead( mBand, mLayerWidth, mLayerHeight, mExtent );
135
136 int iterLeft = 0;
137 int iterTop = 0;
138 int iterCols = 0;
139 int iterRows = 0;
140 bool isNoData = false;
141 std::unique_ptr< QgsRasterBlock > rasterBlock;
142 while ( iter.readNextRasterPart( mBand, iterCols, iterRows, rasterBlock, iterLeft, iterTop ) )
143 {
144 feedback->setProgress( 100 * ( ( iterTop / maxHeight * nbBlocksWidth ) + iterLeft / maxWidth ) / nbBlocks );
145 for ( int row = 0; row < iterRows; row++ )
146 {
147 if ( feedback->isCanceled() )
148 break;
149 for ( int column = 0; column < iterCols; column++ )
150 {
151 const double value = rasterBlock->valueAndNoData( row, column, isNoData );
152 if ( mHasNoDataValue && isNoData )
153 {
154 noDataCount++;
155 }
156 else
157 {
158 uniqueValues[ value ]++;
159 }
160 }
161 }
162 }
163
164 QMap< double, qgssize > sortedUniqueValues;
165 for ( auto it = uniqueValues.constBegin(); it != uniqueValues.constEnd(); ++it )
166 {
167 sortedUniqueValues.insert( it.key(), it.value() );
168 }
169
170 QVariantMap outputs;
171 outputs.insert( QStringLiteral( "EXTENT" ), mExtent.toString() );
172 outputs.insert( QStringLiteral( "CRS_AUTHID" ), mCrs.authid() );
173 outputs.insert( QStringLiteral( "WIDTH_IN_PIXELS" ), mLayerWidth );
174 outputs.insert( QStringLiteral( "HEIGHT_IN_PIXELS" ), mLayerHeight );
175 outputs.insert( QStringLiteral( "TOTAL_PIXEL_COUNT" ), layerSize );
176 outputs.insert( QStringLiteral( "NODATA_PIXEL_COUNT" ), noDataCount );
177
178 const double pixelArea = mRasterUnitsPerPixelX * mRasterUnitsPerPixelY;
179
180 if ( !outputFile.isEmpty() )
181 {
182 QFile file( outputFile );
183 if ( file.open( QIODevice::WriteOnly | QIODevice::Truncate ) )
184 {
185 const QString encodedAreaUnit = QgsStringUtils::ampersandEncode( areaUnit );
186
187 QTextStream out( &file );
188 #if QT_VERSION < QT_VERSION_CHECK(6, 0, 0)
189 out.setCodec( "UTF-8" );
190 #endif
191 out << QStringLiteral( "<html><head><meta http-equiv=\"Content-Type\" content=\"text/html;charset=utf-8\"/></head><body>\n" );
192 out << QStringLiteral( "<p>%1: %2 (%3 %4)</p>\n" ).arg( QObject::tr( "Analyzed file" ), mSource, QObject::tr( "band" ) ).arg( mBand );
193 out << QObject::tr( "<p>%1: %2</p>\n" ).arg( QObject::tr( "Extent" ), mExtent.toString() );
194 out << QObject::tr( "<p>%1: %2</p>\n" ).arg( QObject::tr( "Projection" ), mCrs.userFriendlyIdentifier() );
195 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 );
196 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 );
197 out << QObject::tr( "<p>%1: %2</p>\n" ).arg( QObject::tr( "Total pixel count" ) ).arg( layerSize );
198 if ( mHasNoDataValue )
199 out << QObject::tr( "<p>%1: %2</p>\n" ).arg( QObject::tr( "NODATA pixel count" ) ).arg( noDataCount );
200 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 );
201
202 for ( auto it = sortedUniqueValues.constBegin(); it != sortedUniqueValues.constEnd(); ++it )
203 {
204 const double area = it.value() * pixelArea;
205 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 ) );
206 }
207 out << QStringLiteral( "</table>\n</body></html>" );
208 outputs.insert( QStringLiteral( "OUTPUT_HTML_FILE" ), outputFile );
209 }
210 }
211
212 if ( sink )
213 {
214 for ( auto it = sortedUniqueValues.constBegin(); it != sortedUniqueValues.constEnd(); ++it )
215 {
216 QgsFeature f;
217 const double area = it.value() * pixelArea;
218 f.setAttributes( QgsAttributes() << it.key() << it.value() << area );
219 if ( !sink->addFeature( f, QgsFeatureSink::FastInsert ) )
220 throw QgsProcessingException( writeFeatureError( sink.get(), parameters, QStringLiteral( "OUTPUT_TABLE" ) ) );
221 }
222 outputs.insert( QStringLiteral( "OUTPUT_TABLE" ), tableDest );
223 }
224
225 return outputs;
226 }
227
228
229 ///@endcond
230
231
232
233