QGIS API Documentation  3.22.4-Białowieża (ce8e65e95e)
qgsrastermatrix.cpp
Go to the documentation of this file.
1 /***************************************************************************
2  qgsrastermatrix.cpp
3  -------------------
4  begin : 2010-10-23
5  copyright : (C) 20010 by Marco Hugentobler
6  email : marco dot hugentobler at sourcepole dot ch
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 "qgsrastermatrix.h"
19 #include <cstring>
20 #include <cmath>
21 #include <algorithm>
22 
23 QgsRasterMatrix::QgsRasterMatrix( int nCols, int nRows, double *data, double nodataValue )
24  : mColumns( nCols )
25  , mRows( nRows )
26  , mData( data )
27  , mNodataValue( nodataValue )
28 {
29 }
30 
32 {
33  operator=( m );
34 }
35 
37 {
38  delete[] mData;
39 }
40 
42 {
43  if ( this != &m )
44  {
45  delete[] mData;
46  mColumns = m.nColumns();
47  mRows = m.nRows();
48  const int nEntries = mColumns * mRows;
49  mData = new double[nEntries];
50  memcpy( mData, m.mData, sizeof( double ) * nEntries );
51  mNodataValue = m.nodataValue();
52  }
53  return *this;
54 }
55 
56 void QgsRasterMatrix::setData( int cols, int rows, double *data, double nodataValue )
57 {
58  delete[] mData;
59  mColumns = cols;
60  mRows = rows;
61  mData = data;
62  mNodataValue = nodataValue;
63 }
64 
66 {
67  double *data = mData;
68  mData = nullptr;
69  mColumns = 0;
70  mRows = 0;
71  return data;
72 }
73 
75 {
76  return twoArgumentOperation( opPLUS, other );
77 }
78 
80 {
81  return twoArgumentOperation( opMINUS, other );
82 }
83 
85 {
86  return twoArgumentOperation( opMUL, other );
87 }
88 
90 {
91  return twoArgumentOperation( opDIV, other );
92 }
93 
95 {
96  return twoArgumentOperation( opPOW, other );
97 }
98 
100 {
101  return twoArgumentOperation( opEQ, other );
102 }
103 
105 {
106  return twoArgumentOperation( opNE, other );
107 }
108 
110 {
111  return twoArgumentOperation( opGT, other );
112 }
113 
115 {
116  return twoArgumentOperation( opLT, other );
117 }
118 
120 {
121  return twoArgumentOperation( opGE, other );
122 }
123 
125 {
126  return twoArgumentOperation( opLE, other );
127 }
128 
130 {
131  return twoArgumentOperation( opAND, other );
132 }
133 
135 {
136  return twoArgumentOperation( opOR, other );
137 }
138 
140 {
141  return twoArgumentOperation( opMAX, other );
142 }
143 
145 {
146  return twoArgumentOperation( opMIN, other );
147 }
148 
150 {
151  return oneArgumentOperation( opSQRT );
152 }
153 
155 {
156  return oneArgumentOperation( opSIN );
157 }
158 
160 {
161  return oneArgumentOperation( opASIN );
162 }
163 
165 {
166  return oneArgumentOperation( opCOS );
167 }
168 
170 {
171  return oneArgumentOperation( opACOS );
172 }
173 
175 {
176  return oneArgumentOperation( opTAN );
177 }
178 
180 {
181  return oneArgumentOperation( opATAN );
182 }
183 
185 {
186  return oneArgumentOperation( opSIGN );
187 }
188 
190 {
191  return oneArgumentOperation( opLOG );
192 }
193 
195 {
196  return oneArgumentOperation( opLOG10 );
197 }
198 
200 {
201  return oneArgumentOperation( opABS );
202 }
203 
204 bool QgsRasterMatrix::oneArgumentOperation( OneArgOperator op )
205 {
206  if ( !mData )
207  {
208  return false;
209  }
210 
211  const int nEntries = mColumns * mRows;
212  double value;
213  for ( int i = 0; i < nEntries; ++i )
214  {
215  value = mData[i];
216  if ( value != mNodataValue )
217  {
218  switch ( op )
219  {
220  case opSQRT:
221  if ( value < 0 ) //no complex numbers
222  {
223  mData[i] = mNodataValue;
224  }
225  else
226  {
227  mData[i] = std::sqrt( value );
228  }
229  break;
230  case opSIN:
231  mData[i] = std::sin( value );
232  break;
233  case opCOS:
234  mData[i] = std::cos( value );
235  break;
236  case opTAN:
237  mData[i] = std::tan( value );
238  break;
239  case opASIN:
240  mData[i] = std::asin( value );
241  break;
242  case opACOS:
243  mData[i] = std::acos( value );
244  break;
245  case opATAN:
246  mData[i] = std::atan( value );
247  break;
248  case opSIGN:
249  mData[i] = -value;
250  break;
251  case opLOG:
252  if ( value <= 0 )
253  {
254  mData[i] = mNodataValue;
255  }
256  else
257  {
258  mData[i] = ::log( value );
259  }
260  break;
261  case opLOG10:
262  if ( value <= 0 )
263  {
264  mData[i] = mNodataValue;
265  }
266  else
267  {
268  mData[i] = ::log10( value );
269  }
270  break;
271  case opABS:
272  mData[i] = ::fabs( value );
273  break;
274  }
275  }
276  }
277  return true;
278 }
279 
280 double QgsRasterMatrix::calculateTwoArgumentOp( TwoArgOperator op, double arg1, double arg2 ) const
281 {
282  switch ( op )
283  {
284  case opPLUS:
285  return arg1 + arg2;
286  case opMINUS:
287  return arg1 - arg2;
288  case opMUL:
289  return arg1 * arg2;
290  case opDIV:
291  if ( arg2 == 0 )
292  {
293  return mNodataValue;
294  }
295  else
296  {
297  return arg1 / arg2;
298  }
299  case opPOW:
300  if ( !testPowerValidity( arg1, arg2 ) )
301  {
302  return mNodataValue;
303  }
304  else
305  {
306  return std::pow( arg1, arg2 );
307  }
308  case opEQ:
309  return ( arg1 == arg2 ? 1.0 : 0.0 );
310  case opNE:
311  return ( arg1 == arg2 ? 0.0 : 1.0 );
312  case opGT:
313  return ( arg1 > arg2 ? 1.0 : 0.0 );
314  case opLT:
315  return ( arg1 < arg2 ? 1.0 : 0.0 );
316  case opGE:
317  return ( arg1 >= arg2 ? 1.0 : 0.0 );
318  case opLE:
319  return ( arg1 <= arg2 ? 1.0 : 0.0 );
320  case opAND:
321  return ( arg1 && arg2 ? 1.0 : 0.0 );
322  case opOR:
323  return ( arg1 || arg2 ? 1.0 : 0.0 );
324  case opMAX:
325  return std::max( arg1, arg2 );
326  case opMIN:
327  return std::min( arg1, arg2 );
328  }
329  return mNodataValue;
330 }
331 
332 bool QgsRasterMatrix::twoArgumentOperation( TwoArgOperator op, const QgsRasterMatrix &other )
333 {
334  if ( isNumber() && other.isNumber() ) //operation on two 1x1 matrices
335  {
336  //operations with nodata values always generate nodata
337  if ( mData[0] == mNodataValue || other.number() == other.nodataValue() )
338  {
339  mData[0] = mNodataValue;
340  }
341  else
342  {
343  mData[0] = calculateTwoArgumentOp( op, mData[0], other.number() );
344  }
345  return true;
346  }
347 
348  //two matrices
349  if ( !isNumber() && !other.isNumber() )
350  {
351  double *matrix = other.mData;
352  const int nEntries = mColumns * mRows;
353  double value1, value2;
354 
355  for ( int i = 0; i < nEntries; ++i )
356  {
357  value1 = mData[i];
358  value2 = matrix[i];
359  if ( value1 == mNodataValue || value2 == other.mNodataValue )
360  {
361  mData[i] = mNodataValue;
362  }
363  else
364  {
365  mData[i] = calculateTwoArgumentOp( op, value1, value2 );
366  }
367  }
368  return true;
369  }
370 
371  //this matrix is a single number and the other one a real matrix
372  if ( isNumber() )
373  {
374  double *matrix = other.mData;
375  const int nEntries = other.nColumns() * other.nRows();
376  const double value = mData[0];
377  delete[] mData;
378  mData = new double[nEntries];
379  mColumns = other.nColumns();
380  mRows = other.nRows();
381  mNodataValue = other.nodataValue();
382 
383  if ( value == mNodataValue )
384  {
385  for ( int i = 0; i < nEntries; ++i )
386  {
387  mData[i] = mNodataValue;
388  }
389  return true;
390  }
391 
392  for ( int i = 0; i < nEntries; ++i )
393  {
394  if ( matrix[i] == other.mNodataValue )
395  {
396  mData[i] = mNodataValue;
397  continue;
398  }
399 
400  mData[i] = calculateTwoArgumentOp( op, value, matrix[i] );
401  }
402  return true;
403  }
404  else //this matrix is a real matrix and the other a number
405  {
406  const double value = other.number();
407  const int nEntries = mColumns * mRows;
408 
409  if ( other.number() == other.mNodataValue )
410  {
411  for ( int i = 0; i < nEntries; ++i )
412  {
413  mData[i] = mNodataValue;
414  }
415  return true;
416  }
417 
418  for ( int i = 0; i < nEntries; ++i )
419  {
420  if ( mData[i] == mNodataValue )
421  {
422  continue;
423  }
424 
425  mData[i] = calculateTwoArgumentOp( op, mData[i], value );
426  }
427  return true;
428  }
429 }
430 
431 bool QgsRasterMatrix::testPowerValidity( double base, double power ) const
432 {
433  return !( ( base == 0 && power < 0 ) || ( base < 0 && ( power - std::floor( power ) ) > 0 ) );
434 }
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)
QgsRasterMatrix()=default
Takes ownership of data array.
bool logicalOr(const QgsRasterMatrix &other)
bool lesserEqual(const QgsRasterMatrix &other)
bool absoluteValue()
Calculates the absolute value.
bool logicalAnd(const QgsRasterMatrix &other)
bool isNumber() const
Returns true if matrix is 1x1 (=scalar number)
bool lesserThan(const QgsRasterMatrix &other)
bool add(const QgsRasterMatrix &other)
Adds another matrix to this one.
bool multiply(const QgsRasterMatrix &other)
QgsRasterMatrix & operator=(const QgsRasterMatrix &m)
int nRows() const
bool power(const QgsRasterMatrix &other)
double nodataValue() const
bool min(const QgsRasterMatrix &other)
Calculates the minimum value between two matrices.
double number() const
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.