Quantum GIS API Documentation
1.7.4
|
00001 /*************************************************************************** 00002 qgsderivativefilter.cpp - description 00003 ----------------------- 00004 begin : August 7th, 2009 00005 copyright : (C) 2009 by Marco Hugentobler 00006 email : marco dot hugentobler at karto dot baug dot ethz dot ch 00007 ***************************************************************************/ 00008 00009 /*************************************************************************** 00010 * * 00011 * This program is free software; you can redistribute it and/or modify * 00012 * it under the terms of the GNU General Public License as published by * 00013 * the Free Software Foundation; either version 2 of the License, or * 00014 * (at your option) any later version. * 00015 * * 00016 ***************************************************************************/ 00017 00018 #include "qgsderivativefilter.h" 00019 00020 QgsDerivativeFilter::QgsDerivativeFilter( const QString& inputFile, const QString& outputFile, const QString& outputFormat ): \ 00021 QgsNineCellFilter( inputFile, outputFile, outputFormat ) 00022 { 00023 00024 } 00025 00026 QgsDerivativeFilter::~QgsDerivativeFilter() 00027 { 00028 00029 } 00030 00031 float QgsDerivativeFilter::calcFirstDerX( float* x11, float* x21, float* x31, float* x12, float* x22, float* x32, float* x13, float* x23, float* x33 ) 00032 { 00033 //the basic formula would be simple, but we need to test for nodata values... 00034 //return (( (*x31 - *x11) + 2 * (*x32 - *x12) + (*x33 - *x13) ) / (8 * mCellSizeX)); 00035 00036 int weight = 0; 00037 double sum = 0; 00038 00039 //first row 00040 if ( *x31 != mInputNodataValue && *x11 != mInputNodataValue ) //the normal case 00041 { 00042 sum += ( *x31 - *x11 ); 00043 weight += 2; 00044 } 00045 else if ( *x31 == mInputNodataValue && *x11 != mInputNodataValue && *x21 != mInputNodataValue ) //probably 3x3 window is at the border 00046 { 00047 sum += ( *x21 - *x11 ); 00048 weight += 1; 00049 } 00050 else if ( *x11 == mInputNodataValue && *x31 != mInputNodataValue && *x21 != mInputNodataValue ) //probably 3x3 window is at the border 00051 { 00052 sum += ( *x31 - *x21 ); 00053 weight += 1; 00054 } 00055 00056 //second row 00057 if ( *x32 != mInputNodataValue && *x12 != mInputNodataValue ) //the normal case 00058 { 00059 sum += 2 * ( *x32 - *x12 ); 00060 weight += 4; 00061 } 00062 else if ( *x32 == mInputNodataValue && *x12 != mInputNodataValue && *x22 != mInputNodataValue ) 00063 { 00064 sum += 2 * ( *x22 - *x12 ); 00065 weight += 2; 00066 } 00067 else if ( *x12 == mInputNodataValue && *x32 != mInputNodataValue && *x22 != mInputNodataValue ) 00068 { 00069 sum += 2 * ( *x32 - *x22 ); 00070 weight += 2; 00071 } 00072 00073 //third row 00074 if ( *x33 != mInputNodataValue && *x13 != mInputNodataValue ) //the normal case 00075 { 00076 sum += ( *x33 - *x13 ); 00077 weight += 2; 00078 } 00079 else if ( *x33 == mInputNodataValue && *x13 != mInputNodataValue && *x23 != mInputNodataValue ) 00080 { 00081 sum += ( *x23 - *x13 ); 00082 weight += 1; 00083 } 00084 else if ( *x13 == mInputNodataValue && *x33 != mInputNodataValue && *x23 != mInputNodataValue ) 00085 { 00086 sum += ( *x33 - *x23 ); 00087 weight += 1; 00088 } 00089 00090 if ( weight == 0 ) 00091 { 00092 return mOutputNodataValue; 00093 } 00094 00095 return sum / ( weight * mCellSizeX ); 00096 } 00097 00098 float QgsDerivativeFilter::calcFirstDerY( float* x11, float* x21, float* x31, float* x12, float* x22, float* x32, float* x13, float* x23, float* x33 ) 00099 { 00100 //the basic formula would be simple, but we need to test for nodata values... 00101 //return (((*x11 - *x13) + 2 * (*x21 - *x23) + (*x31 - *x33)) / ( 8 * mCellSizeY)); 00102 00103 double sum = 0; 00104 int weight = 0; 00105 00106 //first row 00107 if ( *x11 != mInputNodataValue && *x13 != mInputNodataValue ) //normal case 00108 { 00109 sum += ( *x11 - *x13 ); 00110 weight += 2; 00111 } 00112 else if ( *x11 == mInputNodataValue && *x13 != mInputNodataValue && *x12 != mInputNodataValue ) 00113 { 00114 sum += ( *x12 - *x13 ); 00115 weight += 1; 00116 } 00117 else if ( *x31 == mInputNodataValue && *x11 != mInputNodataValue && *x12 != mInputNodataValue ) 00118 { 00119 sum += ( *x11 - *x12 ); 00120 weight += 1; 00121 } 00122 00123 //second row 00124 if ( *x21 != mInputNodataValue && *x23 != mInputNodataValue ) 00125 { 00126 sum += 2 * ( *x21 - *x23 ); 00127 weight += 4; 00128 } 00129 else if ( *x21 == mInputNodataValue && *x23 != mInputNodataValue && *x22 != mInputNodataValue ) 00130 { 00131 sum += 2 * ( *x22 - *x23 ); 00132 weight += 2; 00133 } 00134 else if ( *x23 == mInputNodataValue && *x21 != mInputNodataValue && *x22 != mInputNodataValue ) 00135 { 00136 sum += 2 * ( *x21 - *x22 ); 00137 weight += 2; 00138 } 00139 00140 //third row 00141 if ( *x31 != mInputNodataValue && *x33 != mInputNodataValue ) 00142 { 00143 sum += ( *x31 - *x33 ); 00144 weight += 2; 00145 } 00146 else if ( *x31 == mInputNodataValue && *x33 != mInputNodataValue && *x32 != mInputNodataValue ) 00147 { 00148 sum += ( *x32 - *x33 ); 00149 weight += 1; 00150 } 00151 else if ( *x33 == mInputNodataValue && *x31 != mInputNodataValue && *x32 != mInputNodataValue ) 00152 { 00153 sum += ( *x31 - *x32 ); 00154 weight += 1; 00155 } 00156 00157 if ( weight == 0 ) 00158 { 00159 return mOutputNodataValue; 00160 } 00161 00162 return sum / ( weight * mCellSizeY ); 00163 } 00164 00165 00166 00167 00168