![]() |
TerraLib 4.1
|
00001 /************************************************************************************ 00002 TerraLib - a library for developing GIS applications. 00003 Copyright © 2001-2004 INPE and Tecgraf/PUC-Rio. 00004 00005 This code is part of the TerraLib library. 00006 This library is free software; you can redistribute it and/or 00007 modify it under the terms of the GNU Lesser General Public 00008 License as published by the Free Software Foundation; either 00009 version 2.1 of the License, or (at your option) any later version. 00010 00011 You should have received a copy of the GNU Lesser General Public 00012 License along with this library. 00013 00014 The authors reassure the license terms regarding the warranties. 00015 They specifically disclaim any warranties, including, but not limited to, 00016 the implied warranties of merchantability and fitness for a particular purpose. 00017 The library provided hereunder is on an "as is" basis, and the authors have no 00018 obligation to provide maintenance, support, updates, enhancements, or modifications. 00019 In no event shall INPE and Tecgraf / PUC-Rio be held liable to any party for direct, 00020 indirect, special, incidental, or consequential damages arising out of the use 00021 of this library and its documentation. 00022 *************************************************************************************/ 00023 00028 #ifndef __TERRALIB_INTERNAL_INTERPOLATION_H 00029 #define __TERRALIB_INTERNAL_INTERPOLATION_H 00030 00031 #include "TeTheme.h" 00032 00033 #include "TeFunctionsDefines.h" 00034 00035 00036 /* 00037 * NOTE: You should call TeCellInterpolate or TeRasterInterpolate with the required parameters. Don't use 00038 * the class TeInterpolationAlgorithms, it is for internal use only! 00039 * 00040 */ 00041 00042 00057 enum TeInterpolationAlgorithm { TeNNInterpolation, TeAvgInterpolation, TeDistWeightAvgInterpolation, 00058 TeAvgInBoxInterpolation, TeDistWeightAvgInBoxInterpolation }; 00059 00061 00076 template<class ADAPTATIVEKDTREE> class TeInterpolationAlgorithms 00077 { 00078 protected: 00079 00080 typedef typename ADAPTATIVEKDTREE::kdKey kdKey; 00081 typedef typename ADAPTATIVEKDTREE::kdDataItem kdDataItem; 00082 typedef typename ADAPTATIVEKDTREE::kdNode kdNode; 00083 00085 const ADAPTATIVEKDTREE& kd_; 00086 00087 public: 00088 00090 TeInterpolationAlgorithms(const ADAPTATIVEKDTREE& kd) 00091 : kd_(kd) 00092 { 00093 } 00094 00096 ~TeInterpolationAlgorithms() 00097 { 00098 } 00099 00101 double nearestNeighbor(const kdKey& key) 00102 { 00103 vector<kdDataItem> report; 00104 vector<double> sqrDists; 00105 00106 fillNNVector(report, 1); 00107 00108 kd_.nearestNeighborSearch(key, report, sqrDists, 1); 00109 00110 if(sqrDists[0] >= TeMAXFLOAT) 00111 return -TeMAXFLOAT; 00112 else 00113 return report[0].value(); 00114 } 00115 00117 double avgNearestNeighbor(const kdKey& key, const unsigned int& numberOfNeighbors) 00118 { 00119 vector<kdDataItem> report; 00120 vector<double> sqrDists; 00121 00122 fillNNVector(report, numberOfNeighbors); 00123 00124 kd_.nearestNeighborSearch(key, report, sqrDists, numberOfNeighbors); 00125 00126 double numElements = numberOfNeighbors; 00127 00128 double sum = 0.0; 00129 00130 for(unsigned int i = 0; i < numberOfNeighbors; ++i) 00131 { 00132 if(sqrDists[i] >= TeMAXFLOAT) 00133 { 00134 --numElements; 00135 continue; 00136 } 00137 00138 if(sqrDists[i] == 0.0) 00139 return report[i].value(); 00140 00141 sum += report[i].value(); 00142 } 00143 00144 if(numElements > 0.0) 00145 return sum / numElements; 00146 else 00147 return -TeMAXFLOAT; 00148 } 00149 00151 double distWeightAvgNearestNeighbor(const kdKey& key, const unsigned int& numberOfNeighbors, const int& powValue) 00152 { 00153 vector<kdDataItem> report; 00154 vector<double> sqrDists; 00155 00156 fillNNVector(report, numberOfNeighbors); 00157 00158 kd_.nearestNeighborSearch(key, report, sqrDists, numberOfNeighbors); 00159 00160 double num = 0.0; 00161 double den = 0.0; 00162 00163 if(powValue == 1.0) 00164 { 00165 for(unsigned int i = 0; i < numberOfNeighbors; ++i) 00166 { 00167 if(sqrDists[i] >= TeMAXFLOAT) 00168 continue; 00169 00170 if(sqrDists[i] == 0.0) 00171 return report[i].value(); 00172 00173 double wi = 1.0 / (sqrt(sqrDists[i])); 00174 00175 num += (wi * report[i].value()); 00176 den += wi; 00177 } 00178 } 00179 else if(powValue == 2.0) 00180 { 00181 for(unsigned int i = 0; i < numberOfNeighbors; ++i) 00182 { 00183 if(sqrDists[i] >= TeMAXFLOAT) 00184 continue; 00185 00186 if(sqrDists[i] == 0.0) 00187 return report[i].value(); 00188 00189 double wi = 1.0 / (sqrDists[i]); 00190 00191 num += (wi * report[i].value()); 00192 den += wi; 00193 } 00194 } 00195 else 00196 { 00197 for(unsigned int i = 0; i < numberOfNeighbors; ++i) 00198 { 00199 if(sqrDists[i] >= TeMAXFLOAT) 00200 continue; 00201 00202 if(sqrDists[i] == 0.0) 00203 return report[i].value(); 00204 00205 double wi = 1.0 / pow(sqrt(sqrDists[i]), powValue); 00206 00207 num += (wi * report[i].value()); 00208 den += wi; 00209 } 00210 } 00211 00212 if(den != 0.0) 00213 return num / den; 00214 else 00215 return -TeMAXFLOAT; 00216 } 00217 00219 double boxAvg(const TeBox& box) 00220 { 00221 vector<kdNode*> report; 00222 00223 kd_.search(box, report); 00224 00225 unsigned int numberOfNodes = report.size(); 00226 00227 unsigned int numElements = 0; 00228 00229 double sum = 0.0; 00230 00231 for(unsigned int i = 0; i < numberOfNodes; ++i) 00232 { 00233 unsigned int nodeSize = report[i]->getData().size(); 00234 00235 for(unsigned int j = 0; j < nodeSize; ++j) 00236 { 00237 if(TeIntersects(report[i]->getData()[j].location(), box)) 00238 { 00239 sum += report[i]->getData()[j].value(); 00240 ++numElements; 00241 } 00242 } 00243 } 00244 00245 if(numElements > 0) 00246 return sum / double(numElements); 00247 else 00248 return -TeMAXFLOAT; 00249 } 00250 00252 double boxDistWeightAvg(const kdKey& key, const TeBox& box, const int& powValue) 00253 { 00254 vector<kdNode*> report; 00255 00256 kd_.search(box, report); 00257 00258 unsigned int numberOfNodes = report.size(); 00259 00260 double num = 0.0; 00261 double den = 0.0; 00262 00263 for(unsigned int i = 0; i < numberOfNodes; ++i) 00264 { 00265 unsigned int nodeSize = report[i]->getData().size(); 00266 00267 for(unsigned int j = 0; j < nodeSize; ++j) 00268 { 00269 if(TeIntersects(report[i]->getData()[j].location(), box)) 00270 { 00271 double wi = 1 / pow(TeDistance(key, report[i]->getData()[j].location()), powValue); 00272 00273 num += wi * report[i]->getData()[j].value(); 00274 den += wi; 00275 } 00276 } 00277 } 00278 00279 if(den != 0.0) 00280 return num / den; 00281 else 00282 return -TeMAXFLOAT; 00283 } 00284 00285 protected: 00286 00288 void fillNNVector(vector<kdDataItem>& report, const unsigned int& numberOfNeighbors) const 00289 { 00290 for(unsigned int i = 0; i < numberOfNeighbors; ++i) 00291 { 00292 TeCoord2D c(TeMAXFLOAT, TeMAXFLOAT); 00293 kdDataItem item(c, -TeMAXFLOAT); 00294 00295 report.push_back(item); 00296 } 00297 } 00298 00299 private: 00300 00302 TeInterpolationAlgorithms(const TeInterpolationAlgorithms& rhs); 00303 00305 TeInterpolationAlgorithms& operator=(const TeInterpolationAlgorithms& rhs); 00306 }; 00307 00308 00309 00311 /* 00312 \param inputTheme The theme with samples to be used: samples are treated as a layer of points with a column with a value (double or integer) 00313 \param inputAttrTableName Theme table with samples 00314 \param sampleColumnName Table column with sample values: must be int or double values 00315 \param outputTheme The theme where output will be stored: this theme must already exist 00316 \param outputAttrTableName Theme table where output will be stored: this table must already exist 00317 \param outputColumnName Table column where data will be stored: this column must already exist 00318 \param algorithm Type of algorithm to be used: see enum TeInterpolationAlgorithm 00319 \param numberOfNeighbors Number of nearest neighbor to be used in interpolation 00320 \param powValue Pow parameter for inverse distance function 00321 \param boxRatio Ratio of the box used to search samples 00322 \return TRUE if there isn't an error during interpolation otherwise return FALSE 00323 */ 00324 TLFUNCTIONS_DLL bool TeCellInterpolate(TeTheme* inputTheme, const string& inputAttrTableName, const string& sampleColumnName, 00325 TeTheme* outputTheme, const string& outputAttrTableName, const string& outputColumnName, 00326 const TeInterpolationAlgorithm& algorithm = TeNNInterpolation, 00327 const unsigned int& numberOfNeighbors = 1, const int& powValue = 1, const double& boxRatio = 1.0); 00328 00330 /* 00331 \param inputTheme The theme with samples to be used: samples are treated as a layer of points with a column with a value (double or integer) 00332 \param inputAttrTableName Theme table with samples 00333 \param sampleColumnName Table column with sample values 00334 \param outputRaster A pointer to a raster where interpolated data will be stored: this raster mus already exist 00335 \param band The band to be interpolated 00336 \param algorithm Type of algorithm to be used: see enum TeInterpolationAlgorithm 00337 \param numberOfNeighbors Number of nearest neighbor to be used in interpolation 00338 \param powValue Pow parameter for inverse distance function 00339 \param boxRatio Ratio of the box used to search samples 00340 \return TRUE if there isn't an error during interpolation otherwise return FALSE 00341 */ 00342 TLFUNCTIONS_DLL bool TeRasterInterpolate(TeTheme* inputTheme, const string& inputAttrTableName, const string& sampleColumnName, 00343 TeRaster* outputRaster, const int& band = 0, 00344 const TeInterpolationAlgorithm& algorithm = TeNNInterpolation, 00345 const unsigned int& numberOfNeighbors = 1, const int& powValue = 1, const double& boxRatio = 1.0); 00346 00347 00348 00354 #endif // __TERRALIB_INTERNAL_INTERPOLATION_H 00355 00356 00357 00358