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