QGIS API Documentation 3.38.0-Grenoble (exported)
Loading...
Searching...
No Matches
qgsalgorithmrastercalculator.cpp
Go to the documentation of this file.
1/***************************************************************************
2 qgsalgorithmrastercalculator.cpp
3 ---------------------
4 begin : July 2023
5 copyright : (C) 2023 by Alexander Bruy
6 email : alexander dot bruy 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 "qgsalgorithmrastercalculator.h"
19#include "qgsrasterfilewriter.h"
20#include "qgsrastercalculator.h"
21
23
24Qgis::ProcessingAlgorithmFlags QgsRasterCalculatorAlgorithm::flags() const
25{
26 return QgsProcessingAlgorithm::flags() | Qgis::ProcessingAlgorithmFlag::HideFromModeler;
27}
28
29QString QgsRasterCalculatorAlgorithm::name() const
30{
31 return QStringLiteral( "rastercalc" );
32}
33
34QString QgsRasterCalculatorAlgorithm::displayName() const
35{
36 return QObject::tr( "Raster calculator" );
37}
38
39QStringList QgsRasterCalculatorAlgorithm::tags() const
40{
41 return QObject::tr( "raster,calculator" ).split( ',' );
42}
43
44QString QgsRasterCalculatorAlgorithm::group() const
45{
46 return QObject::tr( "Raster analysis" );
47}
48
49QString QgsRasterCalculatorAlgorithm::groupId() const
50{
51 return QStringLiteral( "rasteranalysis" );
52}
53
54QString QgsRasterCalculatorAlgorithm::shortHelpString() const
55{
56 return QObject::tr( "Performing algebraic operations using raster layers." );
57}
58
59QgsRasterCalculatorAlgorithm *QgsRasterCalculatorAlgorithm::createInstance() const
60{
61 return new QgsRasterCalculatorAlgorithm();
62}
63
64void QgsRasterCalculatorAlgorithm::initAlgorithm( const QVariantMap & )
65{
66
67 addParameter( new QgsProcessingParameterMultipleLayers( QStringLiteral( "LAYERS" ), QObject::tr( "Input layers" ), Qgis::ProcessingSourceType::Raster ) );
68 addParameter( new QgsProcessingParameterExpression( QStringLiteral( "EXPRESSION" ), QObject::tr( "Expression" ), QVariant(), QStringLiteral( "LAYERS" ), false, Qgis::ExpressionType::RasterCalculator ) );
69 std::unique_ptr<QgsProcessingParameterExtent> extentParam = std::make_unique<QgsProcessingParameterExtent>( QStringLiteral( "EXTENT" ), QObject::tr( "Output extent" ), QVariant(), true );
70 extentParam->setHelp( QObject::tr( "Extent of the output layer. If not specified, the extent will be the overall extent of all input layers" ) );
71 addParameter( extentParam.release() );
72 std::unique_ptr<QgsProcessingParameterNumber> cellSizeParam = std::make_unique<QgsProcessingParameterNumber>( QStringLiteral( "CELL_SIZE" ), QObject::tr( "Output cell size (leave empty to set automatically)" ), Qgis::ProcessingNumberParameterType::Double, QVariant(), true, 0.0 );
73 cellSizeParam->setHelp( QObject::tr( "Cell size of the output layer. If not specified, the smallest cell size from the input layers will be used" ) );
74 addParameter( cellSizeParam.release() );
75 std::unique_ptr<QgsProcessingParameterCrs> crsParam = std::make_unique<QgsProcessingParameterCrs>( QStringLiteral( "CRS" ), QObject::tr( "Output CRS" ), QVariant(), true );
76 crsParam->setHelp( QObject::tr( "CRS of the output layer. If not specified, the CRS of the first input layer will be used" ) );
77 addParameter( crsParam.release() );
78 addParameter( new QgsProcessingParameterRasterDestination( QStringLiteral( "OUTPUT" ), QObject::tr( "Calculated" ) ) );
79}
80
81bool QgsRasterCalculatorAlgorithm::prepareAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
82{
83 const QList< QgsMapLayer * > layers = parameterAsLayerList( parameters, QStringLiteral( "LAYERS" ), context );
84
85 for ( const QgsMapLayer *layer : std::as_const( layers ) )
86 {
87
88 QgsMapLayer *clonedLayer { layer->clone() };
89 clonedLayer->moveToThread( nullptr );
90 mLayers << clonedLayer;
91 }
92
93 if ( mLayers.isEmpty() )
94 {
95 feedback->reportError( QObject::tr( "No layers selected" ), false );
96 return false;
97 }
98
99 return true;
100}
101
102
103QVariantMap QgsRasterCalculatorAlgorithm::processAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
104{
105 for ( QgsMapLayer *layer : std::as_const( mLayers ) )
106 {
107 layer->moveToThread( QThread::currentThread() );
108 }
109
110 QgsCoordinateReferenceSystem crs;
111 if ( parameters.value( QStringLiteral( "CRS" ) ).isValid() )
112 {
113 crs = parameterAsCrs( parameters, QStringLiteral( "CRS" ), context );
114 }
115 else
116 {
117 crs = mLayers.at( 0 )->crs();
118 }
119
120 QgsRectangle bbox;
121 if ( parameters.value( QStringLiteral( "EXTENT" ) ).isValid() )
122 {
123 bbox = parameterAsExtent( parameters, QStringLiteral( "EXTENT" ), context, crs );
124 }
125 else
126 {
127 bbox = QgsProcessingUtils::combineLayerExtents( mLayers, crs, context );
128 }
129
130 double minCellSize = 1e9;
131
132 QVector< QgsRasterCalculatorEntry > entries;
133 for ( QgsMapLayer *layer : mLayers )
134 {
135 QgsRasterLayer *rLayer = qobject_cast<QgsRasterLayer *>( layer );
136 if ( !rLayer )
137 {
138 continue;
139 }
140
141 const int nBands = rLayer->dataProvider()->bandCount();
142 for ( int i = 0; i < nBands; ++i )
143 {
144 QgsRasterCalculatorEntry entry;
145 entry.ref = QStringLiteral( "%1@%2" ).arg( rLayer->name() ).arg( i + 1 );
146 entry.raster = rLayer;
147 entry.bandNumber = i + 1;
148 entries << entry;
149 }
150
151 QgsRectangle ext = rLayer->extent();
152 if ( rLayer->crs() != crs )
153 {
154 QgsCoordinateTransform ct( rLayer->crs(), crs, context.transformContext() );
155 ext = ct.transformBoundingBox( ext );
156 }
157
158 double cellSize = ( ext.xMaximum() - ext.xMinimum() ) / rLayer->width();
159 if ( cellSize < minCellSize )
160 {
161 minCellSize = cellSize;
162 }
163 }
164
165 double cellSize = parameterAsDouble( parameters, QStringLiteral( "CELL_SIZE" ), context );
166 if ( cellSize == 0 )
167 {
168 cellSize = minCellSize;
169 }
170
171 const QString expression = parameterAsExpression( parameters, QStringLiteral( "EXPRESSION" ), context );
172 const QString outputFile = parameterAsOutputLayer( parameters, QStringLiteral( "OUTPUT" ), context );
173 const QFileInfo fi( outputFile );
174 const QString outputFormat = QgsRasterFileWriter::driverForExtension( fi.suffix() );
175
176 double width = std::round( ( bbox.xMaximum() - bbox.xMinimum() ) / cellSize );
177 double height = std::round( ( bbox.yMaximum() - bbox.yMinimum() ) / cellSize );
178
179 QgsRasterCalculator calc( expression, outputFile, outputFormat, bbox, crs, width, height, entries, context.transformContext() );
180 QgsRasterCalculator::Result result = calc.processCalculation( feedback );
181 qDeleteAll( mLayers );
182 mLayers.clear();
183 switch ( result )
184 {
185 case QgsRasterCalculator::CreateOutputError:
186 throw QgsProcessingException( QObject::tr( "Error creating output file." ) );
187 case QgsRasterCalculator::InputLayerError:
188 throw QgsProcessingException( QObject::tr( "Error reading input layer." ) );
189 case QgsRasterCalculator::ParserError:
190 throw QgsProcessingException( QObject::tr( "Error parsing formula." ) );
191 case QgsRasterCalculator::MemoryError:
192 throw QgsProcessingException( QObject::tr( "Error allocating memory for result." ) );
193 case QgsRasterCalculator::BandError:
194 throw QgsProcessingException( QObject::tr( "Invalid band number for input." ) );
195 case QgsRasterCalculator::CalculationError:
196 throw QgsProcessingException( QObject::tr( "Error occurred while performing calculation." ) );
197 default:
198 break;
199 }
200
201 QVariantMap outputs;
202 outputs.insert( QStringLiteral( "OUTPUT" ), outputFile );
203 return outputs;
204}
205
206Qgis::ProcessingAlgorithmFlags QgsRasterCalculatorModelerAlgorithm::flags() const
207{
208 return QgsProcessingAlgorithm::flags() | Qgis::ProcessingAlgorithmFlag::HideFromToolbox;
209}
210
211QString QgsRasterCalculatorModelerAlgorithm::name() const
212{
213 return QStringLiteral( "modelerrastercalc" );
214}
215
216QString QgsRasterCalculatorModelerAlgorithm::displayName() const
217{
218 return QObject::tr( "Raster calculator" );
219}
220
221QStringList QgsRasterCalculatorModelerAlgorithm::tags() const
222{
223 return QObject::tr( "raster,calculator" ).split( ',' );
224}
225
226QString QgsRasterCalculatorModelerAlgorithm::group() const
227{
228 return QObject::tr( "Raster analysis" );
229}
230
231QString QgsRasterCalculatorModelerAlgorithm::groupId() const
232{
233 return QStringLiteral( "rasteranalysis" );
234}
235
236QgsRasterCalculatorModelerAlgorithm *QgsRasterCalculatorModelerAlgorithm::createInstance() const
237{
238 return new QgsRasterCalculatorModelerAlgorithm();
239}
240
241QVariantMap QgsRasterCalculatorModelerAlgorithm::processAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
242{
243 for ( QgsMapLayer *layer : std::as_const( mLayers ) )
244 {
245 layer->moveToThread( QThread::currentThread() );
246 }
247
248 QgsCoordinateReferenceSystem crs;
249 if ( parameters.value( QStringLiteral( "CRS" ) ).isValid() )
250 {
251 crs = parameterAsCrs( parameters, QStringLiteral( "CRS" ), context );
252 }
253 else
254 {
255 crs = mLayers.at( 0 )->crs();
256 }
257
258 QgsRectangle bbox;
259 if ( parameters.value( QStringLiteral( "EXTENT" ) ).isValid() )
260 {
261 bbox = parameterAsExtent( parameters, QStringLiteral( "EXTENT" ), context, crs );
262 }
263 else
264 {
265 bbox = QgsProcessingUtils::combineLayerExtents( mLayers, crs, context );
266 }
267
268 double minCellSize = 1e9;
269
270 QVector< QgsRasterCalculatorEntry > entries;
271 int n = 0;
272 for ( QgsMapLayer *layer : mLayers )
273 {
274 QgsRasterLayer *rLayer = qobject_cast<QgsRasterLayer *>( layer );
275 if ( !rLayer )
276 {
277 continue;
278 }
279
280 n++;
281 const int nBands = rLayer->dataProvider()->bandCount();
282 for ( int i = 0; i < nBands; ++i )
283 {
284 QgsRasterCalculatorEntry entry;
285 entry.ref = QStringLiteral( "%1@%2" ).arg( indexToName( n ) ).arg( i + 1 );
286 entry.raster = rLayer;
287 entry.bandNumber = i + 1;
288 entries << entry;
289 }
290
291 QgsRectangle ext = rLayer->extent();
292 if ( rLayer->crs() != crs )
293 {
294 QgsCoordinateTransform ct( rLayer->crs(), crs, context.transformContext() );
295 ext = ct.transformBoundingBox( ext );
296 }
297
298 double cellSize = ( ext.xMaximum() - ext.xMinimum() ) / rLayer->width();
299 if ( cellSize < minCellSize )
300 {
301 minCellSize = cellSize;
302 }
303 }
304
305 double cellSize = parameterAsDouble( parameters, QStringLiteral( "CELL_SIZE" ), context );
306 if ( cellSize == 0 )
307 {
308 cellSize = minCellSize;
309 }
310
311 const QString expression = parameterAsExpression( parameters, QStringLiteral( "EXPRESSION" ), context );
312 const QString outputFile = parameterAsOutputLayer( parameters, QStringLiteral( "OUTPUT" ), context );
313 const QFileInfo fi( outputFile );
314 const QString outputFormat = QgsRasterFileWriter::driverForExtension( fi.suffix() );
315
316 double width = std::round( ( bbox.xMaximum() - bbox.xMinimum() ) / cellSize );
317 double height = std::round( ( bbox.yMaximum() - bbox.yMinimum() ) / cellSize );
318
319 QgsRasterCalculator calc( expression, outputFile, outputFormat, bbox, crs, width, height, entries, context.transformContext() );
320 QgsRasterCalculator::Result result = calc.processCalculation( feedback );
321 qDeleteAll( mLayers );
322 mLayers.clear();
323 switch ( result )
324 {
325 case QgsRasterCalculator::CreateOutputError:
326 throw QgsProcessingException( QObject::tr( "Error creating output file." ) );
327 case QgsRasterCalculator::InputLayerError:
328 throw QgsProcessingException( QObject::tr( "Error reading input layer." ) );
329 case QgsRasterCalculator::ParserError:
330 throw QgsProcessingException( QObject::tr( "Error parsing formula." ) );
331 case QgsRasterCalculator::MemoryError:
332 throw QgsProcessingException( QObject::tr( "Error allocating memory for result." ) );
333 case QgsRasterCalculator::BandError:
334 throw QgsProcessingException( QObject::tr( "Invalid band number for input." ) );
335 case QgsRasterCalculator::CalculationError:
336 throw QgsProcessingException( QObject::tr( "Error occurred while performing calculation." ) );
337 default:
338 break;
339 }
340
341 QVariantMap outputs;
342 outputs.insert( QStringLiteral( "OUTPUT" ), outputFile );
343 return outputs;
344}
345
346QString QgsRasterCalculatorModelerAlgorithm::indexToName( int index ) const
347{
348 QString name;
349 int div = index;
350 int mod = 0;
351
352 while ( div > 0 )
353 {
354 mod = ( div - 1 ) % 26;
355 name = static_cast<char>( 65 + mod ) + name;
356 div = ( int )( ( div - mod ) / 26 );
357 }
358 return name;
359}
360
362
Qgis::ExpressionType::RasterCalculator
@ RasterCalculator
Raster calculator expression (since QGIS 3.34)
Qgis::ProcessingAlgorithmFlags
QFlags< ProcessingAlgorithmFlag > ProcessingAlgorithmFlags
Flags indicating how and when an algorithm operates and should be exposed to users.
Definition qgis.h:3071
Qgis::ProcessingAlgorithmFlag::HideFromToolbox
@ HideFromToolbox
Algorithm should be hidden from the toolbox.
Qgis::ProcessingAlgorithmFlag::HideFromModeler
@ HideFromModeler
Algorithm should be hidden from the modeler.
QgsCoordinateReferenceSystem
This class represents a coordinate reference system (CRS).
Definition qgscoordinatereferencesystem.h:212
QgsCoordinateTransform
Class for doing transforms between two map coordinate systems.
Definition qgscoordinatetransform.h:58
QgsMapLayer
Base class for all map layer types.
Definition qgsmaplayer.h:75
QgsMapLayer::name
QString name
Definition qgsmaplayer.h:79
QgsMapLayer::extent
virtual QgsRectangle extent() const
Returns the extent of the layer.
Definition qgsmaplayer.cpp:486
QgsMapLayer::crs
QgsCoordinateReferenceSystem crs
Definition qgsmaplayer.h:82
QgsMapLayer::clone
virtual QgsMapLayer * clone() const =0
Returns a new instance equivalent to this one except for the id which is still unique.
QgsProcessingAlgorithm::flags
virtual Qgis::ProcessingAlgorithmFlags flags() const
Returns the flags indicating how and when the algorithm operates and should be exposed to users.
Definition qgsprocessingalgorithm.cpp:90
QgsProcessingContext
Contains information about the context in which a processing algorithm is executed.
Definition qgsprocessingcontext.h:45
QgsProcessingContext::transformContext
QgsCoordinateTransformContext transformContext() const
Returns the coordinate transform context.
Definition qgsprocessingcontext.h:156
QgsProcessingException
Custom exception class for processing related exceptions.
Definition qgsexception.h:83
QgsProcessingFeedback
Base class for providing feedback from a processing algorithm.
Definition qgsprocessingfeedback.h:38
QgsProcessingFeedback::reportError
virtual void reportError(const QString &error, bool fatalError=false)
Reports that the algorithm encountered an error while executing.
Definition qgsprocessingfeedback.cpp:62
QgsProcessingParameterExpression
An expression parameter for processing algorithms.
Definition qgsprocessingparameters.h:2793
QgsProcessingParameterMultipleLayers
A parameter for processing algorithms which accepts multiple map layers.
Definition qgsprocessingparameters.h:2160
QgsProcessingParameterRasterDestination
A raster layer destination parameter, for specifying the destination path for a raster layer created ...
Definition qgsprocessingparameters.h:3450
QgsProcessingUtils::combineLayerExtents
static QgsRectangle combineLayerExtents(const QList< QgsMapLayer * > &layers, const QgsCoordinateReferenceSystem &crs, QgsProcessingContext &context)
Combines the extent of several map layers.
Definition qgsprocessingutils.cpp:1164
QgsRasterCalculatorEntry
Represents an individual raster layer/band number entry within a raster calculation.
Definition qgsrastercalculator.h:40
QgsRasterCalculatorEntry::raster
QgsRasterLayer * raster
Raster layer associated with entry.
Definition qgsrastercalculator.h:64
QgsRasterCalculatorEntry::bandNumber
int bandNumber
Band number for entry.
Definition qgsrastercalculator.h:69
QgsRasterCalculatorEntry::ref
QString ref
Name of entry.
Definition qgsrastercalculator.h:59
QgsRasterCalculator
Performs raster layer calculations.
Definition qgsrastercalculator.h:77
QgsRasterCalculator::Result
Result
Result of the calculation.
Definition qgsrastercalculator.h:82
QgsRasterCalculator::InputLayerError
@ InputLayerError
Error reading input layer.
Definition qgsrastercalculator.h:85
QgsRasterCalculator::ParserError
@ ParserError
Error parsing formula.
Definition qgsrastercalculator.h:87
QgsRasterCalculator::CalculationError
@ CalculationError
Error occurred while performing calculation.
Definition qgsrastercalculator.h:90
QgsRasterCalculator::MemoryError
@ MemoryError
Error allocating memory for result.
Definition qgsrastercalculator.h:88
QgsRasterCalculator::BandError
@ BandError
Invalid band number for input.
Definition qgsrastercalculator.h:89
QgsRasterCalculator::CreateOutputError
@ CreateOutputError
Error creating output data file.
Definition qgsrastercalculator.h:84
QgsRasterFileWriter::driverForExtension
static QString driverForExtension(const QString &extension)
Returns the GDAL driver name for a specified file extension.
Definition qgsrasterfilewriter.cpp:1066
QgsRasterInterface::bandCount
virtual int bandCount() const =0
Gets number of bands.
QgsRasterLayer
Represents a raster layer.
Definition qgsrasterlayer.h:78
QgsRasterLayer::dataProvider
QgsRasterDataProvider * dataProvider() override
Returns the source data provider.
Definition qgsrasterlayer.cpp:314
QgsRasterLayer::width
int width() const
Returns the width of the (unclipped) raster.
Definition qgsrasterlayer.cpp:2629
QgsRectangle
A rectangle specified with double values.
Definition qgsrectangle.h:42
QgsRectangle::xMinimum
double xMinimum() const
Returns the x minimum value (left side of rectangle).
Definition qgsrectangle.h:201
QgsRectangle::yMinimum
double yMinimum() const
Returns the y minimum value (bottom side of rectangle).
Definition qgsrectangle.h:211
QgsRectangle::xMaximum
double xMaximum() const
Returns the x maximum value (right side of rectangle).
Definition qgsrectangle.h:196
QgsRectangle::yMaximum
double yMaximum() const
Returns the y maximum value (top side of rectangle).
Definition qgsrectangle.h:206
crs
const QgsCoordinateReferenceSystem & crs
Definition qgswfsgetfeature.cpp:53
Generated on Thu Jun 27 2024 20:32:07 for QGIS API Documentation by doxygen 1.9.8