QGIS API Documentation  3.37.0-Master (a5b4d9743e8)
qgsderivativefilter.cpp
Go to the documentation of this file.
1 /***************************************************************************
2  qgsderivativefilter.cpp - description
3  -----------------------
4  begin : August 7th, 2009
5  copyright : (C) 2009 by Marco Hugentobler
6  email : marco dot hugentobler at karto dot baug dot ethz 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 "qgsderivativefilter.h"
19 
20 QgsDerivativeFilter::QgsDerivativeFilter( const QString &inputFile, const QString &outputFile, const QString &outputFormat )
21  : QgsNineCellFilter( inputFile, outputFile, outputFormat )
22 {
23 
24 }
25 
26 float QgsDerivativeFilter::calcFirstDerX( float *x11, float *x21, float *x31, float *x12, float *x22, float *x32, float *x13, float *x23, float *x33 ) const
27 {
28  //the basic formula would be simple, but we need to test for nodata values...
29  //return (( (*x31 - *x11) + 2 * (*x32 - *x12) + (*x33 - *x13) ) / (8 * mCellSizeX));
30 
31  int weight = 0;
32  double sum = 0;
33 
34  //first row
35  if ( *x31 != mInputNodataValue && *x11 != mInputNodataValue ) //the normal case
36  {
37  sum += ( *x31 - *x11 );
38  weight += 2;
39  }
40  else if ( *x31 == mInputNodataValue && *x11 != mInputNodataValue && *x21 != mInputNodataValue ) //probably 3x3 window is at the border
41  {
42  sum += ( *x21 - *x11 );
43  weight += 1;
44  }
45  else if ( *x11 == mInputNodataValue && *x31 != mInputNodataValue && *x21 != mInputNodataValue ) //probably 3x3 window is at the border
46  {
47  sum += ( *x31 - *x21 );
48  weight += 1;
49  }
50 
51  //second row
52  if ( *x32 != mInputNodataValue && *x12 != mInputNodataValue ) //the normal case
53  {
54  sum += 2 * ( *x32 - *x12 );
55  weight += 4;
56  }
57  else if ( *x32 == mInputNodataValue && *x12 != mInputNodataValue && *x22 != mInputNodataValue )
58  {
59  sum += 2 * ( *x22 - *x12 );
60  weight += 2;
61  }
62  else if ( *x12 == mInputNodataValue && *x32 != mInputNodataValue && *x22 != mInputNodataValue )
63  {
64  sum += 2 * ( *x32 - *x22 );
65  weight += 2;
66  }
67 
68  //third row
69  if ( *x33 != mInputNodataValue && *x13 != mInputNodataValue ) //the normal case
70  {
71  sum += ( *x33 - *x13 );
72  weight += 2;
73  }
74  else if ( *x33 == mInputNodataValue && *x13 != mInputNodataValue && *x23 != mInputNodataValue )
75  {
76  sum += ( *x23 - *x13 );
77  weight += 1;
78  }
79  else if ( *x13 == mInputNodataValue && *x33 != mInputNodataValue && *x23 != mInputNodataValue )
80  {
81  sum += ( *x33 - *x23 );
82  weight += 1;
83  }
84 
85  if ( weight == 0 )
86  {
87  return mOutputNodataValue;
88  }
89 
90  return sum / ( weight * mCellSizeX ) * mZFactor;
91 }
92 
93 float QgsDerivativeFilter::calcFirstDerY( float *x11, float *x21, float *x31, float *x12, float *x22, float *x32, float *x13, float *x23, float *x33 ) const
94 {
95  //the basic formula would be simple, but we need to test for nodata values...
96  //return (((*x11 - *x13) + 2 * (*x21 - *x23) + (*x31 - *x33)) / ( 8 * mCellSizeY));
97 
98  double sum = 0;
99  int weight = 0;
100 
101  //first row
102  if ( *x11 != mInputNodataValue && *x13 != mInputNodataValue ) //normal case
103  {
104  sum += ( *x11 - *x13 );
105  weight += 2;
106  }
107  else if ( *x11 == mInputNodataValue && *x13 != mInputNodataValue && *x12 != mInputNodataValue )
108  {
109  sum += ( *x12 - *x13 );
110  weight += 1;
111  }
112  else if ( *x31 == mInputNodataValue && *x11 != mInputNodataValue && *x12 != mInputNodataValue )
113  {
114  sum += ( *x11 - *x12 );
115  weight += 1;
116  }
117 
118  //second row
119  if ( *x21 != mInputNodataValue && *x23 != mInputNodataValue )
120  {
121  sum += 2 * ( *x21 - *x23 );
122  weight += 4;
123  }
124  else if ( *x21 == mInputNodataValue && *x23 != mInputNodataValue && *x22 != mInputNodataValue )
125  {
126  sum += 2 * ( *x22 - *x23 );
127  weight += 2;
128  }
129  else if ( *x23 == mInputNodataValue && *x21 != mInputNodataValue && *x22 != mInputNodataValue )
130  {
131  sum += 2 * ( *x21 - *x22 );
132  weight += 2;
133  }
134 
135  //third row
136  if ( *x31 != mInputNodataValue && *x33 != mInputNodataValue )
137  {
138  sum += ( *x31 - *x33 );
139  weight += 2;
140  }
141  else if ( *x31 == mInputNodataValue && *x33 != mInputNodataValue && *x32 != mInputNodataValue )
142  {
143  sum += ( *x32 - *x33 );
144  weight += 1;
145  }
146  else if ( *x33 == mInputNodataValue && *x31 != mInputNodataValue && *x32 != mInputNodataValue )
147  {
148  sum += ( *x31 - *x32 );
149  weight += 1;
150  }
151 
152  if ( weight == 0 )
153  {
154  return mOutputNodataValue;
155  }
156 
157  return sum / ( weight * mCellSizeY ) * mZFactor;
158 }
159 
160 
161 
162 
163 
QgsDerivativeFilter(const QString &inputFile, const QString &outputFile, const QString &outputFormat)
float calcFirstDerX(float *x11, float *x21, float *x31, float *x12, float *x22, float *x32, float *x13, float *x23, float *x33) const
Calculates the first order derivative in x-direction according to Horn (1981)
float calcFirstDerY(float *x11, float *x21, float *x31, float *x12, float *x22, float *x32, float *x13, float *x23, float *x33) const
Calculates the first order derivative in y-direction according to Horn (1981)
Base class for raster analysis methods that work with a 3x3 cell filter and calculate the value of ea...
float mInputNodataValue
The nodata value of the input layer.
float mOutputNodataValue
The nodata value of the output layer.
double mZFactor
Scale factor for z-value if x-/y- units are different to z-units (111120 for degree->meters and 37040...