QGIS API Documentation 3.99.0-Master (d270888f95f)
Loading...
Searching...
No Matches
qgsrastercalcnode.cpp
Go to the documentation of this file.
1/***************************************************************************
2 qgsrastercalcnode.cpp
3 ---------------------
4 begin : October 2010
5 copyright : (C) 2010 by Marco Hugentobler
6 email : marco dot hugentobler at sourcepole dot ch
7 ***************************************************************************
8 * *
9 * This program is free software; you can redistribute it and/or modify *
10 * it under the terms of the GNU General Public License as published by *
11 * the Free Software Foundation; either version 2 of the License, or *
12 * (at your option) any later version. *
13 * *
14 ***************************************************************************/
15#include "qgsrastercalcnode.h"
16
17#include "qgsrasterblock.h"
18#include "qgsrastermatrix.h"
19
20#include <QString>
21
22using namespace Qt::StringLiterals;
23
25 : mNumber( number )
26{
27}
28
30 : mType( tMatrix )
31 , mMatrix( matrix )
32{
33}
34
36 : mType( tOperator )
37 , mLeft( left )
38 , mRight( right )
39 , mOperator( op )
40{
41}
42
43QgsRasterCalcNode::QgsRasterCalcNode( QString functionName, QVector<QgsRasterCalcNode *> functionArgs )
44 : mType( tFunction )
45 , mFunctionName( functionName )
46 , mFunctionArgs( functionArgs )
47{
48}
49
50QgsRasterCalcNode::QgsRasterCalcNode( const QString &rasterName )
51 : mType( tRasterRef )
52 , mRasterName( rasterName )
53{
54 if ( mRasterName.startsWith( '"' ) && mRasterName.endsWith( '"' ) )
55 mRasterName = mRasterName.mid( 1, mRasterName.size() - 2 );
56}
57
59{
60 delete mLeft;
61 delete mRight;
62 for ( int i = 0; i < mFunctionArgs.size(); ++i )
63 {
64 if ( mFunctionArgs.at( i ) )
65 delete mFunctionArgs.at( i );
66 }
67}
68
69bool QgsRasterCalcNode::calculate( QMap<QString, QgsRasterBlock *> &rasterData, QgsRasterMatrix &result, int row ) const
70{
71 //if type is raster ref: return a copy of the corresponding matrix
72
73 //if type is operator, call the proper matrix operations
74 if ( mType == tRasterRef )
75 {
76 const QMap<QString, QgsRasterBlock *>::iterator it = rasterData.find( mRasterName );
77 if ( it == rasterData.end() )
78 {
79 QgsDebugError( u"Error: could not find raster data for \"%1\""_s.arg( mRasterName ) );
80 return false;
81 }
82
83 const int nRows = ( row >= 0 ? 1 : ( *it )->height() );
84 const int startRow = ( row >= 0 ? row : 0 );
85 const int endRow = startRow + nRows;
86 const int nCols = ( *it )->width();
87 const int nEntries = nCols * nRows;
88 double *data = new double[nEntries];
89
90 //convert input raster values to double, also convert input no data to result no data
91
92 int outRow = 0;
93 bool isNoData = false;
94 for ( int dataRow = startRow; dataRow < endRow; ++dataRow, ++outRow )
95 {
96 for ( int dataCol = 0; dataCol < nCols; ++dataCol )
97 {
98 const double value = ( *it )->valueAndNoData( dataRow, dataCol, isNoData );
99 data[dataCol + nCols * outRow] = isNoData ? result.nodataValue() : value;
100 }
101 }
102 result.setData( nCols, nRows, data, result.nodataValue() );
103 return true;
104 }
105 else if ( mType == tOperator )
106 {
107 QgsRasterMatrix leftMatrix( result.nColumns(), result.nRows(), nullptr, result.nodataValue() );
108 QgsRasterMatrix rightMatrix( result.nColumns(), result.nRows(), nullptr, result.nodataValue() );
109
110 if ( !mLeft || !mLeft->calculate( rasterData, leftMatrix, row ) )
111 {
112 return false;
113 }
114 if ( mRight && !mRight->calculate( rasterData, rightMatrix, row ) )
115 {
116 return false;
117 }
118
119 switch ( mOperator )
120 {
121 case opPLUS:
122 leftMatrix.add( rightMatrix );
123 break;
124 case opMINUS:
125 leftMatrix.subtract( rightMatrix );
126 break;
127 case opMUL:
128 leftMatrix.multiply( rightMatrix );
129 break;
130 case opDIV:
131 leftMatrix.divide( rightMatrix );
132 break;
133 case opPOW:
134 leftMatrix.power( rightMatrix );
135 break;
136 case opEQ:
137 leftMatrix.equal( rightMatrix );
138 break;
139 case opNE:
140 leftMatrix.notEqual( rightMatrix );
141 break;
142 case opGT:
143 leftMatrix.greaterThan( rightMatrix );
144 break;
145 case opLT:
146 leftMatrix.lesserThan( rightMatrix );
147 break;
148 case opGE:
149 leftMatrix.greaterEqual( rightMatrix );
150 break;
151 case opLE:
152 leftMatrix.lesserEqual( rightMatrix );
153 break;
154 case opAND:
155 leftMatrix.logicalAnd( rightMatrix );
156 break;
157 case opOR:
158 leftMatrix.logicalOr( rightMatrix );
159 break;
160 case opMIN:
161 leftMatrix.min( rightMatrix );
162 break;
163 case opMAX:
164 leftMatrix.max( rightMatrix );
165 break;
166 case opSQRT:
167 leftMatrix.squareRoot();
168 break;
169 case opSIN:
170 leftMatrix.sinus();
171 break;
172 case opCOS:
173 leftMatrix.cosinus();
174 break;
175 case opTAN:
176 leftMatrix.tangens();
177 break;
178 case opASIN:
179 leftMatrix.asinus();
180 break;
181 case opACOS:
182 leftMatrix.acosinus();
183 break;
184 case opATAN:
185 leftMatrix.atangens();
186 break;
187 case opSIGN:
188 leftMatrix.changeSign();
189 break;
190 case opLOG:
191 leftMatrix.log();
192 break;
193 case opLOG10:
194 leftMatrix.log10();
195 break;
196 case opABS:
197 leftMatrix.absoluteValue();
198 break;
199 default:
200 return false;
201 }
202 const int newNColumns = leftMatrix.nColumns();
203 const int newNRows = leftMatrix.nRows();
204 result.setData( newNColumns, newNRows, leftMatrix.takeData(), leftMatrix.nodataValue() );
205 return true;
206 }
207 else if ( mType == tNumber )
208 {
209 const size_t nEntries = static_cast<size_t>( result.nColumns() * result.nRows() );
210 double *data = new double[nEntries];
211 std::fill( data, data + nEntries, mNumber );
212 result.setData( result.nColumns(), 1, data, result.nodataValue() );
213
214 return true;
215 }
216 else if ( mType == tMatrix )
217 {
218 const int nEntries = mMatrix->nColumns() * mMatrix->nRows();
219 double *data = new double[nEntries];
220 for ( int i = 0; i < nEntries; ++i )
221 {
222 data[i] = mMatrix->data()[i] == mMatrix->nodataValue() ? result.nodataValue() : mMatrix->data()[i];
223 }
224 result.setData( mMatrix->nColumns(), mMatrix->nRows(), data, result.nodataValue() );
225 return true;
226 }
227 else if ( mType == tFunction )
228 {
229 std::vector<std::unique_ptr<QgsRasterMatrix>> matrixContainer;
230 for ( int i = 0; i < mFunctionArgs.size(); ++i )
231 {
232 auto singleMatrix = std::make_unique<QgsRasterMatrix>( result.nColumns(), result.nRows(), nullptr, result.nodataValue() );
233 if ( !mFunctionArgs.at( i ) || !mFunctionArgs.at( i )->calculate( rasterData, *singleMatrix, row ) )
234 {
235 return false;
236 }
237 matrixContainer.emplace_back( std::move( singleMatrix ) );
238 }
239 evaluateFunction( matrixContainer, result );
240 return true;
241 }
242 return false;
243}
244
245QString QgsRasterCalcNode::toString( bool cStyle ) const
246{
247 QString result;
248 QString left;
249 QString right;
250 if ( mLeft )
251 left = mLeft->toString( cStyle );
252 if ( mRight )
253 right = mRight->toString( cStyle );
254
255 switch ( mType )
256 {
257 case tOperator:
258 switch ( mOperator )
259 {
260 case opPLUS:
261 result = u"( %1 + %2 )"_s.arg( left ).arg( right );
262 break;
263 case opMINUS:
264 result = u"( %1 - %2 )"_s.arg( left ).arg( right );
265 break;
266 case opSIGN:
267 result = u"-%1"_s.arg( left );
268 break;
269 case opMUL:
270 result = u"%1 * %2"_s.arg( left ).arg( right );
271 break;
272 case opDIV:
273 result = u"%1 / %2"_s.arg( left ).arg( right );
274 break;
275 case opPOW:
276 if ( cStyle )
277 result = u"pow( %1, %2 )"_s.arg( left ).arg( right );
278 else
279 result = u"%1^%2"_s.arg( left ).arg( right );
280 break;
281 case opEQ:
282 if ( cStyle )
283 result = u"( float ) ( %1 == %2 )"_s.arg( left ).arg( right );
284 else
285 result = u"%1 = %2"_s.arg( left ).arg( right );
286 break;
287 case opNE:
288 if ( cStyle )
289 result = u"( float ) ( %1 != %2 )"_s.arg( left ).arg( right );
290 else
291 result = u"%1 != %2"_s.arg( left ).arg( right );
292 break;
293 case opGT:
294 if ( cStyle )
295 result = u"( float ) ( %1 > %2 )"_s.arg( left ).arg( right );
296 else
297 result = u"%1 > %2"_s.arg( left ).arg( right );
298 break;
299 case opLT:
300 if ( cStyle )
301 result = u"( float ) ( %1 < %2 )"_s.arg( left ).arg( right );
302 else
303 result = u"%1 < %2"_s.arg( left ).arg( right );
304 break;
305 case opGE:
306 if ( cStyle )
307 result = u"( float ) ( %1 >= %2 )"_s.arg( left ).arg( right );
308 else
309 result = u"%1 >= %2"_s.arg( left ).arg( right );
310 break;
311 case opLE:
312 if ( cStyle )
313 result = u"( float ) ( %1 <= %2 )"_s.arg( left ).arg( right );
314 else
315 result = u"%1 <= %2"_s.arg( left ).arg( right );
316 break;
317 case opAND:
318 if ( cStyle )
319 result = u"( float ) ( %1 && %2 )"_s.arg( left ).arg( right );
320 else
321 result = u"%1 AND %2"_s.arg( left ).arg( right );
322 break;
323 case opOR:
324 if ( cStyle )
325 result = u"( float ) ( %1 || %2 )"_s.arg( left ).arg( right );
326 else
327 result = u"%1 OR %2"_s.arg( left ).arg( right );
328 break;
329 case opSQRT:
330 result = u"sqrt( %1 )"_s.arg( left );
331 break;
332 case opSIN:
333 result = u"sin( %1 )"_s.arg( left );
334 break;
335 case opCOS:
336 result = u"cos( %1 )"_s.arg( left );
337 break;
338 case opTAN:
339 result = u"tan( %1 )"_s.arg( left );
340 break;
341 case opASIN:
342 result = u"asin( %1 )"_s.arg( left );
343 break;
344 case opACOS:
345 result = u"acos( %1 )"_s.arg( left );
346 break;
347 case opATAN:
348 result = u"atan( %1 )"_s.arg( left );
349 break;
350 case opLOG:
351 result = u"log( %1 )"_s.arg( left );
352 break;
353 case opLOG10:
354 result = u"log10( %1 )"_s.arg( left );
355 break;
356 case opABS:
357 if ( cStyle )
358 result = u"fabs( %1 )"_s.arg( left );
359 else
360 // Call the floating point version
361 result = u"abs( %1 )"_s.arg( left );
362 break;
363 case opMIN:
364 if ( cStyle )
365 result = u"min( ( float ) ( %1 ), ( float ) ( %2 ) )"_s.arg( left ).arg( right );
366 else
367 result = u"min( %1, %2 )"_s.arg( left ).arg( right );
368 break;
369 case opMAX:
370 if ( cStyle )
371 result = u"max( ( float ) ( %1 ), ( float ) ( %2 ) )"_s.arg( left ).arg( right );
372 else
373 result = u"max( %1, %2 )"_s.arg( left ).arg( right );
374 break;
375 case opNONE:
376 break;
377 }
378 break;
379 case tRasterRef:
380 if ( cStyle )
381 result = u"( float ) \"%1\""_s.arg( mRasterName );
382 else
383 result = u"\"%1\""_s.arg( mRasterName );
384 break;
385 case tNumber:
386 result = QString::number( mNumber );
387 if ( cStyle )
388 {
389 result = u"( float ) %1"_s.arg( result );
390 }
391 break;
392 case tMatrix:
393 break;
394 case tFunction:
395 if ( mFunctionName == "if" )
396 {
397 const QString argOne = mFunctionArgs.at( 0 )->toString( cStyle );
398 const QString argTwo = mFunctionArgs.at( 1 )->toString( cStyle );
399 const QString argThree = mFunctionArgs.at( 2 )->toString( cStyle );
400 if ( cStyle )
401 result = u" ( %1 ) ? ( %2 ) : ( %3 ) "_s.arg( argOne, argTwo, argThree );
402 else
403 result = u"if( %1 , %2 , %3 )"_s.arg( argOne, argTwo, argThree );
404 }
405 break;
406 }
407 return result;
408}
409
410QList<const QgsRasterCalcNode *> QgsRasterCalcNode::findNodes( const QgsRasterCalcNode::Type type ) const
411{
412 QList<const QgsRasterCalcNode *> nodeList;
413 if ( mType == type )
414 nodeList.push_back( this );
415 if ( mLeft )
416 nodeList.append( mLeft->findNodes( type ) );
417 if ( mRight )
418 nodeList.append( mRight->findNodes( type ) );
419
420 for ( QgsRasterCalcNode *node : mFunctionArgs )
421 nodeList.append( node->findNodes( type ) );
422
423 return nodeList;
424}
425
426QgsRasterCalcNode *QgsRasterCalcNode::parseRasterCalcString( const QString &str, QString &parserErrorMsg )
427{
428 extern QgsRasterCalcNode *localParseRasterCalcString( const QString &str, QString &parserErrorMsg );
429 return localParseRasterCalcString( str, parserErrorMsg );
430}
431
433{
434 QStringList referencedRasters;
435
436 QStringList rasterRef = this->cleanRasterReferences();
437 for ( const auto &i : rasterRef )
438 {
439 if ( referencedRasters.contains( i.mid( 0, i.lastIndexOf( "@" ) ) ) )
440 continue;
441 referencedRasters << i.mid( 0, i.lastIndexOf( "@" ) );
442 }
443
444 return referencedRasters;
445}
446
448{
449 QStringList rasterReferences;
450 const QList<const QgsRasterCalcNode *> rasterRefNodes = this->findNodes( QgsRasterCalcNode::Type::tRasterRef );
451
452 for ( const QgsRasterCalcNode *r : rasterRefNodes )
453 {
454 QString layerRef( r->toString() );
455 if ( layerRef.at( 0 ) == "\""_L1 && layerRef.at( layerRef.size() - 1 ) == "\""_L1 )
456 {
457 layerRef.remove( 0, 1 );
458 layerRef.chop( 1 );
459 }
460 layerRef.remove( QChar( '\\' ), Qt::CaseInsensitive );
461 rasterReferences << layerRef;
462 }
463
464 return rasterReferences;
465}
466
467QgsRasterMatrix QgsRasterCalcNode::evaluateFunction( const std::vector<std::unique_ptr<QgsRasterMatrix>> &matrixVector, QgsRasterMatrix &result ) const
468{
469 if ( mFunctionName == "if" )
470 {
471 //scalar condition
472 if ( matrixVector.at( 0 )->isNumber() )
473 {
474 result = ( matrixVector.at( 0 )->data() ? *matrixVector.at( 1 ) : *matrixVector.at( 2 ) );
475 return result;
476 }
477 int nCols = matrixVector.at( 0 )->nColumns();
478 int nRows = matrixVector.at( 0 )->nRows();
479 int nEntries = nCols * nRows;
480 auto dataResult = std::make_unique<double[]>( nEntries );
481 double *dataResultRawPtr = dataResult.get();
482
483 double *condition = matrixVector.at( 0 )->data();
484 double *firstOption = matrixVector.at( 1 )->data();
485 double *secondOption = matrixVector.at( 2 )->data();
486
487 bool isFirstOptionNumber = matrixVector.at( 1 )->isNumber();
488 bool isSecondCOptionNumber = matrixVector.at( 2 )->isNumber();
489 double noDataValueCondition = matrixVector.at( 0 )->nodataValue();
490
491 for ( int i = 0; i < nEntries; ++i )
492 {
493 if ( condition[i] == noDataValueCondition )
494 {
495 dataResultRawPtr[i] = result.nodataValue();
496 continue;
497 }
498 else if ( condition[i] != 0 )
499 {
500 dataResultRawPtr[i] = isFirstOptionNumber ? firstOption[0] : firstOption[i];
501 continue;
502 }
503 dataResultRawPtr[i] = isSecondCOptionNumber ? secondOption[0] : secondOption[i];
504 }
505
506 result.setData( nCols, nRows, dataResult.release(), result.nodataValue() );
507 }
508 return result;
509}
QString toString(bool cStyle=false) const
Returns a string representation of the expression.
Operator
possible operators
QgsRasterCalcNode()=default
QList< const QgsRasterCalcNode * > findNodes(const QgsRasterCalcNode::Type type) const
Returns a list of nodes of a specific type.
static QgsRasterCalcNode * parseRasterCalcString(const QString &str, QString &parserErrorMsg)
bool calculate(QMap< QString, QgsRasterBlock * > &rasterData, QgsRasterMatrix &result, int row=-1) const
Calculates result of raster calculation (might be real matrix or single number).
QStringList referencedLayerNames() const
Returns a list of raster layer names that are referenced in the formula without the quotation marks.
QStringList cleanRasterReferences() const
Returns a list of raster layer references that are addressed in the formula, without quotation marks.
Type
defines possible types of node
Represents a matrix in a raster calculator operation.
bool equal(const QgsRasterMatrix &other)
bool max(const QgsRasterMatrix &other)
Calculates the maximum value between two matrices.
bool greaterEqual(const QgsRasterMatrix &other)
double * takeData()
Returns data and ownership.
bool notEqual(const QgsRasterMatrix &other)
bool greaterThan(const QgsRasterMatrix &other)
bool logicalOr(const QgsRasterMatrix &other)
bool lesserEqual(const QgsRasterMatrix &other)
bool absoluteValue()
Calculates the absolute value.
bool logicalAnd(const QgsRasterMatrix &other)
bool lesserThan(const QgsRasterMatrix &other)
bool add(const QgsRasterMatrix &other)
Adds another matrix to this one.
bool multiply(const QgsRasterMatrix &other)
int nRows() const
bool power(const QgsRasterMatrix &other)
double nodataValue() const
bool min(const QgsRasterMatrix &other)
Calculates the minimum value between two matrices.
bool divide(const QgsRasterMatrix &other)
int nColumns() const
void setData(int cols, int rows, double *data, double nodataValue)
bool subtract(const QgsRasterMatrix &other)
Subtracts another matrix from this one.
#define QgsDebugError(str)
Definition qgslogger.h:59