src/terralib/rp/Functions.cpp
Go to the documentation of this file.
1 /* Copyright (C) 2008 National Institute For Space Research (INPE) - Brazil.
2 
3  This file is part of the TerraLib - a Framework for building GIS enabled applications.
4 
5  TerraLib is free software: you can redistribute it and/or modify
6  it under the terms of the GNU Lesser General Public License as published by
7  the Free Software Foundation, either version 3 of the License,
8  or (at your option) any later version.
9 
10  TerraLib is distributed in the hope that it will be useful,
11  but WITHOUT ANY WARRANTY; without even the implied warranty of
12  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13  GNU Lesser General Public License for more details.
14 
15  You should have received a copy of the GNU Lesser General Public License
16  along with TerraLib. See COPYING. If not, write to
17  TerraLib Team at <terralib-team@terralib.org>.
18  */
19 
20 /*!
21  \file terralib/rp/Functions.cpp
22 
23  \brief Raster Processing functions.
24 */
25 
26 // TerraLib
27 
28 #ifdef WIN32
29 #define NOMINMAX
30 #endif
31 
32 #include "Functions.h"
33 
34 #include "../common/SystemApplicationSettings.h"
35 #include "../common/Version.h"
36 #include "../core/utils/Platform.h"
37 #include "../dataaccess/dataset/DataSetType.h"
38 #include "../dataaccess/datasource/DataSourceFactory.h"
39 #include "../dataaccess/utils/Utils.h"
40 #include "../datatype/Enums.h"
41 #include "../raster/BandProperty.h"
42 #include "../raster/Grid.h"
43 #include "../raster/RasterFactory.h"
44 #include "../raster/RasterProperty.h"
45 #include "../raster/RasterSummaryManager.h"
46 #include "../raster/RasterIterator.h"
47 #include "../raster/Utils.h"
48 #include "../geometry/Point.h"
49 #include "../common/MatrixUtils.h"
50 #include "../common/progress/TaskProgress.h"
51 #include "Exception.h"
52 #include "Macros.h"
53 #include "RasterHandler.h"
54 #include "Contrast.h"
55 
56 // Boost
57 #include <boost/filesystem.hpp>
58 #include <boost/numeric/ublas/io.hpp>
59 #include <boost/numeric/ublas/matrix.hpp>
60 #include <boost/property_tree/ptree.hpp>
61 #include <boost/property_tree/json_parser.hpp>
62 #include <boost/shared_ptr.hpp>
63 #include <boost/shared_array.hpp>
64 #include <boost/lexical_cast.hpp>
65 #include <boost/thread.hpp>
66 #include <boost/numeric/ublas/lu.hpp>
67 #include <boost/numeric/ublas/matrix.hpp>
68 
69 // STL
70 #include <cstring>
71 #include <string>
72 #include <limits>
73 #include <map>
74 #include <memory>
75 #include <algorithm>
76 
77 #ifndef M_PI
78  #define M_PI 3.14159265358979323846
79 #endif
80 
81 namespace te
82 {
83  namespace rp
84  {
85  struct GetMeanValueThreadParams
86  {
87  te::rst::Band const * m_inputBandPtr;
88  Matrix< bool >* m_rasterBlocksStatusPtr;
89  double* m_meanValuePtr;
90  double* m_pixelsNumberValuePtr;
91  bool m_returnStatus;
92  boost::mutex* m_mutexPtr;
93  };
94 
95  struct GetStdDevValueThreadParams
96  {
97  te::rst::Band const * m_inputBandPtr;
98  Matrix< bool >* m_rasterBlocksStatusPtr;
99  double* m_stdDevValuePtr;
100  double* m_pixelsNumberValuePtr;
101  double m_meanValue;
102  bool m_returnStatus;
103  boost::mutex* m_mutexPtr;
104  };
105 
106  struct GetCovarianceValueThreadParams
107  {
108  te::rst::Band const * m_inputBand1Ptr;
109  te::rst::Band const * m_inputBand2Ptr;
110  Matrix< bool >* m_rasterBlocksStatusPtr;
111  double* m_covarianceValuePtr;
112  double* m_pixelsNumberValuePtr;
113  bool m_returnStatus;
114  double m_mean1Value;
115  double m_mean2Value;
116  boost::mutex* m_mutexPtr;
117  };
118 
119  struct RemapValuesThreadParams
120  {
121  te::rst::Raster const * m_inputRasterPtr;
122  std::vector< unsigned int > const * m_inputRasterBandsPtr;
123  te::rst::Raster* m_outputRasterPtr;
124  std::vector< unsigned int > const * m_outputRasterBandsPtr;
125  boost::numeric::ublas::matrix< double > const * m_remapMatrixPtr;
126  Matrix< bool > * m_rasterBlocksStatusPtr;
127  bool m_returnStatus;
128  boost::mutex* m_mutexPtr;
129  };
130 
131  bool CreateNewRaster( const te::rst::Grid& rasterGrid,
132  const std::vector< te::rst::BandProperty* >& bandsProperties,
133  const std::string& outDataSetName,
134  const std::string& dataSourceType,
135  RasterHandler& outRasterHandler )
136  {
137  // Creating a new datasource
138 
139  std::unique_ptr< te::da::DataSource > dataSourcePtr(
140  te::da::DataSourceFactory::make( dataSourceType, ("file://" + outDataSetName)));
141  if( dataSourcePtr.get() == nullptr ) return false;
142 
143  RasterHandler internalRasterHandler;
144 
145  if( CreateNewRaster( rasterGrid, bandsProperties,
146  outDataSetName, *dataSourcePtr, internalRasterHandler ) )
147  {
148  std::unique_ptr< te::da::DataSource > dummyDataSourcePtr;
149  std::unique_ptr< te::da::DataSourceTransactor > transactorPtr;
150  std::unique_ptr< te::da::DataSet > dataSetPtr;
151  std::unique_ptr< te::rst::Raster > rasterPtr;
152 
153  internalRasterHandler.release( dummyDataSourcePtr, transactorPtr,
154  dataSetPtr, rasterPtr );
155 
156  outRasterHandler.reset( dataSourcePtr.release(), transactorPtr.release(),
157  dataSetPtr.release(), rasterPtr.release() );
158 
159  return true;
160  }
161  else
162  {
163  return false;
164  }
165  }
166 
167  bool CreateNewRaster( const te::rst::Grid& rasterGrid,
168  const std::vector< te::rst::BandProperty* >& bandsProperties,
169  const std::string& outDataSetName,
170  te::da::DataSource& outDataSource,
171  RasterHandler& outRasterHandler )
172  {
173  // Defining the raster properties
174 
175  std::unique_ptr< te::rst::RasterProperty > rasterPropertyPtr(
176  new te::rst::RasterProperty( new te::rst::Grid( rasterGrid ),
177  bandsProperties, std::map< std::string, std::string >(),
178  false, 0, nullptr ) );
179 
180  // acquiring a transactor instance
181 
182  std::unique_ptr< te::da::DataSourceTransactor > transactorPtr(
183  outDataSource.getTransactor() );
184 
185  if( transactorPtr.get() == nullptr )
186  {
187  return false;
188  }
189 
190  // Creating a data set instance
191 
192  if( transactorPtr->dataSetExists( outDataSetName ) )
193  {
194  transactorPtr->dropDataSet( outDataSetName );
195  }
196 
197  std::unique_ptr< te::da::DataSetType > dataSetTypePtr(
198  new te::da::DataSetType( outDataSetName ) );
199  dataSetTypePtr->add( rasterPropertyPtr.release() );
200 
201  try
202  {
203  transactorPtr->createDataSet( dataSetTypePtr.get(),
204  std::map< std::string, std::string >() );
205  }
206  catch( ... )
207  {
208  return false;
209  }
210 
211  if( ! transactorPtr->dataSetExists( outDataSetName ) ) return false;
212 
213  std::unique_ptr< te::da::DataSet > dataSetPtr( transactorPtr->getDataSet(
214  outDataSetName, te::common::FORWARDONLY, true, te::common::RWAccess ) );
215 
216  if( dataSetPtr.get() == nullptr )
217  {
218  return false;
219  }
220 
221  // Creating a raster instance
222 
223  std::unique_ptr< te::rst::Raster > rasterPtr( dataSetPtr->getRaster( 0 ) );
224 
225  if( rasterPtr.get() )
226  {
227  outRasterHandler.reset( transactorPtr.release(), dataSetPtr.release(), rasterPtr.release() );
228  return true;
229  }
230  else
231  {
232  return false;
233  }
234  }
235 
236  bool TERPEXPORT CreateNewRaster( const te::rst::Grid& rasterGrid,
237  const std::vector< te::rst::BandProperty* >& bandsProperties,
238  const std::map< std::string, std::string>& rasterInfo,
239  const std::string& rasterType,
240  std::unique_ptr< te::rst::Raster >& outRasterPtr )
241  {
242  try
243  {
244  outRasterPtr.reset( te::rst::RasterFactory::make( rasterType,
245  new te::rst::Grid( rasterGrid ), bandsProperties, rasterInfo, nullptr, nullptr ) );
246  }
247  catch( const te::common::Exception&/* exc */)
248  {
249  return false;
250  }
251 
252  if( outRasterPtr.get() == nullptr )
253  {
254  return false;
255  }
256 
257  return true;
258  }
259 
260  bool CreateNewMemRaster( const te::rst::Grid& rasterGrid,
261  std::vector< te::rst::BandProperty* > bandsProperties,
262  RasterHandler& outRasterHandler )
263  {
264  std::string dataSetName = std::string( "createNewMemRaster" ) +
265  boost::lexical_cast< std::string >( &outRasterHandler );
266 
267  return CreateNewRaster( rasterGrid, bandsProperties, dataSetName,
268  "MEM", outRasterHandler );
269  }
270 
271  bool CreateNewGdalRaster( const te::rst::Grid& rasterGrid,
272  std::vector< te::rst::BandProperty* > bandsProperties,
273  const std::string& fileName,
274  RasterHandler& outRasterHandler )
275  {
276  std::unique_ptr< te::rst::Raster > outRasterPtr;
277  if( CreateNewGdalRaster( rasterGrid, bandsProperties, fileName, outRasterPtr ) )
278  {
279  outRasterHandler.reset( outRasterPtr.release() );
280  return true;
281  }
282  else
283  {
284  outRasterHandler.reset();
285  return false;
286  }
287  }
288 
289  bool CreateNewGdalRaster( const te::rst::Grid& rasterGrid,
290  std::vector< te::rst::BandProperty* > bandsProperties,
291  const std::string& fileName,
292  std::unique_ptr< te::rst::Raster >& outRasterPtr )
293  {
294  outRasterPtr.reset();
295 
296  std::map< std::string, std::string > rInfo;
297  rInfo[ "URI" ] = fileName;
298 
299  outRasterPtr.reset( te::rst::RasterFactory::make(
300  "GDAL",
301  new te::rst::Grid( rasterGrid ),
302  bandsProperties,
303  rInfo,
304  nullptr,
305  nullptr ) );
306 
307  if( outRasterPtr.get() )
308  {
309  return true;
310  }
311  else
312  {
313  return false;
314  }
315  }
316 
317  bool TERPEXPORT Copy2DiskRaster( const te::rst::Raster& inputRaster,
318  const std::string& fileName )
319  {
320  if( !(inputRaster.getAccessPolicy() & te::common::RAccess ) )
321  {
322  return false;
323  };
324 
325  const unsigned int nBands = static_cast<unsigned int>(inputRaster.getNumberOfBands());
326  const unsigned int nCols = inputRaster.getNumberOfColumns();
327  const unsigned int nRows = inputRaster.getNumberOfRows();
328  unsigned int bandIdx =0;
329  unsigned int col = 0;
330  unsigned int row = 0;
331 
332  std::vector< te::rst::BandProperty* > bandsProperties;
333  for( bandIdx = 0 ; bandIdx < nBands ; ++bandIdx )
334  {
335  bandsProperties.push_back( new te::rst::BandProperty(
336  *( inputRaster.getBand( bandIdx )->getProperty() ) ) );
337  }
338 
339  RasterHandler outRasterHandler;
340 
341  if( ! CreateNewGdalRaster( *( inputRaster.getGrid() ), bandsProperties,
342  fileName, outRasterHandler ) )
343  {
344  return false;
345  }
346 
347  double value = 0;
348 
349  for( bandIdx = 0 ; bandIdx < nBands ; ++bandIdx )
350  {
351  const te::rst::Band& inBand = *inputRaster.getBand( bandIdx );
352  te::rst::Band& outBand = *outRasterHandler.getRasterPtr()->getBand( bandIdx );
353 
354  for( row = 0 ; row < nRows ; ++row )
355  {
356  for( col = 0 ; col < nCols ; ++col )
357  {
358  inBand.getValue( col, row, value );
359  outBand.setValue( col, row, value );
360  }
361  }
362  }
363 
364  return true;
365  }
366 
367  void Convert2DoublesVector( void* inputVector, const int inputVectorDataType,
368  unsigned int inputVectorSize, double* outputVector )
369  {
370  switch( inputVectorDataType )
371  {
372  case te::dt::CHAR_TYPE :
373  {
374  char* vPtr = static_cast<char*>(inputVector);
375  for( unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
376  outputVector[ idx ] = static_cast<double>(vPtr[ idx ]);
377  break;
378  }
379  case te::dt::BIT_TYPE :
380  case te::dt::UCHAR_TYPE :
381  {
382  unsigned char* vPtr = static_cast<unsigned char*>(inputVector);
383  for( unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
384  outputVector[ idx ] = static_cast<double>(vPtr[ idx ]);
385  break;
386  }
387  case te::dt::INT16_TYPE :
388  {
389  short int* vPtr = static_cast<short int*>(inputVector);
390  for( unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
391  outputVector[ idx ] = static_cast<double>(vPtr[ idx ]);
392  break;
393  }
394  case te::dt::CINT16_TYPE :
395  {
396  std::complex< short int >* vPtr = static_cast<std::complex< short int >*>(inputVector);
397  for( unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
398  outputVector[ idx ] = static_cast<double>(vPtr[ idx ].real());
399  break;
400  }
401  case te::dt::UINT16_TYPE :
402  {
403  unsigned short int* vPtr = static_cast<unsigned short int*>(inputVector);
404  for( unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
405  outputVector[ idx ] = static_cast<double>(vPtr[ idx ]);
406  break;
407  }
408  case te::dt::INT32_TYPE :
409  {
410  int* vPtr = static_cast<int*>(inputVector);
411  for( unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
412  outputVector[ idx ] = static_cast<double>(vPtr[ idx ]);
413  break;
414  }
415  case te::dt::CINT32_TYPE :
416  {
417  std::complex< int >* vPtr = static_cast<std::complex< int >*>(inputVector);
418  for( unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
419  outputVector[ idx ] = static_cast<double>(vPtr[ idx ].real());
420  break;
421  }
422  case te::dt::UINT32_TYPE :
423  {
424  unsigned int* vPtr = static_cast<unsigned int*>(inputVector);
425  for( unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
426  outputVector[ idx ] = static_cast<double>(vPtr[ idx ]);
427  break;
428  }
429  case te::dt::INT64_TYPE :
430  {
431  long int* vPtr = static_cast<long int*>(inputVector);
432  for( unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
433  outputVector[ idx ] = static_cast<double>(vPtr[ idx ]);
434  break;
435  }
436  case te::dt::UINT64_TYPE :
437  {
438  unsigned long int* vPtr = static_cast<unsigned long int*>(inputVector);
439  for( unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
440  outputVector[ idx ] = static_cast<double>(vPtr[ idx ]);
441  break;
442  }
443  case te::dt::FLOAT_TYPE :
444  {
445  float* vPtr = static_cast<float*>(inputVector);
446  for( unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
447  outputVector[ idx ] = static_cast<double>(vPtr[ idx ]);
448  break;
449  }
450  case te::dt::CFLOAT_TYPE :
451  {
452  std::complex< float >* vPtr = static_cast<std::complex< float >*>(inputVector);
453  for( unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
454  outputVector[ idx ] = static_cast<double>(vPtr[ idx ].real());
455  break;
456  }
457  case te::dt::DOUBLE_TYPE :
458  {
459  memcpy( outputVector, inputVector, inputVectorSize * sizeof( double ) );
460  break;
461  }
462  case te::dt::CDOUBLE_TYPE :
463  {
464  std::complex< double >* vPtr = static_cast<std::complex< double >*>(inputVector);
465  for( unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
466  outputVector[ idx ] = static_cast<double>(vPtr[ idx ].real());
467  break;
468  }
469  default :
470  throw te::rp::Exception( "Invalid data type" );
471  };
472  }
473 
474  void ConvertDoublesVector( double* inputVector,
475  unsigned int inputVectorSize, const int outputVectorDataType,
476  void* outputVector )
477  {
478  switch( outputVectorDataType )
479  {
480  case te::dt::CHAR_TYPE :
481  {
482  char* vPtr = static_cast<char*>(outputVector);
483  for(unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
484  vPtr[ idx ] = static_cast<char>(inputVector[ idx ]);
485  break;
486  }
487  case te::dt::BIT_TYPE :
488  case te::dt::UCHAR_TYPE :
489  {
490  unsigned char* vPtr = static_cast<unsigned char*>(outputVector);
491  for(unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
492  vPtr[ idx ] = static_cast<unsigned char>(inputVector[ idx ]);
493  break;
494  }
495  case te::dt::INT16_TYPE :
496  {
497  short int* vPtr = static_cast<short int*>(outputVector);
498  for(unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
499  vPtr[ idx ] = static_cast<short int>(inputVector[ idx ]);
500  break;
501  }
502  case te::dt::CINT16_TYPE :
503  {
504  std::complex< short int >* vPtr = static_cast<std::complex< short int >*>(outputVector);
505  for(unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
506  vPtr[ idx ]= ( static_cast<short int>(inputVector[ idx ]) );
507  break;
508  }
509  case te::dt::UINT16_TYPE :
510  {
511  unsigned short int* vPtr = static_cast<unsigned short int*>(outputVector);
512  for(unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
513  vPtr[ idx ] = static_cast<unsigned short int>(inputVector[ idx ]);
514  break;
515  }
516  case te::dt::INT32_TYPE :
517  {
518  int* vPtr = static_cast<int*>(outputVector);
519  for(unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
520  vPtr[ idx ] = static_cast<int>(inputVector[ idx ]);
521  break;
522  }
523  case te::dt::CINT32_TYPE :
524  {
525  std::complex< int >* vPtr = static_cast<std::complex< int >*>(outputVector);
526  for(unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
527  vPtr[ idx ] = ( static_cast<int>(inputVector[ idx ]) );
528  break;
529  }
530  case te::dt::UINT32_TYPE :
531  {
532  unsigned int* vPtr = static_cast<unsigned int*>(outputVector);
533  for(unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
534  vPtr[ idx ] = static_cast<unsigned int>(inputVector[ idx ]);
535  break;
536  }
537  case te::dt::INT64_TYPE :
538  {
539  long int* vPtr = static_cast<long int*>(outputVector);
540  for(unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
541  vPtr[ idx ] = static_cast<long int>(inputVector[ idx ]);
542  break;
543  }
544  case te::dt::UINT64_TYPE :
545  {
546  unsigned long int* vPtr = static_cast<unsigned long int*>(outputVector);
547  for(unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
548  vPtr[ idx ] = static_cast<unsigned long int>(inputVector[ idx ]);
549  break;
550  }
551  case te::dt::FLOAT_TYPE :
552  {
553  float* vPtr = static_cast<float*>(outputVector);
554  for(unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
555  vPtr[ idx ] = static_cast<float>(inputVector[ idx ]);
556  break;
557  }
558  case te::dt::CFLOAT_TYPE :
559  {
560  std::complex< float >* vPtr = static_cast<std::complex< float >*>(outputVector);
561  for(unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
562  vPtr[ idx ] = ( static_cast<float>(inputVector[ idx ]) );
563  break;
564  }
565  case te::dt::DOUBLE_TYPE :
566  {
567  memcpy( outputVector, inputVector, inputVectorSize * sizeof( double ) );
568  break;
569  }
570  case te::dt::CDOUBLE_TYPE :
571  {
572  std::complex< double >* vPtr = static_cast<std::complex< double >*>(outputVector);
573  for( unsigned int idx = 0 ; idx < inputVectorSize ; ++idx )
574  vPtr[ idx ] = ( static_cast<double>(inputVector[ idx ]) );
575  break;
576  }
577  default :
578  throw te::rp::Exception( "Invalid data type" );
579  };
580  }
581 
582  std::string GetSensorFilename()
583  {
584  std::string datalocalPath = te::core::GetAppLocalDataLocation();
585  std::string app = te::common::SystemApplicationSettings::getInstance().getValue("Application.Name");
586  std::string Organization = te::common::SystemApplicationSettings::getInstance().getValue("Application.Organization");
587  std::string vers = te::common::Version::asString();
588 
589  static std::string sens_file = datalocalPath + "/" + Organization + "/" + app + "-" + vers + "/SpectralSensor.json";
590 
591  return sens_file;
592  }
593 
594  std::vector<std::string> GetBandNames()
595  {
596  static std::vector<std::string> bandNames;
597  bandNames.clear();
598 
599  bandNames.push_back("CBERS2_CCD_1_BLUE");
600  bandNames.push_back("CBERS2_CCD_2_GREEN");
601  bandNames.push_back("CBERS2_CCD_3_RED");
602  bandNames.push_back("CBERS2_CCD_4_NIR");
603  bandNames.push_back("CBERS2_CCD_5_PAN");
604 
605  bandNames.push_back("LANDSAT5_TM_1_BLUE");
606  bandNames.push_back("LANDSAT5_TM_2_GREEN");
607  bandNames.push_back("LANDSAT5_TM_3_RED");
608  bandNames.push_back("LANDSAT5_TM_4_NIR");
609  bandNames.push_back("LANDSAT5_TM_5_SWIR");
610  bandNames.push_back("LANDSAT5_TM_6_THERMAL");
611  bandNames.push_back("LANDSAT5_TM_7_MIR");
612 
613  bandNames.push_back("LANDSAT7_ETM+_1_BLUE");
614  bandNames.push_back("LANDSAT7_ETM+_2_GREEN");
615  bandNames.push_back("LANDSAT7_ETM+_3_RED");
616  bandNames.push_back("LANDSAT7_ETM+_4_NIR");
617  bandNames.push_back("LANDSAT7_ETM+_5_SWIR");
618  bandNames.push_back("LANDSAT7_ETM+_6_THERMAL");
619  bandNames.push_back("LANDSAT7_ETM+_7_MIR");
620  bandNames.push_back("LANDSAT7_ETM+_8_PAN");
621 
622  bandNames.push_back("WV2_MUL_1_COASTAL");
623  bandNames.push_back("WV2_MUL_2_BLUE");
624  bandNames.push_back("WV2_MUL_3_GREEN");
625  bandNames.push_back("WV2_MUL_4_YELLOW");
626  bandNames.push_back("WV2_MUL_5_RED");
627  bandNames.push_back("WV2_MUL_6_REDEDGE");
628  bandNames.push_back("WV2_MUL_7_NIR1");
629  bandNames.push_back("WV2_MUL_8_NIR2");
630  bandNames.push_back("WV2_PAN");
631 
632  bandNames.push_back("AWIFS_RED");
633  bandNames.push_back("AWIFS_GREEN");
634  bandNames.push_back("AWIFS_NIR");
635  bandNames.push_back("AWIFS_MIR");
636 
637  bandNames.push_back("RAPID_EYE_BLUE");
638  bandNames.push_back("RAPID_EYE_GREEN");
639  bandNames.push_back("RAPID_EYE_RED");
640  bandNames.push_back("RAPID_EYE_REDEDGE");
641  bandNames.push_back("RAPID_EYE_NIR");
642 
643  bandNames.push_back("LISS3_GREEN");
644  bandNames.push_back("LISS3_RED");
645  bandNames.push_back("LISS3_NEAR_INFRARED");
646  bandNames.push_back("LISS3_SHORTWAVE_INFRARED");
647 
648  bandNames.push_back("DMC_NIR");
649  bandNames.push_back("DMC_RED");
650  bandNames.push_back("DMC_GREEN");
651 
652  bandNames.push_back("LANDSAT8-OLI1_COASTAL_AEROSOL");
653  bandNames.push_back("LANDSAT8-OLI2_BLUE");
654  bandNames.push_back("LANDSAT8-OLI3_GREEN");
655  bandNames.push_back("LANDSAT8-OLI4_RED");
656  bandNames.push_back("LANDSAT8-OLI5_NIR");
657  bandNames.push_back("LANDSAT8-OLI6_SWIR1");
658  bandNames.push_back("LANDSAT8-OLI7_SWIR2");
659  bandNames.push_back("LANDSAT8-OLI8_PAN");
660  bandNames.push_back("LANDSAT8-OLI9_CIRRUS");
661 
662  bandNames.push_back("LANDSAT8-TIRS1");
663  bandNames.push_back("LANDSAT8-TIRS2");
664 
665  bandNames.push_back("CBERS4_PAN_PAN");
666  bandNames.push_back("CBERS4_PAN_GREEN");
667  bandNames.push_back("CBERS4_PAN_RED");
668  bandNames.push_back("CBERS4_PAN_NIR");
669 
670  bandNames.push_back("CBERS4_MUX_BLUE");
671  bandNames.push_back("CBERS4_MUX_GREEN");
672  bandNames.push_back("CBERS4_MUX_RED");
673  bandNames.push_back("CBERS4_MUX_NIR");
674 
675  bandNames.push_back("CBERS4_IRS_PAN");
676  bandNames.push_back("CBERS4_IRS_SWIR_B10");
677  bandNames.push_back("CBERS4_IRS_SWIR_B11");
678  bandNames.push_back("CBERS4_IRS_TH");
679 
680  bandNames.push_back("CBERS4_WFI_BLUE");
681  bandNames.push_back("CBERS4_WFI_GREEN");
682  bandNames.push_back("CBERS4_WFI_RED");
683  bandNames.push_back("CBERS4_WFI_NIR");
684 
685  return bandNames;
686  }
687 
688  std::map<std::string, SpectralSensorParams > getSensorParams()
689  {
690  std::string sens_file = GetSensorFilename();
691 
692  static std::map<std::string, SpectralSensorParams > SensorParams;
693  SensorParams.clear();
694 
695  if (!boost::filesystem::exists(sens_file)) //SpectralSensor.json doesn't exists, then uses the original sensors
696  {
697  std::vector<std::string> bn = GetBandNames();
698  for (size_t i = 0; i < bn.size(); i++)
699  SensorParams[bn[i]] = GetSpectralBandInfo(bn[i]);
700  }
701  else
702  {
703  boost::property_tree::ptree pt;
704  boost::property_tree::json_parser::read_json(sens_file, pt);
705 
706  for(boost::property_tree::ptree::value_type &v: pt.get_child("Spectral Sensors"))
707  {
708  std::string name = v.second.get<std::string>("name");
709  std::string band = v.second.get<std::string>("band");
710  std::string lower = v.second.get<std::string>("lower");
711  std::string upper = v.second.get<std::string>("upper");
712  std::string min = v.second.get<std::string>("min");
713  std::string max = v.second.get<std::string>("max");
714  SensorParams[name] = SpectralSensorParams(std::atoi(band.c_str()), std::atof(lower.c_str()), std::atof(upper.c_str()), std::atof(min.c_str()), std::atof(max.c_str()));
715  }
716  }
717 
718  return SensorParams;
719  }
720 
721  void SaveSensorParams(std::map<std::string, SpectralSensorParams > &params)
722  {
723  std::string filename = GetSensorFilename();
724 
725  boost::property_tree::ptree pt;
726  boost::property_tree::ptree sensors;
727 
728  for (std::map<std::string, te::rp::SpectralSensorParams > ::iterator it = params.begin(); it != params.end(); ++it)
729  {
730  boost::property_tree::ptree sensor;
731 
732  sensor.put("name", it->first);
733  sensor.put("band", it->second.m_band);
734  sensor.put("lower", it->second.m_lower);
735  sensor.put("upper", it->second.m_upper);
736  sensor.put("min", it->second.m_min);
737  sensor.put("max", it->second.m_max);
738 
739  sensors.push_back(std::make_pair("Spectral Sensor", sensor));
740  }
741 
742  pt.add_child("Spectral Sensors", sensors);
743 
744  boost::property_tree::json_parser::write_json(filename, pt);
745 
746  }
747 
748  SpectralSensorParams GetSpectralBandInfo(std::string bandName)
749  {
750  static std::map<std::string, struct SpectralSensorParams > BandInfo;
751  BandInfo.clear();
752 
753  BandInfo["CBERS2_CCD_1_BLUE"] = SpectralSensorParams(1, 0.45, 0.52, 0.0, 255.0);
754  BandInfo["CBERS2_CCD_2_GREEN"] = SpectralSensorParams(2, 0.52, 0.59, 0.0, 255.0);
755  BandInfo["CBERS2_CCD_3_RED"] = SpectralSensorParams(3, 0.63, 0.69, 0.0, 255.0);
756  BandInfo["CBERS2_CCD_4_NIR"] = SpectralSensorParams(4, 0.77, 0.89, 0.0, 255.0);
757  BandInfo["CBERS2_CCD_5_PAN"] = SpectralSensorParams(5, 0.51, 0.73, 0.0, 255.0);
758 
759  BandInfo["LANDSAT5_TM_1_BLUE"] = SpectralSensorParams(1, 0.45, 0.52, 0.0, 255.0);
760  BandInfo["LANDSAT5_TM_2_GREEN"] = SpectralSensorParams(2, 0.52, 0.60, 0.0, 255.0);
761  BandInfo["LANDSAT5_TM_3_RED"] = SpectralSensorParams(3, 0.63, 0.69, 0.0, 255.0);
762  BandInfo["LANDSAT5_TM_4_NIR"] = SpectralSensorParams(4, 0.76, 0.90, 0.0, 255.0);
763  BandInfo["LANDSAT5_TM_5_SWIR"] = SpectralSensorParams(5, 1.55, 1.75, 0.0, 255.0);
764  BandInfo["LANDSAT5_TM_6_THERMAL"] = SpectralSensorParams(6, 10.40, 12.50, 0.0, 255.0);
765  BandInfo["LANDSAT5_TM_7_MIR"] = SpectralSensorParams(7, 2.08, 2.35, 0.0, 255.0);
766 
767  BandInfo["LANDSAT7_ETM+_1_BLUE"] = SpectralSensorParams(1, 0.45, 0.515, 0.0, 255.0);
768  BandInfo["LANDSAT7_ETM+_2_GREEN"] = SpectralSensorParams(2, 0.525, 0.605, 0.0, 255.0);
769  BandInfo["LANDSAT7_ETM+_3_RED"] = SpectralSensorParams(3, 0.63, 0.69, 0.0, 255.0);
770  BandInfo["LANDSAT7_ETM+_4_NIR"] = SpectralSensorParams(4, 0.75, 0.90, 0.0, 255.0);
771  BandInfo["LANDSAT7_ETM+_5_SWIR"] = SpectralSensorParams(5, 1.55, 1.75, 0.0, 255.0);
772  BandInfo["LANDSAT7_ETM+_6_THERMAL"] = SpectralSensorParams(6, 10.40, 12.5, 0.0, 255.0);
773  BandInfo["LANDSAT7_ETM+_7_MIR"] = SpectralSensorParams(7, 2.09, 2.35, 0.0, 255.0);
774  BandInfo["LANDSAT7_ETM+_8_PAN"] = SpectralSensorParams(8, 0.52, 0.90, 0.0, 255.0);
775 
776  BandInfo["WV2_MUL_1_COASTAL"] = SpectralSensorParams(1, 0.4, 0.45, 0.0, 2048.0);
777  BandInfo["WV2_MUL_2_BLUE"] = SpectralSensorParams(2, 0.45, 0.51, 0.0, 2048.0);
778  BandInfo["WV2_MUL_3_GREEN"] = SpectralSensorParams(3, 0.51, 0.58, 0.0, 2048.0);
779  BandInfo["WV2_MUL_4_YELLOW"] = SpectralSensorParams(4, 0.585, 0.625, 0.0, 2048.0);
780  BandInfo["WV2_MUL_5_RED"] = SpectralSensorParams(5, 0.66, 0.69, 0.0, 2048.0);
781  BandInfo["WV2_MUL_6_REDEDGE"] = SpectralSensorParams(6, 0.705, 0.745, 0.0, 2048.0);
782  BandInfo["WV2_MUL_7_NIR1"] = SpectralSensorParams(7, 0.77, 0.895, 0.0, 2048.0);
783  BandInfo["WV2_MUL_8_NIR2"] = SpectralSensorParams(8, 0.86, 0.104, 0.0, 2048.0);
784  BandInfo["WV2_PAN"] = SpectralSensorParams(9, 0.45, 0.8, 0.0, 2048.0);
785 
786  BandInfo["AWIFS_RED"] = SpectralSensorParams(1, 0.62, 0.68, 0.0, 4096.0);
787  BandInfo["AWIFS_GREEN"] = SpectralSensorParams(2, 0.52, 0.59, 0.0, 4096.0);
788  BandInfo["AWIFS_NIR"] = SpectralSensorParams(3, 0.77, 0.86, 0.0, 4096.0);
789  BandInfo["AWIFS_MIR"] = SpectralSensorParams(4, 1.55, 1.70, 0.0, 4096.0);
790 
791  BandInfo["RAPID_EYE_BLUE"] = SpectralSensorParams(1, 440.00, 510.00, 0.0, 4096.0);
792  BandInfo["RAPID_EYE_GREEN"] = SpectralSensorParams(2, 520.00, 590.00, 0.0, 4096.0);
793  BandInfo["RAPID_EYE_RED"] = SpectralSensorParams(3, 630.00, 685.00, 0.0, 4096.0);
794  BandInfo["RAPID_EYE_REDEDGE"] = SpectralSensorParams(4, 690.00, 730.00, 0.0, 4096.0);
795  BandInfo["RAPID_EYE_NIR"] = SpectralSensorParams(5, 760.00, 850.00, 0.0, 4096.0);
796 
797  BandInfo["LISS3_GREEN"] = SpectralSensorParams(1, 0.52, 0.59, 0.0, 128.0);
798  BandInfo["LISS3_RED"] = SpectralSensorParams(2, 0.62, 0.68, 0.0, 128.0);
799  BandInfo["LISS3_NEAR_INFRARED"] = SpectralSensorParams(3, 0.77, 0.86, 0.0, 128.0);
800  BandInfo["LISS3_SHORTWAVE_INFRARED"] = SpectralSensorParams(4, 1.55, 1.70, 0.0, 128.0);
801 
802  BandInfo["DMC_NIR"] = SpectralSensorParams(1, 0.77, 0.90, 0.0, 1024.0);
803  BandInfo["DMC_RED"] = SpectralSensorParams(2, 0.63, 0.69, 0.0, 1024.0);
804  BandInfo["DMC_GREEN"] = SpectralSensorParams(3, 0.52, 0.60, 0.0, 1024.0);
805 
806  BandInfo["LANDSAT8-OLI1_COASTAL_AEROSOL"] = SpectralSensorParams(1, 0.43, 0.45, 0.0, 65536.0);
807  BandInfo["LANDSAT8-OLI2_BLUE"] = SpectralSensorParams(2, 0.45, 0.51, 0.0, 65536.0);
808  BandInfo["LANDSAT8-OLI3_GREEN"] = SpectralSensorParams(3, 0.53, 0.59, 0.0, 65536.0);
809  BandInfo["LANDSAT8-OLI4_RED"] = SpectralSensorParams(4, 0.64, 0.67, 0.0, 65536.0);
810  BandInfo["LANDSAT8-OLI5_NIR"] = SpectralSensorParams(5, 0.85, 0.88, 0.0, 65536.0);
811  BandInfo["LANDSAT8-OLI6_SWIR1"] = SpectralSensorParams(6, 1.57, 1.65, 0.0, 65536.0);
812  BandInfo["LANDSAT8-OLI7_SWIR2"] = SpectralSensorParams(7, 2.11, 2.29, 0.0, 65536.0);
813  BandInfo["LANDSAT8-OLI8_PAN"] = SpectralSensorParams(8, 0.50, 0.68, 0.0, 65536.0);
814  BandInfo["LANDSAT8-OLI9_CIRRUS"] = SpectralSensorParams(9, 1.36, 1.38, 0.0, 65536.0);
815 
816  BandInfo["LANDSAT8-TIRS1"] = SpectralSensorParams(10, 10.60, 11.19, 0.0, 65536.0);
817  BandInfo["LANDSAT8-TIRS2"] = SpectralSensorParams(11, 11.50, 12.51, 0.0, 65536.0);
818 
819  BandInfo["CBERS4_PAN_PAN"] = SpectralSensorParams(1, 0.51, 0.85, 0.0, 255.0);
820  BandInfo["CBERS4_PAN_GREEN"] = SpectralSensorParams(2, 0.52, 0.59, 0.0, 255.0);
821  BandInfo["CBERS4_PAN_RED"] = SpectralSensorParams(3, 0.63, 0.69, 0.0, 255.0);
822  BandInfo["CBERS4_PAN_NIR"] = SpectralSensorParams(4, 0.77, 0.89, 0.0, 255.0);
823 
824  BandInfo["CBERS4_MUX_BLUE"] = SpectralSensorParams(5, 0.45, 0.52, 0.0, 255.0);
825  BandInfo["CBERS4_MUX_GREEN"] = SpectralSensorParams(6, 0.52, 0.59, 0.0, 255.0);
826  BandInfo["CBERS4_MUX_RED"] = SpectralSensorParams(7, 0.63, 0.69, 0.0, 255.0);
827  BandInfo["CBERS4_MUX_NIR"] = SpectralSensorParams(8, 0.77, 0.89, 0.0, 255.0);
828 
829  BandInfo["CBERS4_IRS_PAN"] = SpectralSensorParams(9, 0.50, 0.90, 0.0, 255.0);
830  BandInfo["CBERS4_IRS_SWIR_B10"] = SpectralSensorParams(10, 1.55, 1.75, 0.0, 255.0);
831  BandInfo["CBERS4_IRS_SWIR_B11"] = SpectralSensorParams(11, 2.08, 2.35, 0.0, 255.0);
832  BandInfo["CBERS4_IRS_TH"] = SpectralSensorParams(12, 10.40, 12.50, 0.0, 255.0);
833 
834  BandInfo["CBERS4_WFI_BLUE"] = SpectralSensorParams(13, 0.45, 0.52, 0.0, 1024.0);
835  BandInfo["CBERS4_WFI_GREEN"] = SpectralSensorParams(14, 0.52, 0.59, 0.0, 1024.0);
836  BandInfo["CBERS4_WFI_RED"] = SpectralSensorParams(15, 0.63, 0.69, 0.0, 1024.0);
837  BandInfo["CBERS4_WFI_NIR"] = SpectralSensorParams(16, 0.77, 0.89, 0.0, 1024.0);
838 
839 
840  // needs more Band Info from other sensors/bands
841  if (BandInfo.find(bandName) == BandInfo.end())
842  return SpectralSensorParams(0, 0.0, 1.0, 0.0, 255.0);
843 
844  return BandInfo[bandName];
845  }
846 
847  std::pair<double, double> GetDigitalNumberBandInfo(std::string bandName)
848  {
849  std::map<std::string, SpectralSensorParams > params = getSensorParams();
850  std::map<std::string, SpectralSensorParams > ::iterator it = params.find(bandName);
851  if (it == params.end())
852  return std::pair<double, double>(0.0, 255.0);
853 
854  return std::pair<double, double>(it->second.m_min, it->second.m_max);
855  }
856 
857  double GetDigitalNumberBandMax(std::string bandName)
858  {
859  return GetDigitalNumberBandInfo(bandName).second;
860  }
861 
862  double GetDigitalNumberBandMin(std::string bandName)
863  {
864  return GetDigitalNumberBandInfo(bandName).first;
865  }
866 
867  bool NormalizeRaster(te::rst::Raster& inraster, double nmin, double nmax)
868  {
869  if (nmin > nmax)
870  return false;
871 
872  // computing old maximuns and minimums for each band
873  std::vector<double> omins;
874  std::vector<double> omaxs;
875  std::vector<unsigned int> bands;
876  std::vector<double> dummy;
877 
878  for (unsigned int b = 0; b < inraster.getNumberOfBands(); b++)
879  {
880  const te::rst::RasterSummary* rsMin = te::rst::RasterSummaryManager::getInstance().get(&inraster, te::rst::SUMMARY_MIN, true);
881  const te::rst::RasterSummary* rsMax = te::rst::RasterSummaryManager::getInstance().get(&inraster, te::rst::SUMMARY_MAX, true);
882  const std::complex<double>* cmin = rsMin->at(b).m_minVal;
883  const std::complex<double>* cmax = rsMax->at(b).m_maxVal;
884 
885  omins.push_back(cmin->real());
886  omaxs.push_back(cmax->real());
887  dummy.push_back(inraster.getBand(b)->getProperty()->m_noDataValue);
888 
889  bands.push_back(b);
890  }
891 
892  // calculating amplitudes to avoid multiple subtractions
893  double value;
894  const double namplitude = nmax - nmin;
895  std::vector<double> oamplitude;
896  for (unsigned int b = 0; b < inraster.getNumberOfBands(); b++)
897  oamplitude.push_back(omaxs[b] - omins[b]);
898 
899  // iterating over raster to normalize pixel values
900  te::rst::RasterIterator<double> it = te::rst::RasterIterator<double>::begin(&inraster, bands);
901  te::rst::RasterIterator<double> itend = te::rst::RasterIterator<double>::end(&inraster, bands);
902 
903  while (it != itend)
904  {
905  for (unsigned int b = 0; b < inraster.getNumberOfBands(); b++)
906  {
907  value = (*it)[b];
908  if (value != dummy[b])
909  value = ((*it)[b] - omins[b]) * namplitude / oamplitude[b] + nmin;
910 
911  inraster.setValue(it.getColumn(), it.getRow(), value, b);
912  }
913 
914  ++it;
915  }
916 
917  return true;
918  }
919 
920  te::rst::Raster* NormalizeRaster(const te::rst::Raster* inputRaster)
921  {
922  double nmin = 0.0;
923  double nmax = 255.0;
924 
925  // computing old maximuns and minimums for each band
926  std::vector<double> omins;
927  std::vector<double> omaxs;
928  std::vector<unsigned int> bands;
929  std::vector<double> dummy;
930 
931  for (unsigned int b = 0; b < inputRaster->getNumberOfBands(); b++)
932  {
933  const te::rst::RasterSummary* rsMin = te::rst::RasterSummaryManager::getInstance().get(inputRaster, te::rst::SUMMARY_MIN, true);
934  const te::rst::RasterSummary* rsMax = te::rst::RasterSummaryManager::getInstance().get(inputRaster, te::rst::SUMMARY_MAX, true);
935  const std::complex<double>* cmin = rsMin->at(b).m_minVal;
936  const std::complex<double>* cmax = rsMax->at(b).m_maxVal;
937 
938  omins.push_back(cmin->real());
939  omaxs.push_back(cmax->real());
940  dummy.push_back(inputRaster->getBand(b)->getProperty()->m_noDataValue);
941 
942  bands.push_back(b);
943  }
944 
945  // calculating amplitudes to avoid multiple subtractions
946  double value;
947  const double namplitude = nmax - nmin;
948  std::vector<double> oamplitude;
949  for (unsigned int b = 0; b < inputRaster->getNumberOfBands(); b++)
950  oamplitude.push_back(omaxs[b] - omins[b]);
951 
952  // iterating over raster to normalize pixel values
953  te::rst::RasterIterator<double> it = te::rst::RasterIterator<double>::begin((te::rst::Raster*)inputRaster, bands);
954  te::rst::RasterIterator<double> itend = te::rst::RasterIterator<double>::end((te::rst::Raster*)inputRaster, bands);
955 
956  std::vector<te::rst::BandProperty*> bandsProperties;
957 
958  for (unsigned i = 0; i < inputRaster->getNumberOfBands(); i++)
959  {
960  te::rst::BandProperty* bprop = new te::rst::BandProperty(*(inputRaster->getBand(i)->getProperty()));
961  bprop->m_type = te::dt::UCHAR_TYPE;
962  bprop->m_noDataValue = inputRaster->getBand(i)->getProperty()->m_noDataValue;
963  bandsProperties.push_back(bprop);
964  }
965 
966  std::map< std::string, std::string > rInfo;
967 
968  te::rst::Raster* outputRaster = te::rst::RasterFactory::make("MEM", new te::rst::Grid(*inputRaster->getGrid()), bandsProperties, rInfo);
969 
970  while (it != itend)
971  {
972  for (unsigned int b = 0; b < inputRaster->getNumberOfBands(); b++)
973  {
974  value = (*it)[b];
975  if (value != dummy[b])
976  value = ((*it)[b] - omins[b]) * namplitude / oamplitude[b] + nmin;
977 
978  outputRaster->setValue(it.getColumn(), it.getRow(), value, b);
979  }
980 
981  ++it;
982  }
983 
984  double rasterValue;
985 
986  for (unsigned r = 0; r < inputRaster->getNumberOfRows(); r++)
987  {
988  for (unsigned c = 0; c < inputRaster->getNumberOfColumns(); c++)
989  {
990  for (unsigned b = 0; b < inputRaster->getNumberOfBands(); b++)
991  {
992  inputRaster->getValue(c, r, rasterValue, b);
993 
994  if (rasterValue != inputRaster->getBand(b)->getProperty()->m_noDataValue)
995  rasterValue = 255 * (rasterValue - omins[b]) / (omaxs[b] - omins[b]);
996 
997  outputRaster->setValue(c, r, rasterValue, b);
998  }
999  }
1000  }
1001 
1002  return outputRaster;
1003  }
1004 
1005  bool ConvertRGB2IHS( const te::rst::Raster& inputRGBRaster,
1006  const unsigned int redBandIdx, const unsigned int greenBandIdx,
1007  const unsigned int blueBandIdx, const double rgbRangeMin,
1008  const double rgbRangeMax, te::rst::Raster& outputIHSRaster )
1009  {
1010  if( ( inputRGBRaster.getAccessPolicy() & te::common::RAccess ) == 0 )
1011  return false;
1012  if( redBandIdx > inputRGBRaster.getNumberOfBands() ) return false;
1013  if( greenBandIdx > inputRGBRaster.getNumberOfBands() ) return false;
1014  if( blueBandIdx > inputRGBRaster.getNumberOfBands() ) return false;
1015  if( ( outputIHSRaster.getAccessPolicy() & te::common::WAccess ) == 0 )
1016  return false;
1017  if( outputIHSRaster.getNumberOfBands() < 3 ) return false;
1018  if( inputRGBRaster.getNumberOfColumns() != outputIHSRaster.getNumberOfColumns() )
1019  return false;
1020  if( inputRGBRaster.getNumberOfRows() != outputIHSRaster.getNumberOfRows() )
1021  return false;
1022 
1023  const te::rst::Band& redBand = *inputRGBRaster.getBand( redBandIdx );
1024  const te::rst::Band& greenBand = *inputRGBRaster.getBand( greenBandIdx );
1025  const te::rst::Band& blueBand = *inputRGBRaster.getBand( blueBandIdx );
1026 
1027  te::rst::Band& iBand = *outputIHSRaster.getBand( 0 );
1028  te::rst::Band& hBand = *outputIHSRaster.getBand( 1 );
1029  te::rst::Band& sBand = *outputIHSRaster.getBand( 2 );
1030 
1031  const unsigned int outNCols = outputIHSRaster.getNumberOfColumns();
1032  const unsigned int outNRows = outputIHSRaster.getNumberOfRows();
1033  const double redNoData = redBand.getProperty()->m_noDataValue;
1034  const double greenNoData = greenBand.getProperty()->m_noDataValue;
1035  const double blueNoData = inputRGBRaster.getBand(
1036  blueBandIdx )->getProperty()->m_noDataValue;
1037  const double iNoData = iBand.getProperty()->m_noDataValue;
1038  const double hNoData = hBand.getProperty()->m_noDataValue;
1039  const double sNoData = sBand.getProperty()->m_noDataValue;
1040 
1041  unsigned int outRow = 0;
1042  unsigned int outCol = 0;
1043  double red = 0;
1044  double green = 0;
1045  double blue = 0;
1046  double teta = 0;
1047  double redNorm = 0, greenNorm = 0, blueNorm = 0;
1048  double rMinusG = 0, rMinusB = 0;
1049  double rgbSum = 0;
1050  double cosValue = 0;
1051  const double twoPi = 2.0 * (static_cast<double>(M_PI));
1052  const double pi2 = (static_cast<double>(M_PI)) / 2.0;
1053  const double rgbNormFac = ( rgbRangeMax == rgbRangeMin ) ? 0.0 :
1054  ( 1.0 / ( rgbRangeMax - rgbRangeMin ) );
1055  double intensity = 0;
1056  double hue = 0;
1057  double saturation = 0;
1058 
1059  for( outRow = 0 ; outRow < outNRows ; ++outRow )
1060  {
1061  for( outCol = 0 ; outCol < outNCols ; ++outCol )
1062  {
1063  redBand.getValue( outCol, outRow, red );
1064  greenBand.getValue( outCol, outRow, green );
1065  blueBand.getValue( outCol, outRow, blue );
1066 
1067  if( ( red == redNoData ) || ( green == greenNoData ) ||
1068  ( blue == blueNoData ) )
1069  {
1070  iBand.setValue( outCol, outRow, iNoData );
1071  hBand.setValue( outCol, outRow, hNoData );
1072  sBand.setValue( outCol, outRow, sNoData );
1073  }
1074  else
1075  {
1076  if( ( red == green ) && ( green == blue ) )
1077  { // Gray scale case
1078  // From Wikipedia:
1079  // h = 0 is used for grays though the hue has no geometric
1080  // meaning there, where the saturation s = 0. Similarly,
1081  // the choice of 0 as the value for s when l is equal to 0 or 1
1082  // is arbitrary.
1083 
1084  hue = 0.0;
1085  saturation = 0.0;
1086  intensity = ( red * rgbNormFac ); // or green or blue since they all are the same.
1087  }
1088  else
1089  { // Color case
1090  redNorm = ( red - rgbRangeMin ) * rgbNormFac;
1091  greenNorm = ( green - rgbRangeMin ) * rgbNormFac;
1092  blueNorm = ( blue - rgbRangeMin ) * rgbNormFac;
1093 
1094  rMinusG = redNorm - greenNorm;
1095  rMinusB = redNorm - blueNorm;
1096 
1097  cosValue = std::sqrt( ( rMinusG * rMinusG ) + ( rMinusB *
1098  ( greenNorm - blueNorm ) ) );
1099 
1100  if( cosValue == 0.0 )
1101  {
1102  teta = pi2;
1103  }
1104  else
1105  {
1106  cosValue = ( 0.5 * ( rMinusG + rMinusB ) ) /
1107  cosValue;
1108  teta = std::acos( cosValue );
1109  }
1110 
1111  assert( ( cosValue >= (-1.0) ) && ( cosValue <= (1.0) ) );
1112 
1113  if( blueNorm > greenNorm )
1114  {
1115  hue = ( twoPi - teta );
1116  }
1117  else
1118  {
1119  hue = teta;
1120  }
1121 
1122  rgbSum = ( redNorm + greenNorm + blueNorm );
1123 
1124  saturation = ( 1.0 - ( 3 * std::min( std::min( redNorm, greenNorm ), blueNorm ) /
1125  rgbSum ) );
1126 
1127  intensity = ( rgbSum / 3.0 );
1128  }
1129 
1130  iBand.setValue( outCol, outRow, intensity );
1131  hBand.setValue( outCol, outRow, hue );
1132  sBand.setValue( outCol, outRow, saturation );
1133  }
1134  }
1135  }
1136 
1137  return true;
1138  }
1139 
1140  bool ConvertRGB2IHS(const te::rst::Raster& inputRedRaster, const unsigned int redBandIdx,
1141  const te::rst::Raster& inputGreenRaster, const unsigned int greenBandIdx,
1142  const te::rst::Raster& inputBlueRaster, const unsigned int blueBandIdx,
1143  const double rgbRangeMin, const double rgbRangeMax, te::rst::Raster& outputIHSRaster)
1144  {
1145  if ((inputRedRaster.getAccessPolicy() & te::common::RAccess) == 0 ||
1146  (inputGreenRaster.getAccessPolicy() & te::common::RAccess) == 0 ||
1147  (inputBlueRaster.getAccessPolicy() & te::common::RAccess) == 0)
1148  return false;
1149  if (redBandIdx > inputRedRaster.getNumberOfBands()) return false;
1150  if (greenBandIdx > inputGreenRaster.getNumberOfBands()) return false;
1151  if (blueBandIdx > inputBlueRaster.getNumberOfBands()) return false;
1152  if ((outputIHSRaster.getAccessPolicy() & te::common::WAccess) == 0)
1153  return false;
1154  if (outputIHSRaster.getNumberOfBands() < 3) return false;
1155  if (inputRedRaster.getNumberOfColumns() != inputGreenRaster.getNumberOfColumns() ||
1156  inputRedRaster.getNumberOfColumns() != inputBlueRaster.getNumberOfColumns())
1157  return false;
1158  if (inputRedRaster.getNumberOfRows() != inputGreenRaster.getNumberOfRows() ||
1159  inputRedRaster.getNumberOfRows() != inputBlueRaster.getNumberOfRows())
1160  return false;
1161  if (inputRedRaster.getNumberOfColumns() != outputIHSRaster.getNumberOfColumns())
1162  return false;
1163  if (inputRedRaster.getNumberOfRows() != outputIHSRaster.getNumberOfRows())
1164  return false;
1165 
1166  const te::rst::Band& redBand = *inputRedRaster.getBand(redBandIdx);
1167  const te::rst::Band& greenBand = *inputGreenRaster.getBand(greenBandIdx);
1168  const te::rst::Band& blueBand = *inputBlueRaster.getBand(blueBandIdx);
1169 
1170  te::rst::Band& iBand = *outputIHSRaster.getBand(0);
1171  te::rst::Band& hBand = *outputIHSRaster.getBand(1);
1172  te::rst::Band& sBand = *outputIHSRaster.getBand(2);
1173 
1174  const unsigned int outNCols = outputIHSRaster.getNumberOfColumns();
1175  const unsigned int outNRows = outputIHSRaster.getNumberOfRows();
1176  const double redNoData = redBand.getProperty()->m_noDataValue;
1177  const double greenNoData = greenBand.getProperty()->m_noDataValue;
1178  const double blueNoData = blueBand.getProperty()->m_noDataValue;
1179  const double iNoData = iBand.getProperty()->m_noDataValue;
1180  const double hNoData = hBand.getProperty()->m_noDataValue;
1181  const double sNoData = sBand.getProperty()->m_noDataValue;
1182 
1183  unsigned int outRow = 0;
1184  unsigned int outCol = 0;
1185  double red = 0;
1186  double green = 0;
1187  double blue = 0;
1188  double teta = 0;
1189  double redNorm = 0, greenNorm = 0, blueNorm = 0;
1190  double rMinusG = 0, rMinusB = 0;
1191  double rgbSum = 0;
1192  double cosValue = 0;
1193  const double twoPi = 2.0 * (static_cast<double>(M_PI));
1194  const double pi2 = (static_cast<double>(M_PI)) / 2.0;
1195  const double rgbNormFac = (rgbRangeMax == rgbRangeMin) ? 0.0 :
1196  (1.0 / (rgbRangeMax - rgbRangeMin));
1197  double intensity = 0;
1198  double hue = 0;
1199  double saturation = 0;
1200 
1201  for (outRow = 0; outRow < outNRows; ++outRow)
1202  {
1203  for (outCol = 0; outCol < outNCols; ++outCol)
1204  {
1205  redBand.getValue(outCol, outRow, red);
1206  greenBand.getValue(outCol, outRow, green);
1207  blueBand.getValue(outCol, outRow, blue);
1208 
1209  if ((red == redNoData) || (green == greenNoData) ||
1210  (blue == blueNoData))
1211  {
1212  iBand.setValue(outCol, outRow, iNoData);
1213  hBand.setValue(outCol, outRow, hNoData);
1214  sBand.setValue(outCol, outRow, sNoData);
1215  }
1216  else
1217  {
1218  if ((red == green) && (green == blue))
1219  { // Gray scale case
1220  // From Wikipedia:
1221  // h = 0 is used for grays though the hue has no geometric
1222  // meaning there, where the saturation s = 0. Similarly,
1223  // the choice of 0 as the value for s when l is equal to 0 or 1
1224  // is arbitrary.
1225 
1226  hue = 0.0;
1227  saturation = 0.0;
1228  intensity = (red * rgbNormFac); // or green or blue since they all are the same.
1229  }
1230  else
1231  { // Color case
1232  redNorm = (red - rgbRangeMin) * rgbNormFac;
1233  greenNorm = (green - rgbRangeMin) * rgbNormFac;
1234  blueNorm = (blue - rgbRangeMin) * rgbNormFac;
1235 
1236  rMinusG = redNorm - greenNorm;
1237  rMinusB = redNorm - blueNorm;
1238 
1239  cosValue = std::sqrt((rMinusG * rMinusG) + (rMinusB *
1240  (greenNorm - blueNorm)));
1241 
1242  if (cosValue == 0.0)
1243  {
1244  teta = pi2;
1245  }
1246  else
1247  {
1248  cosValue = (0.5 * (rMinusG + rMinusB)) /
1249  cosValue;
1250  teta = std::acos(cosValue);
1251  }
1252 
1253  assert((cosValue >= (-1.0)) && (cosValue <= (1.0)));
1254 
1255  if (blueNorm > greenNorm)
1256  {
1257  hue = (twoPi - teta);
1258  }
1259  else
1260  {
1261  hue = teta;
1262  }
1263 
1264  rgbSum = (redNorm + greenNorm + blueNorm);
1265 
1266  saturation = (1.0 - (3 * std::min(std::min(redNorm, greenNorm), blueNorm) /
1267  rgbSum));
1268 
1269  intensity = (rgbSum / 3.0);
1270  }
1271 
1272  iBand.setValue(outCol, outRow, intensity);
1273  hBand.setValue(outCol, outRow, hue);
1274  sBand.setValue(outCol, outRow, saturation);
1275  }
1276  }
1277  }
1278 
1279  return true;
1280  }
1281 
1282  bool ConvertIHS2RGB( const te::rst::Raster& inputIHSRaster,
1283  const unsigned int intensityBandIdx, const unsigned int hueBandIdx,
1284  const unsigned int saturationBandIdx, const double rgbRangeMin,
1285  const double rgbRangeMax, te::rst::Raster& outputRGBRaster )
1286  {
1287  if( ( inputIHSRaster.getAccessPolicy() & te::common::RAccess ) == 0 )
1288  return false;
1289  if( intensityBandIdx > inputIHSRaster.getNumberOfBands() ) return false;
1290  if( hueBandIdx > inputIHSRaster.getNumberOfBands() ) return false;
1291  if( saturationBandIdx > inputIHSRaster.getNumberOfBands() ) return false;
1292  if( ( outputRGBRaster.getAccessPolicy() & te::common::WAccess ) == 0 )
1293  return false;
1294  if( outputRGBRaster.getNumberOfBands() < 3 ) return false;
1295  if( inputIHSRaster.getNumberOfColumns() != outputRGBRaster.getNumberOfColumns() )
1296  return false;
1297  if( inputIHSRaster.getNumberOfRows() != outputRGBRaster.getNumberOfRows() )
1298  return false;
1299 
1300  const te::rst::Band& intensityBand = *inputIHSRaster.getBand( intensityBandIdx );
1301  const te::rst::Band& hueBand = *inputIHSRaster.getBand( hueBandIdx );
1302  const te::rst::Band& saturationBand = *inputIHSRaster.getBand( saturationBandIdx );
1303  te::rst::Band& redBand = *outputRGBRaster.getBand( 0 );
1304  te::rst::Band& greenBand = *outputRGBRaster.getBand( 1 );
1305  te::rst::Band& blueBand = *outputRGBRaster.getBand( 2 );
1306 
1307  const unsigned int nCols = outputRGBRaster.getNumberOfColumns();
1308  const unsigned int nRows = outputRGBRaster.getNumberOfRows();
1309 
1310  const double intensityNoData = intensityBand.getProperty()->m_noDataValue;
1311  const double hueNoData = hueBand.getProperty()->m_noDataValue;
1312  const double saturationNoData = saturationBand.getProperty()->m_noDataValue;
1313 
1314  const double rNoData = redBand.getProperty()->m_noDataValue;
1315  const double gNoData = greenBand.getProperty()->m_noDataValue;
1316  const double bNoData = blueBand.getProperty()->m_noDataValue;
1317 
1318  const double rgbNormFac = ( rgbRangeMax == rgbRangeMin ) ? 0.0 :
1319  ( rgbRangeMax - rgbRangeMin );
1320  const double pi3 = M_PI / 3.0; // 60
1321  const double twoPi3 = 2.0 * M_PI / 3.0; // 120
1322  const double fourPi3 = 4.0 * M_PI / 3.0; // 240
1323  unsigned int row = 0;
1324  unsigned int col = 0;
1325  double hue = 0;
1326  double lig = 0;
1327  double sat = 0;
1328  double red = 0;
1329  double green = 0;
1330  double blue = 0;
1331 
1332 
1333  for( row = 0 ; row < nRows ; ++row )
1334  {
1335  for( col = 0 ; col < nCols ; ++col )
1336  {
1337  intensityBand.getValue( col, row, lig );
1338  hueBand.getValue( col, row, hue );
1339  saturationBand.getValue( col, row, sat );
1340 
1341  if( ( lig == intensityNoData ) || ( hue == hueNoData ) ||
1342  ( sat == saturationNoData ) )
1343  {
1344  redBand.setValue( col, row, rNoData );
1345  greenBand.setValue( col, row, gNoData );
1346  blueBand.setValue( col, row, bNoData );
1347  }
1348  else
1349  {
1350  if( ( hue == 0.0 ) && ( sat == 0.0 ) )
1351  { // Gray scale case
1352  red = green = blue = ( lig * rgbNormFac );
1353  }
1354  else
1355  { // color case
1356  /* Hue inside RG sector */
1357  if( hue < twoPi3 )
1358  {
1359  blue = lig * ( 1.0 - sat );
1360  red = lig * ( 1.0 + ( sat * std::cos( hue ) /
1361  std::cos( pi3 - hue ) ) );
1362  green = ( 3.0 * lig ) - ( red + blue );
1363  }
1364  else if( hue < fourPi3 )
1365  { /* Hue inside GB sector */
1366 
1367  hue -= twoPi3;
1368 
1369  red = lig * ( 1.0 - sat );
1370  green = lig * ( 1.0 + ( sat * std::cos( hue ) /
1371  std::cos( pi3 - hue ) ) );
1372  blue = ( 3.0 * lig ) - ( red + green );
1373  }
1374  else
1375  { /* Hue inside BR sector */
1376 
1377  hue -= fourPi3;
1378 
1379  green = lig * ( 1.0 - sat );
1380  blue = lig * ( 1.0 + ( sat * std::cos( hue ) /
1381  std::cos( pi3 - hue ) ) );
1382  red = ( 3.0 * lig ) - ( green + blue );
1383  }
1384 
1385  red = ( red * rgbNormFac ) + rgbRangeMin;
1386  green = ( green * rgbNormFac ) + rgbRangeMin;
1387  blue = ( blue * rgbNormFac ) + rgbRangeMin;
1388  }
1389 
1390  red = MIN( red, rgbRangeMax );
1391  green = MIN( green, rgbRangeMax );
1392  blue = MIN( blue, rgbRangeMax );
1393 
1394  red = MAX( red, rgbRangeMin );
1395  green = MAX( green, rgbRangeMin );
1396  blue = MAX( blue, rgbRangeMin );
1397 
1398  redBand.setValue( col, row, red );
1399  greenBand.setValue( col, row, green );
1400  blueBand.setValue( col, row, blue );
1401  }
1402  }
1403  }
1404 
1405  return true;
1406  }
1407 
1408  bool ConvertIHS2RGB(const te::rst::Raster& inputIRaster, const unsigned int intensityBandIdx,
1409  const te::rst::Raster& inputHRaster, const unsigned int hueBandIdx,
1410  const te::rst::Raster& inputSRaster, const unsigned int saturationBandIdx,
1411  const double rgbRangeMin, const double rgbRangeMax, te::rst::Raster& outputRGBRaster)
1412  {
1413  if ((inputIRaster.getAccessPolicy() & te::common::RAccess) == 0 ||
1414  (inputHRaster.getAccessPolicy() & te::common::RAccess) == 0 ||
1415  (inputSRaster.getAccessPolicy() & te::common::RAccess) == 0)
1416  return false;
1417  if (intensityBandIdx > inputIRaster.getNumberOfBands()) return false;
1418  if (hueBandIdx > inputHRaster.getNumberOfBands()) return false;
1419  if (saturationBandIdx > inputSRaster.getNumberOfBands()) return false;
1420  if (inputIRaster.getNumberOfColumns() != inputHRaster.getNumberOfColumns() ||
1421  inputIRaster.getNumberOfColumns() != inputSRaster.getNumberOfColumns())
1422  return false;
1423  if (inputIRaster.getNumberOfRows() != inputHRaster.getNumberOfRows() ||
1424  inputIRaster.getNumberOfRows() != inputSRaster.getNumberOfRows())
1425  return false;
1426  if ((outputRGBRaster.getAccessPolicy() & te::common::WAccess) == 0)
1427  return false;
1428  if (outputRGBRaster.getNumberOfBands() < 3) return false;
1429  if (inputIRaster.getNumberOfColumns() != outputRGBRaster.getNumberOfColumns())
1430  return false;
1431  if (inputIRaster.getNumberOfRows() != outputRGBRaster.getNumberOfRows())
1432  return false;
1433 
1434  const te::rst::Band& intensityBand = *inputIRaster.getBand(intensityBandIdx);
1435  const te::rst::Band& hueBand = *inputHRaster.getBand(hueBandIdx);
1436  const te::rst::Band& saturationBand = *inputSRaster.getBand(saturationBandIdx);
1437  te::rst::Band& redBand = *outputRGBRaster.getBand(0);
1438  te::rst::Band& greenBand = *outputRGBRaster.getBand(1);
1439  te::rst::Band& blueBand = *outputRGBRaster.getBand(2);
1440 
1441  const unsigned int nCols = outputRGBRaster.getNumberOfColumns();
1442  const unsigned int nRows = outputRGBRaster.getNumberOfRows();
1443 
1444  const double intensityNoData = intensityBand.getProperty()->m_noDataValue;
1445  const double hueNoData = hueBand.getProperty()->m_noDataValue;
1446  const double saturationNoData = saturationBand.getProperty()->m_noDataValue;
1447 
1448  const double rNoData = redBand.getProperty()->m_noDataValue;
1449  const double gNoData = greenBand.getProperty()->m_noDataValue;
1450  const double bNoData = blueBand.getProperty()->m_noDataValue;
1451 
1452  const double rgbNormFac = (rgbRangeMax == rgbRangeMin) ? 0.0 :
1453  (rgbRangeMax - rgbRangeMin);
1454  const double pi3 = M_PI / 3.0; // 60
1455  const double twoPi3 = 2.0 * M_PI / 3.0; // 120
1456  const double fourPi3 = 4.0 * M_PI / 3.0; // 240
1457  unsigned int row = 0;
1458  unsigned int col = 0;
1459  double hue = 0;
1460  double lig = 0;
1461  double sat = 0;
1462  double red = 0;
1463  double green = 0;
1464  double blue = 0;
1465 
1466 
1467  for (row = 0; row < nRows; ++row)
1468  {
1469  for (col = 0; col < nCols; ++col)
1470  {
1471  intensityBand.getValue(col, row, lig);
1472  hueBand.getValue(col, row, hue);
1473  saturationBand.getValue(col, row, sat);
1474 
1475  if ((lig == intensityNoData) || (hue == hueNoData) ||
1476  (sat == saturationNoData))
1477  {
1478  redBand.setValue(col, row, rNoData);
1479  greenBand.setValue(col, row, gNoData);
1480  blueBand.setValue(col, row, bNoData);
1481  }
1482  else
1483  {
1484  if ((hue == 0.0) && (sat == 0.0))
1485  { // Gray scale case
1486  red = green = blue = (lig * rgbNormFac);
1487  }
1488  else
1489  { // color case
1490  /* Hue inside RG sector */
1491  if (hue < twoPi3)
1492  {
1493  blue = lig * (1.0 - sat);
1494  red = lig * (1.0 + (sat * std::cos(hue) /
1495  std::cos(pi3 - hue)));
1496  green = (3.0 * lig) - (red + blue);
1497  }
1498  else if (hue < fourPi3)
1499  { /* Hue inside GB sector */
1500 
1501  hue -= twoPi3;
1502 
1503  red = lig * (1.0 - sat);
1504  green = lig * (1.0 + (sat * std::cos(hue) /
1505  std::cos(pi3 - hue)));
1506  blue = (3.0 * lig) - (red + green);
1507  }
1508  else
1509  { /* Hue inside BR sector */
1510 
1511  hue -= fourPi3;
1512 
1513  green = lig * (1.0 - sat);
1514  blue = lig * (1.0 + (sat * std::cos(hue) /
1515  std::cos(pi3 - hue)));
1516  red = (3.0 * lig) - (green + blue);
1517  }
1518 
1519  red = (red * rgbNormFac) + rgbRangeMin;
1520  green = (green * rgbNormFac) + rgbRangeMin;
1521  blue = (blue * rgbNormFac) + rgbRangeMin;
1522  }
1523 
1524  red = MIN(red, rgbRangeMax);
1525  green = MIN(green, rgbRangeMax);
1526  blue = MIN(blue, rgbRangeMax);
1527 
1528  red = MAX(red, rgbRangeMin);
1529  green = MAX(green, rgbRangeMin);
1530  blue = MAX(blue, rgbRangeMin);
1531 
1532  redBand.setValue(col, row, red);
1533  greenBand.setValue(col, row, green);
1534  blueBand.setValue(col, row, blue);
1535  }
1536  }
1537  }
1538 
1539  return true;
1540  }
1541 
1542  bool ConvertRGB2HLS(const te::rst::Raster& inputRGBRaster,
1543  const unsigned int redBandIdx, const unsigned int greenBandIdx,
1544  const unsigned int blueBandIdx, const double /*rgbRangeMin*/,
1545  const double rgbRangeMax, te::rst::Raster& outputHLSRaster)
1546  {
1547  if ((inputRGBRaster.getAccessPolicy() & te::common::RAccess) == 0)
1548  return false;
1549  if (redBandIdx > inputRGBRaster.getNumberOfBands()) return false;
1550  if (greenBandIdx > inputRGBRaster.getNumberOfBands()) return false;
1551  if (blueBandIdx > inputRGBRaster.getNumberOfBands()) return false;
1552  if ((outputHLSRaster.getAccessPolicy() & te::common::WAccess) == 0)
1553  return false;
1554  if (outputHLSRaster.getNumberOfBands() < 3) return false;
1555  if (inputRGBRaster.getNumberOfColumns() != outputHLSRaster.getNumberOfColumns())
1556  return false;
1557  if (inputRGBRaster.getNumberOfRows() != outputHLSRaster.getNumberOfRows())
1558  return false;
1559 
1560  const te::rst::Band& redBand = *inputRGBRaster.getBand(redBandIdx);
1561  const te::rst::Band& greenBand = *inputRGBRaster.getBand(greenBandIdx);
1562  const te::rst::Band& blueBand = *inputRGBRaster.getBand(blueBandIdx);
1563 
1564  te::rst::Band& hBand = *outputHLSRaster.getBand(0);
1565  te::rst::Band& lBand = *outputHLSRaster.getBand(1);
1566  te::rst::Band& sBand = *outputHLSRaster.getBand(2);
1567 
1568  const unsigned int outNCols = outputHLSRaster.getNumberOfColumns();
1569  const unsigned int outNRows = outputHLSRaster.getNumberOfRows();
1570 
1571  unsigned int outRow = 0;
1572  unsigned int outCol = 0;
1573  double red = 0;
1574  double green = 0;
1575  double blue = 0;
1576  double light = 0;
1577  double hue = 0;
1578  double saturation = 0;
1579 
1580  for (outRow = 0; outRow < outNRows; ++outRow)
1581  {
1582  for (outCol = 0; outCol < outNCols; ++outCol)
1583  {
1584  redBand.getValue(outCol, outRow, red);
1585  greenBand.getValue(outCol, outRow, green);
1586  blueBand.getValue(outCol, outRow, blue);
1587 
1588  red /= rgbRangeMax;
1589  green /= rgbRangeMax;
1590  blue /= rgbRangeMax;
1591 
1592  double vmax = MAX(red, green);
1593  vmax = MAX(vmax, blue);
1594 
1595  double vmin = MIN(red, green);
1596  vmin = MIN(vmin, blue);
1597 
1598  hue = (vmax + vmin) / 2.;
1599  saturation = (vmax + vmin) / 2.;
1600  light = (vmax + vmin) / 2.;
1601 
1602  if (vmax == vmin)
1603  {// achromatic
1604  hue = 0.;
1605  saturation = 0.;
1606  }
1607  else
1608  {
1609  double diff = vmax - vmin;
1610 
1611  saturation = light > 0.5 ? diff / (2. - vmax - vmin) : diff / (vmax + vmin);
1612 
1613  if (vmax == red)
1614  {
1615  hue = (green - blue) / diff + (green < blue ? 6. : 0.);
1616  }
1617  else if (vmax == green)
1618  {
1619  hue = (blue - red) / diff + 2.;
1620  }
1621  else if (vmax == blue)
1622  {
1623  hue = (red - green) / diff + 4.;
1624  }
1625 
1626  hue *= 60.;
1627  }
1628 
1629  hBand.setValue(outCol, outRow, hue);
1630  lBand.setValue(outCol, outRow, light);
1631  sBand.setValue(outCol, outRow, saturation);
1632  }
1633  }
1634 
1635  return true;
1636  }
1637 
1638  bool ConvertRGB2HLS(const te::rst::Raster& inputRedRaster, const unsigned int redBandIdx,
1639  const te::rst::Raster& inputGreenRaster, const unsigned int greenBandIdx,
1640  const te::rst::Raster& inputBlueRaster, const unsigned int blueBandIdx,
1641  const double /*rgbRangeMin*/, const double rgbRangeMax, te::rst::Raster& outputHLSRaster)
1642  {
1643  if ((inputRedRaster.getAccessPolicy() & te::common::RAccess) == 0 ||
1644  (inputGreenRaster.getAccessPolicy() & te::common::RAccess) == 0 ||
1645  (inputBlueRaster.getAccessPolicy() & te::common::RAccess) == 0)
1646  return false;
1647  if (redBandIdx > inputRedRaster.getNumberOfBands()) return false;
1648  if (greenBandIdx > inputGreenRaster.getNumberOfBands()) return false;
1649  if (blueBandIdx > inputBlueRaster.getNumberOfBands()) return false;
1650  if ((outputHLSRaster.getAccessPolicy() & te::common::WAccess) == 0)
1651  return false;
1652  if (outputHLSRaster.getNumberOfBands() < 3) return false;
1653  if (inputRedRaster.getNumberOfRows() != inputGreenRaster.getNumberOfRows() ||
1654  inputRedRaster.getNumberOfRows() != inputBlueRaster.getNumberOfRows())
1655  return false;
1656  if (inputRedRaster.getNumberOfColumns() != outputHLSRaster.getNumberOfColumns())
1657  return false;
1658  if (inputRedRaster.getNumberOfRows() != outputHLSRaster.getNumberOfRows())
1659  return false;
1660 
1661  const te::rst::Band& redBand = *inputRedRaster.getBand(redBandIdx);
1662  const te::rst::Band& greenBand = *inputGreenRaster.getBand(greenBandIdx);
1663  const te::rst::Band& blueBand = *inputBlueRaster.getBand(blueBandIdx);
1664 
1665  te::rst::Band& hBand = *outputHLSRaster.getBand(0);
1666  te::rst::Band& lBand = *outputHLSRaster.getBand(1);
1667  te::rst::Band& sBand = *outputHLSRaster.getBand(2);
1668 
1669  const unsigned int outNCols = outputHLSRaster.getNumberOfColumns();
1670  const unsigned int outNRows = outputHLSRaster.getNumberOfRows();
1671 
1672  unsigned int outRow = 0;
1673  unsigned int outCol = 0;
1674  double red = 0;
1675  double green = 0;
1676  double blue = 0;
1677  double light = 0;
1678  double hue = 0;
1679  double saturation = 0;
1680 
1681  for (outRow = 0; outRow < outNRows; ++outRow)
1682  {
1683  for (outCol = 0; outCol < outNCols; ++outCol)
1684  {
1685  redBand.getValue(outCol, outRow, red);
1686  greenBand.getValue(outCol, outRow, green);
1687  blueBand.getValue(outCol, outRow, blue);
1688 
1689  red /= rgbRangeMax;
1690  green /= rgbRangeMax;
1691  blue /= rgbRangeMax;
1692 
1693  double vmax = MAX(red, green);
1694  vmax = MAX(vmax, blue);
1695 
1696  double vmin = MIN(red, green);
1697  vmin = MIN(vmin, blue);
1698 
1699  hue = (vmax + vmin) / 2.;
1700  saturation = (vmax + vmin) / 2.;
1701  light = (vmax + vmin) / 2.;
1702 
1703  if (vmax == vmin)
1704  {// achromatic
1705  hue = 0.;
1706  saturation = 0.;
1707  }
1708  else
1709  {
1710  double diff = vmax - vmin;
1711 
1712  saturation = light > 0.5 ? diff / (2. - vmax - vmin) : diff / (vmax + vmin);
1713 
1714  if (vmax == red)
1715  {
1716  hue = (green - blue) / diff + (green < blue ? 6. : 0.);
1717  }
1718  else if (vmax == green)
1719  {
1720  hue = (blue - red) / diff + 2.;
1721  }
1722  else if (vmax == blue)
1723  {
1724  hue = (red - green) / diff + 4.;
1725  }
1726 
1727  hue *= 60.;
1728  }
1729 
1730  hBand.setValue(outCol, outRow, hue);
1731  lBand.setValue(outCol, outRow, light);
1732  sBand.setValue(outCol, outRow, saturation);
1733  }
1734  }
1735 
1736  return true;
1737  }
1738 
1739  bool ConvertHLS2RGB(const te::rst::Raster& inputHLSRaster,
1740  const unsigned int hueBandIdx, const unsigned int lightBandIdx,
1741  const unsigned int saturationBandIdx, const double /*rgbRangeMin*/,
1742  const double rgbRangeMax, te::rst::Raster& outputRGBRaster)
1743  {
1744  if ((inputHLSRaster.getAccessPolicy() & te::common::RAccess) == 0)
1745  return false;
1746  if (hueBandIdx > inputHLSRaster.getNumberOfBands()) return false;
1747  if (lightBandIdx > inputHLSRaster.getNumberOfBands()) return false;
1748  if (saturationBandIdx > inputHLSRaster.getNumberOfBands()) return false;
1749  if ((outputRGBRaster.getAccessPolicy() & te::common::WAccess) == 0)
1750  return false;
1751  if (outputRGBRaster.getNumberOfBands() < 3) return false;
1752  if (inputHLSRaster.getNumberOfColumns() != outputRGBRaster.getNumberOfColumns())
1753  return false;
1754  if (inputHLSRaster.getNumberOfRows() != outputRGBRaster.getNumberOfRows())
1755  return false;
1756 
1757  const te::rst::Band& hueBand = *inputHLSRaster.getBand(hueBandIdx);
1758  const te::rst::Band& lightBand = *inputHLSRaster.getBand(lightBandIdx);
1759  const te::rst::Band& saturationBand = *inputHLSRaster.getBand(saturationBandIdx);
1760  te::rst::Band& redBand = *outputRGBRaster.getBand(0);
1761  te::rst::Band& greenBand = *outputRGBRaster.getBand(1);
1762  te::rst::Band& blueBand = *outputRGBRaster.getBand(2);
1763 
1764  const unsigned int nCols = outputRGBRaster.getNumberOfColumns();
1765  const unsigned int nRows = outputRGBRaster.getNumberOfRows();
1766 
1767  const double hueNoData = hueBand.getProperty()->m_noDataValue;
1768  const double lightNoData = lightBand.getProperty()->m_noDataValue;
1769  const double saturationNoData = saturationBand.getProperty()->m_noDataValue;
1770 
1771  const double rNoData = redBand.getProperty()->m_noDataValue;
1772  const double gNoData = greenBand.getProperty()->m_noDataValue;
1773  const double bNoData = blueBand.getProperty()->m_noDataValue;
1774 
1775  unsigned int row = 0;
1776  unsigned int col = 0;
1777  double hue = 0;
1778  double light = 0;
1779  double saturation = 0;
1780 
1781  for (row = 0; row < nRows; ++row)
1782  {
1783  for (col = 0; col < nCols; ++col)
1784  {
1785  hueBand.getValue(col, row, hue);
1786  lightBand.getValue(col, row, light);
1787  saturationBand.getValue(col, row, saturation);
1788 
1789  if ((hue == hueNoData) || (light == lightNoData) ||
1790  (saturation == saturationNoData))
1791  {
1792  redBand.setValue(col, row, rNoData);
1793  greenBand.setValue(col, row, gNoData);
1794  blueBand.setValue(col, row, bNoData);
1795  }
1796  else
1797  {
1798  std::vector< double > rgb;
1799  hue = hue / 360.;
1800 
1801  if (saturation == 0)
1802  {// achromatic
1803  rgb.push_back(light * rgbRangeMax);
1804  rgb.push_back(light * rgbRangeMax);
1805  rgb.push_back(light * rgbRangeMax);
1806  }
1807  else
1808  {
1809  for (unsigned int i = 0; i < 3; ++i)
1810  {
1811  double q = light < 0.5 ? light * (1. + saturation) : light + saturation - light * saturation;
1812  double p = 2. * light - q;
1813 
1814  double t = 0;
1815 
1816  switch (i)
1817  {
1818  case 0://red
1819  {
1820  t = hue + 1. / 3.;
1821  break;
1822  }
1823  case 1://green
1824  {
1825  t = hue;
1826  break;
1827  }
1828  case 2://blue
1829  {
1830  t = hue - 1. / 3.;
1831  break;
1832  }
1833  }
1834 
1835  if (t < 0)
1836  {
1837  t += 1.;
1838  }
1839  if (t > 1)
1840  {
1841  t -= 1.;
1842  }
1843 
1844  if (t < (1. / 6.))
1845  {
1846  p = (p + (q - p) * 6. * t);
1847  }
1848  else if (t < (1. / 2.))
1849  {
1850  p = q;
1851  }
1852  else if (t < (2. / 3.))
1853  {
1854  p = (p + (q - p) * (2. / 3. - t) * 6.);
1855  }
1856 
1857  rgb.push_back(p * rgbRangeMax);
1858  }
1859  }
1860 
1861  redBand.setValue(col, row, rgb[0]);
1862  greenBand.setValue(col, row, rgb[1]);
1863  blueBand.setValue(col, row, rgb[2]);
1864  }
1865  }
1866  }
1867 
1868  return true;
1869  }
1870 
1871  bool ConvertHLS2RGB(const te::rst::Raster& inputHRaster, const unsigned int hueBandIdx,
1872  const te::rst::Raster& inputLRaster, const unsigned int lightBandIdx,
1873  const te::rst::Raster& inputSRaster, const unsigned int saturationBandIdx,
1874  const double /*rgbRangeMin*/, const double rgbRangeMax, te::rst::Raster& outputRGBRaster)
1875  {
1876  if ((inputHRaster.getAccessPolicy() & te::common::RAccess) == 0 ||
1877  (inputLRaster.getAccessPolicy() & te::common::RAccess) == 0 ||
1878  (inputSRaster.getAccessPolicy() & te::common::RAccess) == 0)
1879  return false;
1880  if (hueBandIdx > inputHRaster.getNumberOfBands()) return false;
1881  if (lightBandIdx > inputLRaster.getNumberOfBands()) return false;
1882  if (saturationBandIdx > inputSRaster.getNumberOfBands()) return false;
1883  if (inputHRaster.getNumberOfColumns() != inputLRaster.getNumberOfColumns() ||
1884  inputHRaster.getNumberOfColumns() != inputSRaster.getNumberOfColumns())
1885  return false;
1886  if (inputHRaster.getNumberOfRows() != inputLRaster.getNumberOfRows() ||
1887  inputHRaster.getNumberOfRows() != inputSRaster.getNumberOfRows())
1888  return false;
1889  if ((outputRGBRaster.getAccessPolicy() & te::common::WAccess) == 0)
1890  return false;
1891  if (outputRGBRaster.getNumberOfBands() < 3) return false;
1892  if (inputHRaster.getNumberOfColumns() != outputRGBRaster.getNumberOfColumns())
1893  return false;
1894  if (inputHRaster.getNumberOfRows() != outputRGBRaster.getNumberOfRows())
1895  return false;
1896 
1897  const te::rst::Band& hueBand = *inputHRaster.getBand(hueBandIdx);
1898  const te::rst::Band& lightBand = *inputLRaster.getBand(lightBandIdx);
1899  const te::rst::Band& saturationBand = *inputSRaster.getBand(saturationBandIdx);
1900  te::rst::Band& redBand = *outputRGBRaster.getBand(0);
1901  te::rst::Band& greenBand = *outputRGBRaster.getBand(1);
1902  te::rst::Band& blueBand = *outputRGBRaster.getBand(2);
1903 
1904  const unsigned int nCols = outputRGBRaster.getNumberOfColumns();
1905  const unsigned int nRows = outputRGBRaster.getNumberOfRows();
1906 
1907  const double hueNoData = hueBand.getProperty()->m_noDataValue;
1908  const double lightNoData = lightBand.getProperty()->m_noDataValue;
1909  const double saturationNoData = saturationBand.getProperty()->m_noDataValue;
1910 
1911  const double rNoData = redBand.getProperty()->m_noDataValue;
1912  const double gNoData = greenBand.getProperty()->m_noDataValue;
1913  const double bNoData = blueBand.getProperty()->m_noDataValue;
1914 
1915  unsigned int row = 0;
1916  unsigned int col = 0;
1917  double hue = 0;
1918  double light = 0;
1919  double saturation = 0;
1920 
1921  for (row = 0; row < nRows; ++row)
1922  {
1923  for (col = 0; col < nCols; ++col)
1924  {
1925  hueBand.getValue(col, row, hue);
1926  lightBand.getValue(col, row, light);
1927  saturationBand.getValue(col, row, saturation);
1928 
1929  if ((hue == hueNoData) || (light == lightNoData) ||
1930  (saturation == saturationNoData))
1931  {
1932  redBand.setValue(col, row, rNoData);
1933  greenBand.setValue(col, row, gNoData);
1934  blueBand.setValue(col, row, bNoData);
1935  }
1936  else
1937  {
1938  std::vector< double > rgb;
1939  hue = hue / 360.;
1940 
1941  if (saturation == 0)
1942  {// achromatic
1943  rgb.push_back(light * rgbRangeMax);
1944  rgb.push_back(light * rgbRangeMax);
1945  rgb.push_back(light * rgbRangeMax);
1946  }
1947  else
1948  {
1949  for (unsigned int i = 0; i < 3; ++i)
1950  {
1951  double q = light < 0.5 ? light * (1. + saturation) : light + saturation - light * saturation;
1952  double p = 2. * light - q;
1953 
1954  double t = 0;
1955 
1956  switch (i)
1957  {
1958  case 0://red
1959  {
1960  t = hue + 1. / 3.;
1961  break;
1962  }
1963  case 1://green
1964  {
1965  t = hue;
1966  break;
1967  }
1968  case 2://blue
1969  {
1970  t = hue - 1. / 3.;
1971  break;
1972  }
1973  }
1974 
1975  if (t < 0)
1976  {
1977  t += 1.;
1978  }
1979  if (t > 1)
1980  {
1981  t -= 1.;
1982  }
1983 
1984  if (t < (1. / 6.))
1985  {
1986  p = (p + (q - p) * 6. * t);
1987  }
1988  else if (t < (1. / 2.))
1989  {
1990  p = q;
1991  }
1992  else if (t < (2. / 3.))
1993  {
1994  p = (p + (q - p) * (2. / 3. - t) * 6.);
1995  }
1996 
1997  rgb.push_back(p * rgbRangeMax);
1998  }
1999  }
2000 
2001  redBand.setValue(col, row, rgb[0]);
2002  greenBand.setValue(col, row, rgb[1]);
2003  blueBand.setValue(col, row, rgb[2]);
2004  }
2005  }
2006  }
2007 
2008  return true;
2009  }
2010 
2011  void GetMeanValueThread( GetMeanValueThreadParams* paramsPtr )
2012  {
2013  paramsPtr->m_mutexPtr->lock();
2014  const unsigned int blockElementsNumber = static_cast<unsigned int>(paramsPtr->m_inputBandPtr->getProperty()->m_blkh *
2015  paramsPtr->m_inputBandPtr->getProperty()->m_blkw);
2016  const int blockDataType = paramsPtr->m_inputBandPtr->getProperty()->getType();
2017  const int blockSizeBytes = paramsPtr->m_inputBandPtr->getBlockSize();
2018  const double noDataValue = paramsPtr->m_inputBandPtr->getProperty()->m_noDataValue;
2019  paramsPtr->m_mutexPtr->unlock();
2020 
2021  boost::scoped_array< unsigned char > blockBuffer( new unsigned char[ blockSizeBytes ] );
2022  boost::scoped_array< double > doubleBuffer( new double[ blockElementsNumber ] );
2023  unsigned blkX = 0;
2024  unsigned int elementIdx = 0;
2025  double mean = 0;
2026  double meanElementsNumber = 0;
2027 
2028  for( unsigned blkY = 0 ; blkY < paramsPtr->m_rasterBlocksStatusPtr->getLinesNumber() ;
2029  ++blkY )
2030  {
2031  for( blkX = 0 ; blkX < paramsPtr->m_rasterBlocksStatusPtr->getColumnsNumber() ;
2032  ++blkX )
2033  {
2034  paramsPtr->m_mutexPtr->lock();
2035 
2036  if( paramsPtr->m_rasterBlocksStatusPtr->operator()( blkY, blkX ) )
2037  {
2038  paramsPtr->m_mutexPtr->unlock();
2039  }
2040  else
2041  {
2042  paramsPtr->m_rasterBlocksStatusPtr->operator()( blkY, blkX ) = true;
2043 
2044  paramsPtr->m_inputBandPtr->read( static_cast<int>(blkX), static_cast<int>(blkY), blockBuffer.get() );
2045 
2046  paramsPtr->m_mutexPtr->unlock();
2047 
2048  Convert2DoublesVector( blockBuffer.get(), blockDataType, blockElementsNumber,
2049  doubleBuffer.get() );
2050 
2051  for( elementIdx = 0 ; elementIdx < blockElementsNumber ; ++elementIdx )
2052  {
2053  if( noDataValue != doubleBuffer[ elementIdx ] )
2054  {
2055  mean =
2056  (
2057  (
2058  mean
2059  *
2060  meanElementsNumber
2061  )
2062  +
2063  doubleBuffer[ elementIdx ]
2064  )
2065  /
2066  (
2067  meanElementsNumber
2068  +
2069  1.0
2070  );
2071 
2072  meanElementsNumber = meanElementsNumber + 1.0;
2073  }
2074  }
2075  }
2076  }
2077  }
2078 
2079  paramsPtr->m_mutexPtr->lock();
2080  ( *(paramsPtr->m_meanValuePtr) ) =
2081  (
2082  (
2083  ( *(paramsPtr->m_meanValuePtr) )
2084  *
2085  ( *( paramsPtr->m_pixelsNumberValuePtr ) )
2086  )
2087  +
2088  (
2089  mean
2090  *
2091  meanElementsNumber
2092  )
2093  )
2094  /
2095  (
2096  ( *( paramsPtr->m_pixelsNumberValuePtr ) )
2097  +
2098  meanElementsNumber
2099  );
2100  ( *( paramsPtr->m_pixelsNumberValuePtr ) ) += meanElementsNumber;
2101  paramsPtr->m_mutexPtr->unlock();
2102  }
2103 
2104  bool GetMeanValue( const te::rst::Band& band,
2105  const unsigned int maxThreads,
2106  double& meanValue )
2107  {
2108  Matrix< bool > rasterBlocksStatus;
2109  if( ! rasterBlocksStatus.reset( static_cast<unsigned int>(band.getProperty()->m_nblocksy),
2110  static_cast<unsigned int>(band.getProperty()->m_nblocksx), Matrix< bool >::RAMMemPol ) )
2111  {
2112  return false;
2113  }
2114  for( unsigned int row = 0 ; row < rasterBlocksStatus.getLinesNumber() ;
2115  ++row )
2116  {
2117  for( unsigned int col = 0 ; col < rasterBlocksStatus.getColumnsNumber() ;
2118  ++col )
2119  {
2120  rasterBlocksStatus( row, col ) = false;
2121  }
2122  }
2123 
2124  meanValue = 0.0;
2125 
2126  double pixelsNumber = 0.0;
2127 
2128  boost::mutex mutex;
2129 
2130  GetMeanValueThreadParams params;
2131  params.m_inputBandPtr = &band;
2132  params.m_rasterBlocksStatusPtr = &rasterBlocksStatus;
2133  params.m_meanValuePtr = &meanValue;
2134  params.m_pixelsNumberValuePtr = &pixelsNumber;
2135  params.m_returnStatus = true;
2136  params.m_mutexPtr = &mutex;
2137 
2138  if( maxThreads == 1 )
2139  {
2140  GetMeanValueThread( &params );
2141  }
2142  else
2143  {
2144  unsigned int threadsNumber = maxThreads ? maxThreads : te::common::GetPhysProcNumber();
2145 
2146  boost::thread_group threads;
2147 
2148  for( unsigned int threadIdx = 0 ; threadIdx < threadsNumber ;
2149  ++threadIdx )
2150  {
2151  threads.add_thread( new boost::thread( GetMeanValueThread,
2152  &params ) );
2153  };
2154 
2155  threads.join_all();
2156  }
2157 
2158  return params.m_returnStatus;
2159  }
2160 
2161  void GetStdDevValueThread( GetStdDevValueThreadParams* paramsPtr )
2162  {
2163  paramsPtr->m_mutexPtr->lock();
2164  const unsigned int blockElementsNumber = static_cast<unsigned int>(paramsPtr->m_inputBandPtr->getProperty()->m_blkh *
2165  paramsPtr->m_inputBandPtr->getProperty()->m_blkw);
2166  const int blockDataType = paramsPtr->m_inputBandPtr->getProperty()->getType();
2167  const int blockSizeBytes = paramsPtr->m_inputBandPtr->getBlockSize();
2168  const double noDataValue = paramsPtr->m_inputBandPtr->getProperty()->m_noDataValue;
2169  const double& meanValue = paramsPtr->m_meanValue;
2170  paramsPtr->m_mutexPtr->unlock();
2171 
2172  boost::scoped_array< unsigned char > blockBuffer( new unsigned char[ blockSizeBytes ] );
2173  boost::scoped_array< double > doubleBuffer( new double[ blockElementsNumber ] );
2174  unsigned blkX = 0;
2175  unsigned int elementIdx = 0;
2176  double diff = 0;
2177  double squaresDifsSum = 0;
2178  double elementsNumber = 0;
2179 
2180 
2181  for( unsigned blkY = 0 ; blkY < paramsPtr->m_rasterBlocksStatusPtr->getLinesNumber() ;
2182  ++blkY )
2183  {
2184  for( blkX = 0 ; blkX < paramsPtr->m_rasterBlocksStatusPtr->getColumnsNumber() ;
2185  ++blkX )
2186  {
2187  paramsPtr->m_mutexPtr->lock();
2188 
2189  if( paramsPtr->m_rasterBlocksStatusPtr->operator()( blkY, blkX ) )
2190  {
2191  paramsPtr->m_mutexPtr->unlock();
2192  }
2193  else
2194  {
2195  paramsPtr->m_rasterBlocksStatusPtr->operator()( blkY, blkX ) = true;
2196 
2197  paramsPtr->m_inputBandPtr->read( static_cast<int>(blkX), static_cast<int>(blkY), blockBuffer.get() );
2198 
2199  paramsPtr->m_mutexPtr->unlock();
2200 
2201  Convert2DoublesVector( blockBuffer.get(), blockDataType, blockElementsNumber,
2202  doubleBuffer.get() );
2203 
2204  for( elementIdx = 0 ; elementIdx < blockElementsNumber ; ++elementIdx )
2205  {
2206  if( noDataValue != doubleBuffer[ elementIdx ] )
2207  {
2208  diff = doubleBuffer[ elementIdx ] - meanValue;
2209  diff *= diff;
2210  squaresDifsSum += diff;
2211  elementsNumber = elementsNumber + 1.0;
2212  }
2213  }
2214  }
2215  }
2216  }
2217 
2218  paramsPtr->m_mutexPtr->lock();
2219  ( *(paramsPtr->m_stdDevValuePtr) ) =
2220  (
2221  ( *(paramsPtr->m_stdDevValuePtr) )
2222  +
2223  squaresDifsSum
2224  );
2225  ( *( paramsPtr->m_pixelsNumberValuePtr ) ) += elementsNumber;
2226  paramsPtr->m_mutexPtr->unlock();
2227  }
2228 
2229  bool GetStdDevValue( const te::rst::Band& band,
2230  const unsigned int maxThreads, double const * const meanValuePtr,
2231  double& stdDevValue )
2232  {
2233  double mean = 0;
2234  if( meanValuePtr )
2235  {
2236  mean = (*meanValuePtr);
2237  }
2238  else
2239  {
2240  if( ! GetMeanValue( band, maxThreads, mean ) )
2241  {
2242  return false;
2243  }
2244  }
2245 
2246  Matrix< bool > rasterBlocksStatus;
2247  if( ! rasterBlocksStatus.reset( static_cast<unsigned int>(band.getProperty()->m_nblocksy),
2248  static_cast<unsigned int>(band.getProperty()->m_nblocksx), Matrix< bool >::RAMMemPol ) )
2249  {
2250  return false;
2251  }
2252  for( unsigned int row = 0 ; row < rasterBlocksStatus.getLinesNumber() ;
2253  ++row )
2254  {
2255  for( unsigned int col = 0 ; col < rasterBlocksStatus.getColumnsNumber() ;
2256  ++col )
2257  {
2258  rasterBlocksStatus( row, col ) = false;
2259  }
2260  }
2261 
2262  stdDevValue = 0.0;
2263 
2264  double pixelsNumber = 0.0;
2265 
2266  boost::mutex mutex;
2267 
2268  GetStdDevValueThreadParams params;
2269  params.m_inputBandPtr = &band;
2270  params.m_rasterBlocksStatusPtr = &rasterBlocksStatus;
2271  params.m_stdDevValuePtr = &stdDevValue;
2272  params.m_pixelsNumberValuePtr = &pixelsNumber;
2273  params.m_meanValue = mean;
2274  params.m_returnStatus = true;
2275  params.m_mutexPtr = &mutex;
2276 
2277  if( maxThreads == 1 )
2278  {
2279  GetStdDevValueThread( &params );
2280  }
2281  else
2282  {
2283  unsigned int threadsNumber = maxThreads ? maxThreads : te::common::GetPhysProcNumber();
2284 
2285  boost::thread_group threads;
2286 
2287  for( unsigned int threadIdx = 0 ; threadIdx < threadsNumber ;
2288  ++threadIdx )
2289  {
2290  threads.add_thread( new boost::thread( GetStdDevValueThread,
2291  &params ) );
2292  };
2293 
2294  threads.join_all();
2295  }
2296 
2297  stdDevValue = std::sqrt( stdDevValue / pixelsNumber );
2298 
2299  return params.m_returnStatus;
2300  }
2301 
2302  void GetCovarianceValueThread( GetCovarianceValueThreadParams* paramsPtr )
2303  {
2304  paramsPtr->m_mutexPtr->lock();
2305  const unsigned int blockElementsNumber = static_cast<unsigned int>(paramsPtr->m_inputBand1Ptr->getProperty()->m_blkh *
2306  paramsPtr->m_inputBand1Ptr->getProperty()->m_blkw);
2307  const int blockDataType1 = paramsPtr->m_inputBand1Ptr->getProperty()->getType();
2308  const int blockDataType2 = paramsPtr->m_inputBand2Ptr->getProperty()->getType();
2309  const int blockSizeBytes1 = paramsPtr->m_inputBand1Ptr->getBlockSize();
2310  const int blockSizeBytes2 = paramsPtr->m_inputBand2Ptr->getBlockSize();
2311  const double noDataValue1 = paramsPtr->m_inputBand1Ptr->getProperty()->m_noDataValue;
2312  const double noDataValue2 = paramsPtr->m_inputBand2Ptr->getProperty()->m_noDataValue;
2313  paramsPtr->m_mutexPtr->unlock();
2314 
2315  boost::scoped_array< unsigned char > blockBuffer1( new unsigned char[ blockSizeBytes1 ] );
2316  boost::scoped_array< unsigned char > blockBuffer2( new unsigned char[ blockSizeBytes2 ] );
2317  boost::scoped_array< double > doubleBuffer1( new double[ blockElementsNumber ] );
2318  boost::scoped_array< double > doubleBuffer2( new double[ blockElementsNumber ] );
2319  unsigned blkX = 0;
2320  unsigned int elementIdx = 0;
2321  double covariance = 0;
2322  double elementsNumber = 0;
2323  double diff1 = 0;
2324  double diff2 = 0;
2325 
2326  for( unsigned blkY = 0 ; blkY < paramsPtr->m_rasterBlocksStatusPtr->getLinesNumber() ;
2327  ++blkY )
2328  {
2329  for( blkX = 0 ; blkX < paramsPtr->m_rasterBlocksStatusPtr->getColumnsNumber() ;
2330  ++blkX )
2331  {
2332  paramsPtr->m_mutexPtr->lock();
2333 
2334  if( paramsPtr->m_rasterBlocksStatusPtr->operator()( blkY, blkX ) )
2335  {
2336  paramsPtr->m_mutexPtr->unlock();
2337  }
2338  else
2339  {
2340  paramsPtr->m_rasterBlocksStatusPtr->operator()( blkY, blkX ) = true;
2341 
2342  paramsPtr->m_inputBand1Ptr->read( static_cast<int>(blkX), static_cast<int>(blkY), blockBuffer1.get() );
2343  paramsPtr->m_inputBand2Ptr->read( static_cast<int>(blkX), static_cast<int>(blkY), blockBuffer2.get() );
2344 
2345  paramsPtr->m_mutexPtr->unlock();
2346 
2347  Convert2DoublesVector( blockBuffer1.get(), blockDataType1, blockElementsNumber,
2348  doubleBuffer1.get() );
2349  Convert2DoublesVector( blockBuffer2.get(), blockDataType2, blockElementsNumber,
2350  doubleBuffer2.get() );
2351 
2352  for( elementIdx = 0 ; elementIdx < blockElementsNumber ; ++elementIdx )
2353  {
2354  if( ( noDataValue1 != doubleBuffer1[ elementIdx ] ) &&
2355  ( noDataValue2 != doubleBuffer2[ elementIdx ] ) )
2356  {
2357  diff1 = doubleBuffer1[ elementIdx ] - paramsPtr->m_mean1Value;
2358  diff2 = doubleBuffer2[ elementIdx ] - paramsPtr->m_mean2Value;
2359 
2360  covariance += ( diff1 * diff2 );
2361 
2362  elementsNumber = elementsNumber + 1.0;
2363  }
2364  }
2365  }
2366  }
2367  }
2368 
2369  paramsPtr->m_mutexPtr->lock();
2370  ( *(paramsPtr->m_covarianceValuePtr) ) += covariance;
2371  ( *( paramsPtr->m_pixelsNumberValuePtr ) ) += elementsNumber;
2372  paramsPtr->m_mutexPtr->unlock();
2373  }
2374 
2375  bool GetCovarianceValue( const te::rst::Band& band1,
2376  const te::rst::Band& band2, const unsigned int maxThreads,
2377  double const * const mean1ValuePtr, double const * const mean2ValuePtr,
2378  double& covarianceValue )
2379  {
2380  if( ( band1.getProperty()->m_blkw * band1.getProperty()->m_nblocksx ) !=
2381  ( band2.getProperty()->m_blkw * band2.getProperty()->m_nblocksx ) )
2382  {
2383  return false;
2384  }
2385 
2386  if( ( band1.getProperty()->m_blkh * band1.getProperty()->m_nblocksy ) !=
2387  ( band2.getProperty()->m_blkh * band2.getProperty()->m_nblocksy ) )
2388  {
2389  return false;
2390  }
2391 
2392  double mean1 = 0;
2393  if( mean1ValuePtr )
2394  {
2395  mean1 = (*mean1ValuePtr);
2396  }
2397  else
2398  {
2399  if( ! GetMeanValue( band1, maxThreads, mean1 ) )
2400  {
2401  return false;
2402  }
2403  }
2404 
2405  double mean2 = 0;
2406  if( mean2ValuePtr )
2407  {
2408  mean2 = (*mean2ValuePtr);
2409  }
2410  else
2411  {
2412  if( ! GetMeanValue( band2, maxThreads, mean2 ) )
2413  {
2414  return false;
2415  }
2416  }
2417 
2418  if( ( band1.getProperty()->m_blkw == band2.getProperty()->m_blkw )
2419  && ( band1.getProperty()->m_blkh == band2.getProperty()->m_blkh )
2420  && ( band1.getProperty()->m_nblocksx == band2.getProperty()->m_nblocksx )
2421  && ( band1.getProperty()->m_nblocksy == band2.getProperty()->m_nblocksy ) )
2422  {
2423  Matrix< bool > rasterBlocksStatus;
2424  if( ! rasterBlocksStatus.reset( static_cast<unsigned int>(band1.getProperty()->m_nblocksy),
2425  static_cast<unsigned int>(band1.getProperty()->m_nblocksx), Matrix< bool >::RAMMemPol ) )
2426  {
2427  return false;
2428  }
2429  for( unsigned int row = 0 ; row < rasterBlocksStatus.getLinesNumber() ;
2430  ++row )
2431  {
2432  for( unsigned int col = 0 ; col < rasterBlocksStatus.getColumnsNumber() ;
2433  ++col )
2434  {
2435  rasterBlocksStatus( row, col ) = false;
2436  }
2437  }
2438 
2439  covarianceValue = 0.0;
2440 
2441  double pixelsNumber = 0.0;
2442 
2443  boost::mutex mutex;
2444 
2445  GetCovarianceValueThreadParams params;
2446  params.m_inputBand1Ptr = &band1;
2447  params.m_inputBand2Ptr = &band2;
2448  params.m_rasterBlocksStatusPtr = &rasterBlocksStatus;
2449  params.m_covarianceValuePtr = &covarianceValue;
2450  params.m_pixelsNumberValuePtr = &pixelsNumber;
2451  params.m_returnStatus = true;
2452  params.m_mutexPtr = &mutex;
2453  params.m_mean1Value = mean1;
2454  params.m_mean2Value = mean2;
2455 
2456  if( maxThreads == 1 )
2457  {
2458  GetCovarianceValueThread( &params );
2459  }
2460  else
2461  {
2462  unsigned int threadsNumber = maxThreads ? maxThreads : te::common::GetPhysProcNumber();
2463 
2464  boost::thread_group threads;
2465 
2466  for( unsigned int threadIdx = 0 ; threadIdx < threadsNumber ;
2467  ++threadIdx )
2468  {
2469  threads.add_thread( new boost::thread( GetCovarianceValueThread,
2470  &params ) );
2471  };
2472 
2473  threads.join_all();
2474  }
2475 
2476  if( pixelsNumber != 0.0 )
2477  {
2478  covarianceValue /= pixelsNumber;
2479  }
2480 
2481  return params.m_returnStatus;
2482  }
2483  else
2484  {
2485  const unsigned int nCols = static_cast<unsigned int>(band1.getProperty()->m_blkw * band1.getProperty()->m_nblocksx);
2486  const unsigned int nRows = static_cast<unsigned int>(band1.getProperty()->m_blkh * band1.getProperty()->m_nblocksy);
2487  double value1 = 0;
2488  double value2 = 0;
2489  const double noDataValue1 = band1.getProperty()->m_noDataValue;
2490  const double noDataValue2 = band2.getProperty()->m_noDataValue;
2491  unsigned int col = 0;
2492  double pixelsNumber = 0.0;
2493 
2494  covarianceValue = 0;
2495 
2496  for( unsigned int row = 0 ; row < nRows ; ++row )
2497  {
2498  for( col = 0 ; col < nCols ; ++col )
2499  {
2500  band1.getValue( col, row, value1 );
2501  band2.getValue( col, row, value2 );
2502 
2503  if( ( noDataValue1 != value1 ) &&
2504  ( noDataValue2 != value2 ) )
2505  {
2506  covarianceValue += ( ( value1 - mean1 ) * ( value2 - mean2 ) );
2507 
2508  pixelsNumber = pixelsNumber + 1.0;
2509  }
2510  }
2511  }
2512 
2513  if( pixelsNumber != 0.0 )
2514  {
2515  covarianceValue /= pixelsNumber;
2516  }
2517 
2518  return true;
2519  }
2520  }
2521 
2522  bool DirectPrincipalComponents(
2523  const te::rst::Raster& inputRaster,
2524  const std::vector< unsigned int >& inputRasterBands,
2525  boost::numeric::ublas::matrix< double >& pcaMatrix,
2526  te::rst::Raster& pcaRaster,
2527  const std::vector< unsigned int >& pcaRasterBands,
2528  const unsigned int maxThreads )
2529  {
2530  if( ( inputRaster.getAccessPolicy() & te::common::RAccess ) == 0 )
2531  {
2532  return false;
2533  }
2534  if( ( pcaRaster.getAccessPolicy() & te::common::WAccess ) == 0 )
2535  {
2536  return false;
2537  }
2538  if( inputRaster.getNumberOfBands() != pcaRaster.getNumberOfBands() )
2539  {
2540  return false;
2541  }
2542  if( pcaRaster.getNumberOfBands() != inputRasterBands.size() )
2543  {
2544  return false;
2545  }
2546  if( inputRaster.getNumberOfRows() != pcaRaster.getNumberOfRows() )
2547  {
2548  return false;
2549  }
2550  if( inputRaster.getNumberOfColumns() != pcaRaster.getNumberOfColumns() )
2551  {
2552  return false;
2553  }
2554 
2555  // Covariance matrix
2556 
2557  boost::numeric::ublas::matrix< double > covarianceMatrix;
2558 
2559  covarianceMatrix.resize( inputRasterBands.size(), inputRasterBands.size() );
2560 
2561  for( unsigned int covMatrixIdx1 = 0 ; covMatrixIdx1 < inputRasterBands.size() ;
2562  ++covMatrixIdx1 )
2563  {
2564  if( inputRasterBands[ covMatrixIdx1 ] >= inputRaster.getNumberOfBands() )
2565  {
2566  return false;
2567  }
2568 
2569  for( unsigned int covMatrixIdx2 = 0 ; covMatrixIdx2 < inputRasterBands.size() ;
2570  ++covMatrixIdx2 )
2571  {
2572  if( inputRasterBands[ covMatrixIdx2 ] >= inputRaster.getNumberOfBands() )
2573  {
2574  return false;
2575  }
2576 
2577  if( covMatrixIdx1 > covMatrixIdx2 )
2578  {
2579  covarianceMatrix( covMatrixIdx1, covMatrixIdx2 ) =
2580  covarianceMatrix( covMatrixIdx2, covMatrixIdx1 );
2581  }
2582  else
2583  {
2584  if( ! GetCovarianceValue( *( inputRaster.getBand( inputRasterBands[ covMatrixIdx1 ] ) ),
2585  *( inputRaster.getBand( inputRasterBands[ covMatrixIdx2 ] ) ),
2586  maxThreads, nullptr, nullptr, covarianceMatrix( covMatrixIdx1, covMatrixIdx2 ) ) )
2587  {
2588  return false;
2589  }
2590  }
2591  }
2592  }
2593 
2594 //std::cout << std::endl << "Covariance matrix:" << covarianceMatrix << std::endl;
2595 
2596  // Eigen stuff
2597 
2598  boost::numeric::ublas::matrix< double > eigenValues;
2599  boost::numeric::ublas::matrix< double > eigenVectors;
2600 
2601  if( ! te::common::EigenVectors( covarianceMatrix, eigenVectors, eigenValues ) )
2602  {
2603  return false;
2604  }
2605 
2606  pcaMatrix = boost::numeric::ublas::trans( eigenVectors );
2607 
2608  return RemapValues( inputRaster, inputRasterBands, pcaMatrix, pcaRaster,
2609  pcaRasterBands, maxThreads );
2610  }
2611 
2612  bool InversePrincipalComponents(
2613  const te::rst::Raster& pcaRaster,
2614  const boost::numeric::ublas::matrix< double >& pcaMatrix,
2615  te::rst::Raster& outputRaster,
2616  const std::vector< unsigned int >& outputRasterBands,
2617  const unsigned int maxThreads )
2618  {
2619  boost::numeric::ublas::matrix< double > inversePcaMatrix;
2620  if( ! te::common::GetInverseMatrix( pcaMatrix, inversePcaMatrix ) )
2621  {
2622  return false;
2623  }
2624 
2625  std::vector< unsigned int > inputRasterBands;
2626 
2627  for( unsigned int bandIdx = 0 ; bandIdx < pcaRaster.getNumberOfBands() ;
2628  ++bandIdx )
2629  {
2630  inputRasterBands.push_back( bandIdx );
2631  }
2632 
2633  return RemapValues( pcaRaster, inputRasterBands, inversePcaMatrix,
2634  outputRaster, outputRasterBands, maxThreads );
2635  }
2636 
2637  void RemapValuesThread( RemapValuesThreadParams* paramsPtr )
2638  {
2639  paramsPtr->m_mutexPtr->lock();
2640 
2641  assert( paramsPtr->m_inputRasterBandsPtr->operator[]( 0 ) <
2642  paramsPtr->m_inputRasterPtr->getNumberOfBands() );
2643  const unsigned int blockElementsNumber = static_cast<unsigned int>(paramsPtr->m_inputRasterPtr->getBand(
2644  paramsPtr->m_inputRasterBandsPtr->operator[]( 0 ) )->getProperty()->m_blkh *
2645  paramsPtr->m_inputRasterPtr->getBand( paramsPtr->m_inputRasterBandsPtr->operator[]( 0 ) )->getProperty()->m_blkw);
2646 
2647  const std::vector< unsigned int > inputRasterBands = *( paramsPtr->m_inputRasterBandsPtr );
2648  const unsigned int inputRasterBandsSize = static_cast<unsigned int>(inputRasterBands.size());
2649  const std::vector< unsigned int > outputRasterBands = *( paramsPtr->m_inputRasterBandsPtr );
2650  assert( inputRasterBandsSize == paramsPtr->m_outputRasterPtr->getNumberOfBands() );
2651  assert( inputRasterBandsSize == outputRasterBands.size() );
2652 
2653  unsigned int maxBlocksSizesBytes = 0;
2654  std::vector< double > inputBandsNoDataValues( inputRasterBandsSize, 0.0 );
2655  std::vector< double > outputBandsNoDataValues( inputRasterBandsSize, 0.0 );
2656  std::vector< boost::shared_array< double > > inDoubleBuffersHandlers( inputRasterBandsSize );
2657  std::vector< boost::shared_array< double > > outDoubleBuffersHandlers( inputRasterBandsSize );
2658  unsigned int inputRasterBandsIdx = 0;
2659  boost::numeric::ublas::matrix< double > remapMatrix = *( paramsPtr->m_remapMatrixPtr );
2660  std::vector< double > outputMinValue( inputRasterBandsSize );
2661  std::vector< double > outputMaxValue( inputRasterBandsSize );
2662  boost::shared_array< double* > inDoubleBuffersPtrsHandler( new double*[ inputRasterBandsSize ] );
2663  boost::shared_array< double* > outDoubleBuffersPtrsHandler( new double*[ inputRasterBandsSize ] );
2664  double** inDoubleBuffers = inDoubleBuffersPtrsHandler.get();
2665  double** outDoubleBuffers = outDoubleBuffersPtrsHandler.get();
2666 
2667  for( inputRasterBandsIdx = 0 ; inputRasterBandsIdx < inputRasterBandsSize ;
2668  ++inputRasterBandsIdx )
2669  {
2670  const unsigned int& inBandIdx = inputRasterBands[ inputRasterBandsIdx ];
2671  assert( inBandIdx < paramsPtr->m_inputRasterPtr->getNumberOfBands() );
2672 
2673  te::rst::Band const * const inBandPtr = paramsPtr->m_inputRasterPtr->getBand( inBandIdx );
2674  te::rst::Band * const outBandPtr = paramsPtr->m_outputRasterPtr->getBand(
2675  outputRasterBands[ inputRasterBandsIdx ] );
2676 
2677  maxBlocksSizesBytes = std::max( maxBlocksSizesBytes, static_cast<unsigned int>(inBandPtr->getBlockSize()) );
2678  maxBlocksSizesBytes = std::max( maxBlocksSizesBytes, static_cast<unsigned int>(outBandPtr->getBlockSize()) );
2679 
2680  inputBandsNoDataValues[ inputRasterBandsIdx ] = inBandPtr->getProperty()->m_noDataValue;
2681  outputBandsNoDataValues[ inputRasterBandsIdx ] = outBandPtr->getProperty()->m_noDataValue;
2682 
2683  inDoubleBuffersHandlers[ inputRasterBandsIdx ].reset( new double[ blockElementsNumber ] );
2684  inDoubleBuffers[ inputRasterBandsIdx ] = inDoubleBuffersHandlers[ inputRasterBandsIdx ].get();
2685 
2686  outDoubleBuffersHandlers[ inputRasterBandsIdx ].reset( new double[ blockElementsNumber ] );
2687  outDoubleBuffers[ inputRasterBandsIdx ] = outDoubleBuffersHandlers[ inputRasterBandsIdx ].get();
2688 
2689  te::rst::GetDataTypeRanges(
2690  outBandPtr->getProperty()->getType(),
2691  outputMinValue[ inputRasterBandsIdx ],
2692  outputMaxValue[ inputRasterBandsIdx ] );
2693  }
2694 
2695  paramsPtr->m_mutexPtr->unlock();
2696 
2697  boost::scoped_array< unsigned char > blockBuffer( new unsigned char[ maxBlocksSizesBytes ] );
2698  unsigned blkX = 0;
2699  unsigned int elementIdx = 0;
2700  boost::numeric::ublas::matrix< double > pixelValues( inputRasterBandsSize, 1 );
2701  boost::numeric::ublas::matrix< double > remappedPixelValues( inputRasterBandsSize, 1 );
2702  bool elementIsValid = false;
2703 
2704  for( unsigned blkY = 0 ; blkY < paramsPtr->m_rasterBlocksStatusPtr->getLinesNumber() ;
2705  ++blkY )
2706  {
2707  for( blkX = 0 ; blkX < paramsPtr->m_rasterBlocksStatusPtr->getColumnsNumber() ;
2708  ++blkX )
2709  {
2710  paramsPtr->m_mutexPtr->lock();
2711 
2712  if( paramsPtr->m_rasterBlocksStatusPtr->operator()( blkY, blkX ) )
2713  {
2714  paramsPtr->m_mutexPtr->unlock();
2715  }
2716  else
2717  {
2718  paramsPtr->m_rasterBlocksStatusPtr->operator()( blkY, blkX ) = true;
2719 
2720  // reading the input blocks
2721 
2722  for( inputRasterBandsIdx = 0 ; inputRasterBandsIdx < inputRasterBandsSize ;
2723  ++inputRasterBandsIdx )
2724  {
2725  const unsigned int& inBandIdx = inputRasterBands[ inputRasterBandsIdx ];
2726  te::rst::Band const * const inBandPtr = paramsPtr->m_inputRasterPtr->getBand( inBandIdx );
2727 
2728  inBandPtr->read( static_cast<int>(blkX), static_cast<int>(blkY), blockBuffer.get() );
2729  Convert2DoublesVector( blockBuffer.get(),
2730  inBandPtr->getProperty()->getType(),
2731  blockElementsNumber,
2732  inDoubleBuffers[ inputRasterBandsIdx ] );
2733  }
2734 
2735  paramsPtr->m_mutexPtr->unlock();
2736 
2737  // Remapping the values using the remap matrix
2738 
2739  for( elementIdx = 0 ; elementIdx < blockElementsNumber ; ++elementIdx )
2740  {
2741  elementIsValid = true;
2742 
2743  for( inputRasterBandsIdx = 0 ; inputRasterBandsIdx < inputRasterBandsSize ;
2744  ++inputRasterBandsIdx )
2745  {
2746  if( inDoubleBuffers[ inputRasterBandsIdx ][ elementIdx ] ==
2747  inputBandsNoDataValues[ inputRasterBandsIdx ] )
2748  {
2749  elementIsValid = false;
2750  break;
2751  }
2752  else
2753  {
2754  pixelValues( inputRasterBandsIdx, 0 ) =
2755  inDoubleBuffers[ inputRasterBandsIdx ][ elementIdx ];
2756  }
2757  }
2758 
2759  if( elementIsValid )
2760  {
2761 //std::cout << std::endl << "pixelValues:" << pixelValues << std::endl;
2762 
2763  remappedPixelValues = boost::numeric::ublas::prod( remapMatrix, pixelValues );
2764 
2765 //std::cout << std::endl << "remappedPixelValues:" << remappedPixelValues << std::endl;
2766 
2767  for( inputRasterBandsIdx = 0 ; inputRasterBandsIdx < inputRasterBandsSize ;
2768  ++inputRasterBandsIdx )
2769  {
2770  outDoubleBuffers[ inputRasterBandsIdx ][ elementIdx ] =
2771  std::max(
2772  outputMinValue[ inputRasterBandsIdx ],
2773  std::min(
2774  outputMaxValue[ inputRasterBandsIdx ],
2775  remappedPixelValues( inputRasterBandsIdx, 0 )
2776  )
2777  );
2778  }
2779  }
2780  else
2781  {
2782  for( inputRasterBandsIdx = 0 ; inputRasterBandsIdx < inputRasterBandsSize ;
2783  ++inputRasterBandsIdx )
2784  {
2785  outDoubleBuffers[ inputRasterBandsIdx ][ elementIdx ] =
2786  outputBandsNoDataValues[ inputRasterBandsIdx ] ;
2787  }
2788  }
2789  }
2790 
2791  // Writing the remmaped blocks
2792 
2793  paramsPtr->m_mutexPtr->lock();
2794 
2795  for( inputRasterBandsIdx = 0 ; inputRasterBandsIdx < inputRasterBandsSize ;
2796  ++inputRasterBandsIdx )
2797  {
2798  te::rst::Band* outBandPtr = paramsPtr->m_outputRasterPtr->getBand(
2799  outputRasterBands[ inputRasterBandsIdx ] );
2800 
2801  ConvertDoublesVector( outDoubleBuffers[ inputRasterBandsIdx ],
2802  blockElementsNumber, outBandPtr->getProperty()->getType(),
2803  blockBuffer.get() );
2804 
2805  outBandPtr->write( static_cast<int>(blkX), static_cast<int>(blkY), blockBuffer.get() );
2806  }
2807 
2808  paramsPtr->m_mutexPtr->unlock();
2809  }
2810  }
2811  }
2812 
2813  }
2814 
2815  bool RemapValues(
2816  const te::rst::Raster& inputRaster,
2817  const std::vector< unsigned int >& inputRasterBands,
2818  const boost::numeric::ublas::matrix< double >& remapMatrix,
2819  te::rst::Raster& outputRaster,
2820  const std::vector< unsigned int >& outputRasterBands,
2821  const unsigned int maxThreads )
2822  {
2823  if( ( inputRaster.getAccessPolicy() & te::common::RAccess ) == 0 )
2824  {
2825  return false;
2826  }
2827  if( ( outputRaster.getAccessPolicy() & te::common::WAccess ) == 0 )
2828  {
2829  return false;
2830  }
2831  if( remapMatrix.size1() != inputRasterBands.size() )
2832  {
2833  return false;
2834  }
2835  if( remapMatrix.size2() != inputRasterBands.size() )
2836  {
2837  return false;
2838  }
2839  if( outputRaster.getNumberOfBands() != inputRasterBands.size() )
2840  {
2841  return false;
2842  }
2843  if( inputRaster.getNumberOfRows() != outputRaster.getNumberOfRows() )
2844  {
2845  return false;
2846  }
2847  if( inputRaster.getNumberOfColumns() != outputRaster.getNumberOfColumns() )
2848  {
2849  return false;
2850  }
2851  if( remapMatrix.size1() != outputRasterBands.size() )
2852  {
2853  return false;
2854  }
2855  if( remapMatrix.size2() != outputRasterBands.size() )
2856  {
2857  return false;
2858  }
2859  for( unsigned int inputRasterBandsIdx = 0 ; inputRasterBandsIdx < inputRasterBands.size() ;
2860  ++inputRasterBandsIdx )
2861  {
2862  if( inputRasterBands[ inputRasterBandsIdx ] >= inputRaster.getNumberOfBands() )
2863  {
2864  return false;
2865  }
2866  }
2867  for( unsigned int outputRasterBandsIdx = 0 ; outputRasterBandsIdx < outputRasterBands.size() ;
2868  ++outputRasterBandsIdx )
2869  {
2870  if( outputRasterBands[ outputRasterBandsIdx ] >= outputRaster.getNumberOfBands() )
2871  {
2872  return false;
2873  }
2874  }
2875 
2876  // Checking if optmized PCA can be executed
2877 
2878  bool useOptimizedRemap = true;
2879 
2880  {
2881  for( unsigned int inputRasterBandsIdx = 0 ; inputRasterBandsIdx <
2882  inputRaster.getNumberOfBands() ; ++inputRasterBandsIdx )
2883  {
2884  if(
2885  (
2886  inputRaster.getBand( inputRasterBands[ inputRasterBandsIdx ] )->getProperty()->m_blkw
2887  !=
2888  inputRaster.getBand( inputRasterBands[ 0 ] )->getProperty()->m_blkw
2889  )
2890  ||
2891  (
2892  inputRaster.getBand( inputRasterBands[ inputRasterBandsIdx ] )->getProperty()->m_blkh
2893  !=
2894  inputRaster.getBand( inputRasterBands[ 0 ] )->getProperty()->m_blkh
2895  )
2896  ||
2897  (
2898  inputRaster.getBand( inputRasterBands[ inputRasterBandsIdx ] )->getProperty()->m_nblocksx
2899  !=
2900  inputRaster.getBand( inputRasterBands[ 0 ] )->getProperty()->m_nblocksx
2901  )
2902  ||
2903  (
2904  inputRaster.getBand( inputRasterBands[ inputRasterBandsIdx ] )->getProperty()->m_nblocksy
2905  !=
2906  inputRaster.getBand( inputRasterBands[ 0 ] )->getProperty()->m_nblocksy
2907  )
2908  ||
2909  (
2910  outputRaster.getBand( inputRasterBandsIdx )->getProperty()->m_blkw
2911  !=
2912  inputRaster.getBand( inputRasterBands[ 0 ] )->getProperty()->m_blkw
2913  )
2914  ||
2915  (
2916  outputRaster.getBand( inputRasterBandsIdx )->getProperty()->m_blkh
2917  !=
2918  inputRaster.getBand( inputRasterBands[ 0 ] )->getProperty()->m_blkh
2919  )
2920  ||
2921  (
2922  outputRaster.getBand( inputRasterBandsIdx )->getProperty()->m_nblocksx
2923  !=
2924  inputRaster.getBand( inputRasterBands[ 0 ] )->getProperty()->m_nblocksx
2925  )
2926  ||
2927  (
2928  outputRaster.getBand( inputRasterBandsIdx )->getProperty()->m_nblocksy
2929  !=
2930  inputRaster.getBand( inputRasterBands[ 0 ] )->getProperty()->m_nblocksy
2931  )
2932  )
2933  {
2934  useOptimizedRemap = false;
2935  break;
2936  }
2937  }
2938  }
2939 
2940  // remapping
2941 
2942  if( useOptimizedRemap )
2943  {
2944  Matrix< bool > rasterBlocksStatus;
2945  if( ! rasterBlocksStatus.reset(
2946  static_cast<unsigned int>(inputRaster.getBand( inputRasterBands[ 0 ] )->getProperty()->m_nblocksy),
2947  static_cast<unsigned int>(inputRaster.getBand( inputRasterBands[ 0 ] )->getProperty()->m_nblocksx),
2948  Matrix< bool >::RAMMemPol ) )
2949  {
2950  return false;
2951  }
2952  for( unsigned int row = 0 ; row < rasterBlocksStatus.getLinesNumber() ;
2953  ++row )
2954  {
2955  for( unsigned int col = 0 ; col < rasterBlocksStatus.getColumnsNumber() ;
2956  ++col )
2957  {
2958  rasterBlocksStatus( row, col ) = false;
2959  }
2960  }
2961 
2962  boost::mutex mutex;
2963 
2964  RemapValuesThreadParams params;
2965  params.m_inputRasterPtr = &inputRaster;
2966  params.m_inputRasterBandsPtr = &inputRasterBands;
2967  params.m_outputRasterPtr = &outputRaster;
2968  params.m_outputRasterBandsPtr = &outputRasterBands;
2969  params.m_remapMatrixPtr = &remapMatrix;
2970  params.m_rasterBlocksStatusPtr = &rasterBlocksStatus;
2971  params.m_returnStatus = true;
2972  params.m_mutexPtr = &mutex;
2973 
2974  if( maxThreads == 1 )
2975  {
2976  RemapValuesThread( &params );
2977  }
2978  else
2979  {
2980  unsigned int threadsNumber = maxThreads ? maxThreads : te::common::GetPhysProcNumber();
2981 
2982  boost::thread_group threads;
2983 
2984  for( unsigned int threadIdx = 0 ; threadIdx < threadsNumber ;
2985  ++threadIdx )
2986  {
2987  threads.add_thread( new boost::thread( RemapValuesThread,
2988  &params ) );
2989  };
2990 
2991  threads.join_all();
2992  }
2993 
2994  return params.m_returnStatus;
2995  }
2996  else
2997  {
2998  const unsigned int inputRasterBandsSize = static_cast<unsigned int>(inputRasterBands.size());
2999  const unsigned int nRows = inputRaster.getNumberOfRows();
3000  const unsigned int nCols = inputRaster.getNumberOfColumns();
3001  boost::numeric::ublas::matrix< double > pixelValues( inputRasterBands.size(), 1 );
3002  boost::numeric::ublas::matrix< double > remappedPixelValues( inputRasterBands.size(), 1 );
3003  std::vector< double > inputNoDataValues( inputRasterBandsSize );
3004  std::vector< double > outputNoDataValues( inputRasterBandsSize );
3005  std::vector< double > outputMinValue( inputRasterBandsSize );
3006  std::vector< double > outputMaxValue( inputRasterBandsSize );
3007  unsigned int inputRasterBandsIdx = 0;
3008  unsigned int row = 0;
3009  unsigned int col = 0;
3010  bool pixelIsValid = false;
3011 
3012  for( inputRasterBandsIdx = 0 ; inputRasterBandsIdx < inputRasterBandsSize ;
3013  ++inputRasterBandsIdx )
3014  {
3015  if( inputRasterBands[ inputRasterBandsIdx ] >= inputRaster.getNumberOfBands() )
3016  {
3017  return false;
3018  }
3019 
3020  inputNoDataValues[ inputRasterBandsIdx ] = inputRaster.getBand(
3021  inputRasterBands[ inputRasterBandsIdx ] )->getProperty()->m_noDataValue;
3022 
3023  outputNoDataValues[ inputRasterBandsIdx ] = outputRaster.getBand(
3024  outputRasterBands[ inputRasterBandsIdx ] )->getProperty()->m_noDataValue;
3025 
3026  te::rst::GetDataTypeRanges(
3027  outputRaster.getBand( outputRasterBands[ inputRasterBandsIdx ] )->getProperty()->getType(),
3028  outputMinValue[ inputRasterBandsIdx ],
3029  outputMaxValue[ inputRasterBandsIdx ] );
3030  }
3031 
3032  for( row = 0 ; row < nRows ; ++row )
3033  {
3034  for( col = 0 ; col < nCols ; ++col )
3035  {
3036  pixelIsValid = true;
3037 
3038  for( inputRasterBandsIdx = 0 ; inputRasterBandsIdx < inputRasterBandsSize ;
3039  ++inputRasterBandsIdx )
3040  {
3041  inputRaster.getValue( col, row, pixelValues( inputRasterBandsIdx, 0 ),
3042  inputRasterBands[ inputRasterBandsIdx ] );
3043 
3044  if( pixelValues( inputRasterBandsIdx, 0 ) == inputNoDataValues[ inputRasterBandsIdx ] )
3045  {
3046  pixelIsValid = false;
3047  break;
3048  }
3049  }
3050 
3051  if( pixelIsValid )
3052  {
3053  remappedPixelValues = boost::numeric::ublas::prod( remapMatrix, pixelValues );
3054 
3055  for( inputRasterBandsIdx = 0 ; inputRasterBandsIdx < inputRasterBandsSize ;
3056  ++inputRasterBandsIdx )
3057  {
3058  outputRaster.setValue(
3059  col,
3060  row,
3061  std::max(
3062  outputMinValue[ inputRasterBandsIdx ]
3063  ,
3064  std::min(
3065  outputMaxValue[ inputRasterBandsIdx ]
3066  ,
3067  remappedPixelValues( inputRasterBandsIdx, 0 )
3068  )
3069  ),
3070  outputRasterBands[ inputRasterBandsIdx ] );
3071  }
3072  }
3073  else
3074  {
3075  for( inputRasterBandsIdx = 0 ; inputRasterBandsIdx < inputRasterBandsSize ;
3076  ++inputRasterBandsIdx )
3077  {
3078  outputRaster.setValue( col, row, outputNoDataValues[ inputRasterBandsIdx ],
3079  outputRasterBands[ inputRasterBandsIdx ] );
3080  }
3081  }
3082  }
3083  }
3084  }
3085 
3086  return true;
3087  }
3088 
3089  bool DecomposeBands(
3090  const te::rst::Raster& inputRaster,
3091  const std::vector< unsigned int >& inputRasterBands,
3092  const std::vector< std::map<std::string, std::string> > & outputRastersInfos,
3093  const std::string& outputDataSourceType,
3094  std::vector< boost::shared_ptr< te::rst::Raster > > & outputRastersPtrs )
3095  {
3096  outputRastersPtrs.clear();
3097 
3098  if( !( inputRaster.getAccessPolicy() & te::common::RAccess ) )
3099  {
3100  return false;
3101  }
3102  if( outputRastersInfos.size() != inputRasterBands.size() )
3103  {
3104  return false;
3105  }
3106  if( outputDataSourceType.empty() )
3107  {
3108  return false;
3109  }
3110 
3111  outputRastersPtrs.resize( inputRasterBands.size() );
3112 
3113  const unsigned int nRows = inputRaster.getNumberOfRows();
3114  const unsigned int nCols = inputRaster.getNumberOfColumns();
3115 
3116  for( unsigned int inputRasterBandsIdx = 0 ; inputRasterBandsIdx <
3117  inputRasterBands.size() ; ++inputRasterBandsIdx )
3118  {
3119  const unsigned int bandIdx = inputRasterBands[ inputRasterBandsIdx ];
3120  if( bandIdx >= inputRaster.getNumberOfBands() )
3121  {
3122  return false;
3123  }
3124 
3125  std::vector< te::rst::BandProperty* > bandsProperties;
3126 
3127  bandsProperties.push_back( new te::rst::BandProperty(
3128  *( inputRaster.getBand( bandIdx )->getProperty()) ) );
3129 
3130  outputRastersPtrs[ inputRasterBandsIdx].reset(
3131  te::rst::RasterFactory::make( outputDataSourceType,
3132  new te::rst::Grid( *inputRaster.getGrid() ), bandsProperties,
3133  outputRastersInfos[ inputRasterBandsIdx ], nullptr, nullptr ) );
3134  if( outputRastersPtrs[ inputRasterBandsIdx].get() == nullptr ) return false;
3135 
3136  unsigned int col = 0;
3137  unsigned int row = 0;
3138  double value = 0;
3139  const te::rst::Band& inBand = *(inputRaster.getBand( bandIdx ));
3140  te::rst::Band& outBand = *(outputRastersPtrs[ inputRasterBandsIdx]->getBand( 0 ));
3141 
3142  for( row = 0 ; row < nRows ; ++row )
3143  {
3144  for( col = 0 ; col < nCols ; ++col )
3145  {
3146  inBand.getValue( col, row, value );
3147  outBand.setValue( col, row, value );
3148  }
3149  }
3150  }
3151 
3152  return true;
3153  }
3154 
3155  bool ComposeBands(
3156  te::rp::FeederConstRaster& feeder,
3157  const std::vector< unsigned int >& inputRasterBands,
3158  const te::rst::Interpolator::Method& interpMethod,
3159  const std::map<std::string, std::string>& outputRasterInfo,
3160  const std::string& outputDataSourceType,
3161  std::unique_ptr< te::rst::Raster >& outputRasterPtr )
3162  {
3163  outputRasterPtr.reset();
3164 
3165  if( inputRasterBands.size() != feeder.getObjsCount() )
3166  {
3167  return false;
3168  }
3169  if( outputDataSourceType.empty() )
3170  {
3171  return false;
3172  }
3173 
3174  // creating the output raster
3175 
3176  {
3177  te::rst::Raster const * inputRasterPtr = nullptr;
3178  std::unique_ptr< te::rst::Grid > outputGridPtr;
3179  std::vector< te::rst::BandProperty* > bandsProperties;
3180 
3181  feeder.reset();
3182  while( (inputRasterPtr = feeder.getCurrentObj()) )
3183  {
3184  if( inputRasterBands[ feeder.getCurrentOffset() ] >=
3185  inputRasterPtr->getNumberOfBands() )
3186  {
3187  return false;
3188  }
3189 
3190  if( feeder.getCurrentOffset() == 0 )
3191  {
3192  outputGridPtr.reset( new te::rst::Grid( *inputRasterPtr->getGrid() ) );
3193  }
3194 
3195  bandsProperties.push_back( new te::rst::BandProperty(
3196  *( inputRasterPtr->getBand( inputRasterBands[
3197  feeder.getCurrentOffset() ] )->getProperty()) ) );
3198 
3199  feeder.moveNext();
3200  }
3201 
3202  outputRasterPtr.reset( te::rst::RasterFactory::make( outputDataSourceType,
3203  outputGridPtr.release(), bandsProperties, outputRasterInfo, nullptr, nullptr ) );
3204  if( outputRasterPtr.get() == nullptr ) return false;
3205  }
3206 
3207  // copying data from each input band
3208 
3209  {
3210  te::rst::Raster const * inputRasterPtr = nullptr;
3211 
3212  feeder.reset();
3213  while( (inputRasterPtr = feeder.getCurrentObj()) )
3214  {
3215  const unsigned int inBandIdx = inputRasterBands[
3216  feeder.getCurrentOffset() ];
3217  te::rst::Interpolator interp( inputRasterPtr, interpMethod );
3218  unsigned int outRow = 0;
3219  unsigned int outCol = 0;
3220  const unsigned int nOutRows = outputRasterPtr->getNumberOfRows();
3221  const unsigned int nOutCols = outputRasterPtr->getNumberOfColumns();
3222  te::rst::Band& outBand = *outputRasterPtr->getBand(
3223  feeder.getCurrentOffset() );
3224  const te::rst::Grid& inGrid = *inputRasterPtr->getGrid();
3225  const te::rst::Grid& outGrid = *outputRasterPtr->getGrid();
3226  double xOutCoord = 0;
3227  double yOutCoord = 0;
3228  double xInCoord = 0;
3229  double yInCoord = 0;
3230  double inRow = 0;
3231  double inCol = 0;
3232  std::complex< double > value = 0;
3233  te::srs::Converter conv( outputRasterPtr->getSRID(),
3234  inputRasterPtr->getSRID() );
3235 
3236  if( inputRasterPtr->getSRID() == outputRasterPtr->getSRID() )
3237  {
3238  for( outRow = 0 ; outRow < nOutRows ; ++outRow )
3239  {
3240  for( outCol = 0 ; outCol < nOutCols ; ++outCol )
3241  {
3242  outGrid.gridToGeo( (double)outCol, (double)outRow, xOutCoord, yOutCoord );
3243  inGrid.geoToGrid( xOutCoord, yOutCoord, inCol, inRow );
3244  interp.getValue( inCol, inRow, value, inBandIdx );
3245  outBand.setValue( outCol, outRow, value );
3246  }
3247  }
3248  }
3249  else
3250  {
3251  for( outRow = 0 ; outRow < nOutRows ; ++outRow )
3252  {
3253  for( outCol = 0 ; outCol < nOutCols ; ++outCol )
3254  {
3255  outGrid.gridToGeo( (double)outCol, (double)outRow, xOutCoord, yOutCoord );
3256  conv.convert( xOutCoord, yOutCoord, xInCoord, yInCoord );
3257  inGrid.geoToGrid( xInCoord, yInCoord, inCol, inRow );
3258  interp.getValue( inCol, inRow, value, inBandIdx );
3259  outBand.setValue( outCol, outRow, value );
3260  }
3261  }
3262  }
3263 
3264  feeder.moveNext();
3265  }
3266  }
3267 
3268  return true;
3269  }
3270 
3271  bool GetDetailedExtent( const te::rst::Grid& grid, te::gm::LinearRing& detailedExtent )
3272  {
3273  const unsigned int nCols = grid.getNumberOfColumns();
3274  const unsigned int nRows = grid.getNumberOfRows();
3275  if( ( nCols == 0 ) || ( nRows == 0 ) )
3276  {
3277  return false;
3278  }
3279 
3280  const unsigned int ringSize = ( 2 * ( nCols + 1 ) ) +
3281  ( 2 * ( nRows - 1 ) ) + 1;
3282 
3283  te::gm::LinearRing ring( ringSize , te::gm::LineStringType,
3284  grid.getSRID(), nullptr );
3285 
3286  const double lLX = grid.getExtent()->getLowerLeftX();
3287  const double lLY = grid.getExtent()->getLowerLeftY();
3288  const double uRX = grid.getExtent()->getUpperRightX();
3289  const double uRY = grid.getExtent()->getUpperRightY();
3290  const double& resX = grid.getResolutionX();
3291  const double& resY = grid.getResolutionY();
3292  unsigned int ringIdx = 0;
3293  unsigned int row = 0;
3294  unsigned int col = 0;
3295 
3296  ring.setPoint( 0, lLX, uRY );
3297 
3298  for( col = 1 ; col <= nCols ; ++col )
3299  {
3300  ring.setPoint( ++ringIdx, lLX + ( (static_cast<double>( col ) ) * resX ), uRY );
3301  }
3302 
3303  for( row = 1 ; row <= nRows ; ++row )
3304  {
3305  ring.setPoint( ++ringIdx, uRX, uRY - ( (static_cast<double>( row ) ) * resY ) );
3306  }
3307 
3308  for( col = 1 ; col <= nCols ; ++col )
3309  {
3310  ring.setPoint( ++ringIdx, uRX - ( (static_cast<double>( col ) ) * resX ), lLY );
3311  }
3312 
3313  for( row = 1 ; row <= nRows ; ++row )
3314  {
3315  ring.setPoint( ++ringIdx, lLX, lLY + ( (static_cast<double>( row ) ) * resY ) );
3316  }
3317 
3318  // Ensuring correct closing
3319  ring.setPoint( ringSize - 1, lLX, uRY );
3320 
3321  detailedExtent = ring;
3322 
3323  return true;
3324  }
3325 
3326  bool GetIndexedDetailedExtent( const te::rst::Grid& grid,
3327  te::gm::LinearRing& indexedDetailedExtent )
3328  {
3329  const unsigned int nCols = grid.getNumberOfColumns();
3330  const unsigned int nRows = grid.getNumberOfRows();
3331  if( ( nCols == 0 ) || ( nRows == 0 ) )
3332  {
3333  return false;
3334  }
3335 
3336  const unsigned int ringSize = ( 2 * ( nCols + 1 ) ) +
3337  ( 2 * ( nRows - 1 ) ) + 1;
3338 
3339  te::gm::LinearRing ring( ringSize , te::gm::LineStringType,
3340  0, nullptr );
3341 
3342  const double lLY = (static_cast<double>(nRows)) - 0.5;
3343  const double uRX = (static_cast<double>(nCols)) - 0.5;
3344  unsigned int ringIdx = 0;
3345  unsigned int row = 0;
3346  unsigned int col = 0;
3347 
3348  ring.setPoint( 0, -0.5, -0.5 );
3349 
3350  for( col = 0 ; col < nCols ; ++col )
3351  {
3352  ring.setPoint( ++ringIdx, 0.5 + (static_cast<double>(col)), (-0.5) );
3353  }
3354 
3355  for( row = 0 ; row < nRows ; ++row )
3356  {
3357  ring.setPoint( ++ringIdx, uRX, 0.5 + (static_cast<double>(row)) );
3358  }
3359 
3360  for( col = 1; col <= nCols ; ++col )
3361  {
3362  ring.setPoint( ++ringIdx, uRX - (static_cast<double>(col)), lLY );
3363  }
3364 
3365  for( row = 1 ; row <= nRows ; ++row )
3366  {
3367  ring.setPoint( ++ringIdx, (-0.5), lLY - (static_cast<double>(row)) );
3368  }
3369 
3370  // Ensuring correct closing
3371  ring.setPoint( ringSize - 1, -0.5, -0.5 );
3372 
3373  indexedDetailedExtent = ring;
3374 
3375  return true;
3376  }
3377 
3378  boost::numeric::ublas::matrix< double > CreateWaveletAtrousFilter(
3379  const WaveletAtrousFilterType& filterType )
3380  {
3381  boost::numeric::ublas::matrix< double > emptyFilter;
3382 
3383  switch( filterType )
3384  {
3385  case B3SplineWAFilter :
3386  {
3387  boost::numeric::ublas::matrix< double > internalFilter( 5, 5 );
3388  const double weight = 256;
3389 
3390  internalFilter(0,0) = 1/weight; internalFilter(0,1) = 4/weight; internalFilter(0,2) = 6/weight; internalFilter(0,3) = 4/weight; internalFilter(0,4) = 1/weight;
3391  internalFilter(1,0) = 4/weight; internalFilter(1,1) = 16/weight; internalFilter(1,2) = 24/weight; internalFilter(1,3) = 16/weight; internalFilter(1,4) = 4/weight;
3392  internalFilter(2,0) = 6/weight; internalFilter(2,1) = 24/weight; internalFilter(2,2) = 36/weight; internalFilter(2,3) = 24/weight; internalFilter(2,4) = 6/weight;
3393  internalFilter(3,0) = 4/weight; internalFilter(3,1) = 16/weight; internalFilter(3,2) = 24/weight; internalFilter(3,3) = 16/weight; internalFilter(3,4) = 4/weight;
3394  internalFilter(4,0) = 1/weight; internalFilter(4,1) = 4/weight; internalFilter(4,2) = 6/weight; internalFilter(4,3) = 4/weight; internalFilter(4,4) = 1/weight;
3395 
3396  return internalFilter;
3397  }
3398  case TriangleWAFilter :
3399  {
3400  boost::numeric::ublas::matrix< double > internalFilter( 3, 3 );
3401  const double weight = 16;
3402 
3403  internalFilter(0,0) = 1/weight; internalFilter(0,1) = 2/weight; internalFilter(0,2) = 1/weight;
3404  internalFilter(1,0) = 2/weight; internalFilter(1,1) = 4/weight; internalFilter(1,2) = 2/weight;
3405  internalFilter(2,0) = 1/weight; internalFilter(2,1) = 2/weight; internalFilter(2,2) = 1/weight;
3406 
3407  return internalFilter;
3408  }
3409  default :
3410  {
3411  throw te::rp::Exception( "Invalid filter type" );
3412  }
3413  }
3414 
3415  return emptyFilter;
3416  }
3417 
3418  bool DirectWaveletAtrous(
3419  const te::rst::Raster& inputRaster,
3420  const std::vector< unsigned int >& inputRasterBands,
3421  te::rst::Raster& waveletRaster,
3422  const unsigned int levelsNumber,
3423  const boost::numeric::ublas::matrix< double >& filter )
3424  {
3425  if( ! ( inputRaster.getAccessPolicy() & te::common::RAccess ) )
3426  {
3427  return false;
3428  }
3429  for( unsigned int inputRasterBandsIdx = 0 ; inputRasterBandsIdx <
3430  inputRasterBands.size() ; ++inputRasterBandsIdx )
3431  {
3432  if( inputRasterBands[ inputRasterBandsIdx ] >= inputRaster.getNumberOfBands() )
3433  {
3434  return false;
3435  }
3436  }
3437  if(
3438  ( ( waveletRaster.getAccessPolicy() & te::common::WAccess ) == 0 )
3439  ||
3440  ( waveletRaster.getNumberOfColumns() != inputRaster.getNumberOfColumns() )
3441  ||
3442  ( waveletRaster.getNumberOfRows() != inputRaster.getNumberOfRows() )
3443  ||
3444  ( waveletRaster.getNumberOfBands() < ( 2 * levelsNumber * inputRasterBands.size() ) )
3445  )
3446  {
3447  return false;
3448  }
3449  if( levelsNumber == 0 )
3450  {
3451  return false;
3452  }
3453  if( ( filter.size1() == 0 ) || ( filter.size2() == 0 ) ||
3454  ( filter.size1() != filter.size2() )
3455  )
3456  {
3457  return false;
3458  }
3459 
3460  // Creating the coeficients and the resitual planes for each required
3461  // wavelet level
3462 
3463  const int filterWidth = static_cast<int>(filter.size1());
3464  const int offset = filterWidth / 2;
3465  const int nLines = static_cast<int>(inputRaster.getNumberOfRows());
3466  const int nCols = static_cast<int>(inputRaster.getNumberOfColumns());
3467 
3468  for( unsigned int inputRasterBandsIdx = 0 ; inputRasterBandsIdx < inputRasterBands.size() ;
3469  ++inputRasterBandsIdx )
3470  {
3471  for( unsigned int levelIndex = 0; levelIndex < levelsNumber; ++levelIndex)
3472  {
3473  const unsigned int currentSmoothBandIdx = ( 2 * levelsNumber *
3474  inputRasterBandsIdx ) + ( 2 * levelIndex );
3475  te::rst::Band& currentSmoothBand = *waveletRaster.getBand( currentSmoothBandIdx );
3476 
3477  const unsigned int currentWaveletBandIdx = currentSmoothBandIdx + 1;
3478  te::rst::Band& currentWaveletBand = *waveletRaster.getBand(
3479  currentWaveletBandIdx );
3480 
3481  const unsigned int prevSmoothBandIdx = currentSmoothBandIdx - 2;
3482  te::rst::Band const* prevSmoothBandPtr = nullptr;
3483  if( levelIndex == 0 )
3484  {
3485  prevSmoothBandPtr = inputRaster.getBand( inputRasterBands[ inputRasterBandsIdx ] );
3486  }
3487  else
3488  {
3489  prevSmoothBandPtr = waveletRaster.getBand( prevSmoothBandIdx );
3490  }
3491  const te::rst::Band& prevSmoothBand = *prevSmoothBandPtr;
3492 
3493  const int filterScale = static_cast<int>(std::pow(2.0, static_cast<double>(levelIndex)));
3494 
3495  int col = 0;
3496  int row = 0;
3497  int convolutionCenterCol = 0;
3498  int convolutionCenterRow = 0;
3499  int filterCol = 0;
3500  int filterRow = 0;
3501  double valueOriginal = 0.0;
3502  double valueNew = 0.0;
3503 
3504  for (convolutionCenterRow = 0; convolutionCenterRow < nLines; convolutionCenterRow++)
3505  {
3506  for (convolutionCenterCol = 0; convolutionCenterCol < nCols; convolutionCenterCol++)
3507  {
3508  valueNew = 0.0;
3509 
3510  for (filterRow = 0; filterRow < filterWidth; filterRow++)
3511  {
3512  for (filterCol = 0; filterCol < filterWidth; filterCol++)
3513  {
3514  col = convolutionCenterCol+(filterCol-offset)*filterScale;
3515  row = convolutionCenterRow+(filterRow-offset)*filterScale;
3516 
3517  if (col < 0)
3518  col = nCols + col;
3519  else if (col >= nCols)
3520  col = col - nCols;
3521  if (row < 0)
3522  row = nLines + row;
3523  else if (row >= nLines)
3524  row = row - nLines;
3525 
3526  prevSmoothBand.getValue( static_cast<unsigned int>(col), static_cast<unsigned int>(row), valueOriginal );
3527 
3528  valueNew += valueOriginal * filter( static_cast<size_t>(filterRow), static_cast<size_t>(filterCol) );
3529  }
3530  }
3531 
3532  currentSmoothBand.setValue( static_cast<unsigned int>(convolutionCenterCol), static_cast<unsigned int>(convolutionCenterRow),
3533  valueNew );
3534  prevSmoothBand.getValue( static_cast<unsigned int>(convolutionCenterCol), static_cast<unsigned int>(convolutionCenterRow),
3535  valueOriginal );
3536  currentWaveletBand.setValue( static_cast<unsigned int>(convolutionCenterCol), static_cast<unsigned int>(convolutionCenterRow),
3537  valueOriginal - valueNew );
3538  }
3539  }
3540  }
3541  }
3542 
3543  return true;
3544  }
3545 
3546  bool InverseWaveletAtrous(
3547  const te::rst::Raster& waveletRaster,
3548  const unsigned int levelsNumber,
3549  te::rst::Raster& outputRaster,
3550  const std::vector< unsigned int >& outputRasterBands )
3551  {
3552  if(
3553  ( ! ( waveletRaster.getAccessPolicy() & te::common::RAccess ) )
3554  ||
3555  ( waveletRaster.getNumberOfBands() < ( 2 * levelsNumber * outputRasterBands.size() ) )
3556  )
3557  {
3558  return false;
3559  }
3560  if( levelsNumber == 0 )
3561  {
3562  return false;
3563  }
3564  if(
3565  ( ( outputRaster.getAccessPolicy() & te::common::WAccess ) == 0 )
3566  ||
3567  ( waveletRaster.getNumberOfColumns() != outputRaster.getNumberOfColumns() )
3568  ||
3569  ( waveletRaster.getNumberOfRows() != outputRaster.getNumberOfRows() )
3570  )
3571  {
3572  return false;
3573  }
3574  for( unsigned int outputRasterBandsIdx = 0 ; outputRasterBandsIdx <
3575  outputRasterBands.size() ; ++outputRasterBandsIdx )
3576  {
3577  if( outputRasterBands[ outputRasterBandsIdx ] >= outputRaster.getNumberOfBands() )
3578  {
3579  return false;
3580  }
3581  }
3582 
3583  for( unsigned int outputRasterBandsIdx = 0 ; outputRasterBandsIdx <
3584  outputRasterBands.size() ; ++outputRasterBandsIdx )
3585  {
3586  const unsigned int outputRasterBandIdx = outputRasterBands[ outputRasterBandsIdx ];
3587  const unsigned int nRows = waveletRaster.getNumberOfRows();
3588  const unsigned int nCols = waveletRaster.getNumberOfColumns();
3589  const unsigned int firstSmoothBandIdx = outputRasterBandsIdx *
3590  levelsNumber * 2;
3591  const unsigned int firstWaveletBandIdx = firstSmoothBandIdx
3592  + 1;
3593  const unsigned int lastSmoothBandIdx = firstSmoothBandIdx
3594  + ( 2 * levelsNumber ) - 2;
3595  const unsigned int lastWaveletBandIdx = lastSmoothBandIdx + 1;
3596  unsigned int col = 0;
3597  unsigned int row = 0;
3598  double value = 0;
3599  double sum = 0;
3600  unsigned int waveletRasterBand = 0;
3601 
3602  double bandAllowedMinValue = 0;
3603  double bandAllowedMaxValue = 0;
3604  te::rst::GetDataTypeRanges( outputRaster.getBand(
3605  outputRasterBandIdx )->getProperty()->m_type, bandAllowedMinValue,
3606  bandAllowedMaxValue );
3607 
3608 
3609  for( row = 0 ; row < nRows ; ++row )
3610  {
3611  for( col = 0 ; col < nCols ; ++col )
3612  {
3613  sum = 0.0;
3614 
3615  for( waveletRasterBand = firstWaveletBandIdx ;
3616  waveletRasterBand <= lastWaveletBandIdx ;
3617  waveletRasterBand += 2)
3618  {
3619  waveletRaster.getValue( col, row, value, waveletRasterBand );
3620  sum += value;
3621  }
3622 
3623  waveletRaster.getValue( col, row, value, lastSmoothBandIdx );
3624  sum += value;
3625 
3626  sum = std::max( sum, bandAllowedMinValue );
3627  sum = std::min( sum, bandAllowedMaxValue );
3628 
3629  outputRaster.setValue( col, row, sum, outputRasterBandIdx );
3630  }
3631  }
3632  }
3633 
3634  return true;
3635  }
3636 
3637  bool RasterResample(
3638  const te::rst::Raster& inputRaster,
3639  const std::vector< unsigned int >& inputRasterBands,
3640  const te::rst::Interpolator::Method interpMethod,
3641  const unsigned int firstRow,
3642  const unsigned int firstColumn,
3643  const unsigned int height,
3644  const unsigned int width,
3645  const unsigned int newheight,
3646  const unsigned int newwidth,
3647  const std::map<std::string, std::string>& rinfo,
3648  const std::string& dataSourceType,
3649  std::unique_ptr< te::rst::Raster >& resampledRasterPtr )
3650  {
3651  if( ( firstRow + height ) > inputRaster.getNumberOfRows() )
3652  {
3653  return false;
3654  }
3655  if( ( firstColumn + width ) > inputRaster.getNumberOfColumns() )
3656  {
3657  return false;
3658  }
3659  for (std::size_t inputRasterBandsIdx = 0; inputRasterBandsIdx <
3660  inputRasterBands.size(); inputRasterBandsIdx++)
3661  {
3662  if( inputRasterBands[ inputRasterBandsIdx ] >= inputRaster.getNumberOfBands() )
3663  {
3664  return false;
3665  }
3666  }
3667 
3668  te::gm::Coord2D ulc = inputRaster.getGrid()->gridToGeo( (static_cast<double>(firstColumn))
3669  - 0.5, (static_cast<double>(firstRow)) - 0.5);
3670  te::gm::Coord2D lrc = inputRaster.getGrid()->gridToGeo( (static_cast<double>(firstColumn
3671  + width)) - 0.5, (static_cast<double>(firstRow + height)) - 0.5);
3672 
3673  std::unique_ptr< te::gm::Envelope > newEnvelopePtr( new te::gm::Envelope(
3674  ulc.x, lrc.y, lrc.x, ulc.y ) );
3675 
3676  // create output parameters and raster
3677 
3678  std::unique_ptr< te::rst::Grid > gridPtr( new te::rst::Grid(
3679  newwidth, newheight, newEnvelopePtr.release(),
3680  inputRaster.getSRID()) );
3681 
3682  std::vector<te::rst::BandProperty*> bands;
3683 
3684  for (std::size_t inputRasterBandsIdx = 0; inputRasterBandsIdx <
3685  inputRasterBands.size(); inputRasterBandsIdx++)
3686  {
3687  bands.push_back( new te::rst::BandProperty(
3688  *(inputRaster.getBand( inputRasterBands[ inputRasterBandsIdx ] )->getProperty())));
3689  bands[ inputRasterBandsIdx ]->m_blkh = 1;
3690  bands[ inputRasterBandsIdx ]->m_blkw = static_cast<int>(newwidth);
3691  bands[ inputRasterBandsIdx ]->m_nblocksx = 1;
3692  bands[ inputRasterBandsIdx ]->m_nblocksy = static_cast<int>(newheight);
3693  }
3694 
3695  resampledRasterPtr.reset( te::rst::RasterFactory::make( dataSourceType,
3696  gridPtr.release(), bands, rinfo, nullptr ) );
3697  if( resampledRasterPtr.get() == nullptr )
3698  {
3699  return false;
3700  }
3701 
3702  // fill output raster
3703 
3704  std::complex<double> interpValue;
3705  te::rst::Interpolator interp( &inputRaster, interpMethod);
3706  const double rowsFactor = (static_cast<double>(height-1)) / (static_cast<double>(newheight-1));
3707  const double colsFactor = (static_cast<double>(width-1)) / (static_cast<double>(newwidth-1));
3708  double inputRow = 0;
3709  double inputCol = 0;
3710  unsigned int outputCol = 0;
3711  unsigned int outputRow = 0;
3712  unsigned int inputBandIdx = 0;
3713  double allowedMin = 0;
3714  double allowedMax = 0;
3715 
3716  for (std::size_t inputRasterBandsIdx = 0; inputRasterBandsIdx <
3717  inputRasterBands.size(); inputRasterBandsIdx++)
3718  {
3719  te::rst::Band& outputBand = *resampledRasterPtr->getBand( inputRasterBandsIdx );
3720  inputBandIdx = inputRasterBands[ inputRasterBandsIdx ];
3721  te::rp::GetDataTypeRange( outputBand.getProperty()->m_type, allowedMin,
3722  allowedMax );
3723 
3724  for ( outputRow = 0; outputRow < newheight; ++outputRow)
3725  {
3726  inputRow = ( (static_cast<double>( outputRow ) ) * rowsFactor ) + (static_cast<double>(firstRow));
3727 
3728  for ( outputCol = 0; outputCol < newwidth; ++outputCol )
3729  {
3730  inputCol = ( (static_cast<double>( outputCol ) ) * colsFactor ) + (static_cast<double>(firstColumn));
3731 
3732  interp.getValue(inputCol, inputRow, interpValue, inputBandIdx);
3733 
3734  interpValue.real( std::max( allowedMin, interpValue.real() ) );
3735  interpValue.real( std::min( allowedMax, interpValue.real() ) );
3736 
3737  outputBand.setValue(outputCol, outputRow, interpValue);
3738  }
3739  }
3740  }
3741 
3742  return true;
3743  }
3744 
3745  void CreateFixedStepPalette(
3746  const unsigned int paletteSize, const bool randomize,
3747  std::vector< te::rst::BandProperty::ColorEntry >& palette )
3748  {
3749  palette.resize( paletteSize );
3750  if( paletteSize == 0 ) return;
3751 
3752  const short int step = static_cast<short int>( 256.0 /
3753  std::cbrt( static_cast<double>(paletteSize) ) );
3754  short int green = 0;
3755  short int blue = 0;
3756  unsigned int paletteIdx = 0;
3757 
3758  for( short int red = 0 ; red < 256 ; red += step )
3759  {
3760  for( green = 0 ; green < 256 ; green += step )
3761  {
3762  for( blue = 0 ; blue < 256 ; blue += step )
3763  {
3764  if( paletteIdx >= paletteSize )
3765  {
3766  green = 255;
3767  red = 255;
3768  break;
3769  }
3770 
3771  palette[ paletteIdx ].c1 = red;
3772  palette[ paletteIdx ].c2 = green;
3773  palette[ paletteIdx ].c3 = blue;
3774  palette[ paletteIdx ].c4 = 255;
3775 
3776  ++paletteIdx;
3777  }
3778  }
3779  }
3780 
3781  if( randomize )
3782  {
3783  std::random_shuffle( palette.begin(), palette.end() );
3784  }
3785  }
3786 
3787  bool RasterSlicing(
3788  const te::rst::Raster& inputRaster,
3789  const unsigned int inputRasterBand,
3790  const bool createPaletteRaster,
3791  const unsigned int slicesNumber,
3792  const bool eqHistogram,
3793  const std::map< std::string, std::string >& rasterInfo,
3794  const std::string& rasterType,
3795  const bool enableProgress,
3796  std::vector< te::rst::BandProperty::ColorEntry > const * const palettePtr,
3797  std::unique_ptr< te::rst::Raster >& outRasterPtr ,
3798  std::map< double, double > limits)
3799  {
3800  TERP_TRUE_OR_RETURN_FALSE( inputRasterBand < inputRaster.getNumberOfBands(),
3801  "Invalid raster band index" );
3802 
3803  te::rst::Raster const * internalInRasterPtr = &inputRaster;
3804  unsigned int internalInRasterBand = inputRasterBand;
3805  std::unique_ptr< te::rst::Raster > eqRasterPtr;
3806  if( eqHistogram )
3807  {
3808  double min = 0;
3809  double max = 0;
3810  GetDataTypeRange( inputRaster.getBand( inputRasterBand )->getProperty()->m_type,
3811  min, max );
3812 
3813  Contrast::InputParameters contInputPars;
3814  contInputPars.m_type = Contrast::InputParameters::HistogramEqualizationContrastT;
3815  contInputPars.m_inRasterPtr = &inputRaster;
3816  contInputPars.m_inRasterBands.push_back( inputRasterBand );
3817  contInputPars.m_hECMaxInput.push_back( max );
3818  contInputPars.m_enableProgress = enableProgress;
3819 
3820  Contrast contInstance;
3821  TERP_TRUE_OR_RETURN_FALSE( contInstance.initialize( contInputPars ),
3822  "Contrast init error" );
3823 
3824  Contrast::OutputParameters contOutPars;
3825  contOutPars.m_createdOutRasterDSType = "EXPANSIBLE";
3826  TERP_TRUE_OR_RETURN_FALSE( contInstance.execute( contOutPars ),
3827  "Contrast exec error" );
3828 
3829  eqRasterPtr.reset( contOutPars.m_createdOutRasterPtr.release() );
3830  internalInRasterPtr = eqRasterPtr.get();
3831  internalInRasterBand = 0;
3832  }
3833 
3834  std::vector< te::rst::BandProperty::ColorEntry > const * internalPalettePtr =
3835  palettePtr;
3836  std::vector< te::rst::BandProperty::ColorEntry > internalPalette;
3837  if( palettePtr == nullptr )
3838  {
3839  CreateFixedStepPalette( slicesNumber, true, internalPalette );
3840  internalPalettePtr = &internalPalette;
3841  }
3842 
3843  std::unique_ptr< te::common::TaskProgress > progressPtr;
3844  if( enableProgress )
3845  {
3846  progressPtr.reset( new te::common::TaskProgress );
3847  progressPtr->setMessage( "Slicing" );
3848  progressPtr->setTotalSteps( static_cast<int>(2 * internalInRasterPtr->getNumberOfRows()) );
3849  }
3850 
3851  if( createPaletteRaster )
3852  {
3853  std::vector< te::rst::BandProperty* > bandsProperties;
3854  bandsProperties.push_back( new te::rst::BandProperty(
3855  * internalInRasterPtr->getBand( internalInRasterBand )->getProperty() ) );
3856  bandsProperties[ 0 ]->m_colorInterp = te::rst::PaletteIdxCInt;
3857  bandsProperties[ 0 ]->m_paletteInterp = te::rst::RGBPalInt;
3858  bandsProperties[ 0 ]->m_palette = (*internalPalettePtr);
3859  bandsProperties[ 0 ]->m_noDataValue = static_cast<double>(internalPalettePtr->size());
3860  if( internalPalettePtr->size() <= std::numeric_limits< unsigned char >::max() )
3861  {
3862  bandsProperties[ 0 ]->m_type = te::dt::UCHAR_TYPE;
3863  }
3864  else if( internalPalettePtr->size() <= std::numeric_limits< unsigned short int >::max() )
3865  {
3866  bandsProperties[ 0 ]->m_type = te::dt::UINT16_TYPE;
3867  }
3868  else if( internalPalettePtr->size() <= std::numeric_limits< unsigned int >::max() )
3869  {
3870  bandsProperties[ 0 ]->m_type = te::dt::UINT32_TYPE;
3871  }
3872  else
3873  {
3874  TERP_LOG_AND_RETURN_FALSE( "Invalid palette size" );
3875  }
3876 
3877  TERP_TRUE_OR_RETURN_FALSE( CreateNewRaster( *internalInRasterPtr->getGrid(),
3878  bandsProperties, rasterInfo, rasterType, outRasterPtr ),
3879  "Output raster creation error" );
3880 
3881  const te::rst::Band& inBand = (*internalInRasterPtr->getBand( internalInRasterBand ));
3882  te::rst::Band& outBand = (*outRasterPtr->getBand( 0 ));
3883  const unsigned int nRows = internalInRasterPtr->getNumberOfRows();
3884  const unsigned int nCols = internalInRasterPtr->getNumberOfColumns();
3885  const double inputNoDataValue = inBand.getProperty()->m_noDataValue;
3886  const double outputNoDataValue = outBand.getProperty()->m_noDataValue;
3887  unsigned int col = 0;
3888  double value = 0;
3889 
3890  for( unsigned int row = 0 ; row < nRows ; ++row )
3891  {
3892  for( col = 0 ; col < nCols ; ++col )
3893  {
3894  inBand.getValue( col, row, value );
3895 
3896  if( value == inputNoDataValue )
3897  {
3898  outBand.setValue( col, row, outputNoDataValue );
3899  }
3900  else
3901  {
3902  size_t index = 0;
3903  for (std::map<double, double>::iterator it = limits.begin(); it != limits.end(); ++it, ++index)
3904  {
3905  if (value >= (*it).first && value < (*it).second)
3906  {
3907  outBand.setValue(col, row, (double)index);
3908  break;
3909  }
3910  }
3911  if (index == limits.size())
3912  outBand.setValue(col, row, outputNoDataValue);
3913  }
3914  }
3915 
3916  if( enableProgress )
3917  {
3918  progressPtr->pulse();
3919  if( ! progressPtr->isActive() ) return false;
3920  }
3921  }
3922  }
3923  else
3924  {
3925  short internalPaletteMaxValue = 0;
3926  for( unsigned int internalPaletteIdx = 0 ; internalPaletteIdx <
3927  internalPalettePtr->size() ; ++internalPaletteIdx )
3928  {
3929  if( internalPalettePtr->operator[]( internalPaletteIdx ).c1 >
3930  internalPaletteMaxValue )
3931  {
3932  internalPaletteMaxValue = internalPalettePtr->operator[]( internalPaletteIdx ).c1;
3933  }
3934  if( internalPalettePtr->operator[]( internalPaletteIdx ).c2 >
3935  internalPaletteMaxValue )
3936  {
3937  internalPaletteMaxValue = internalPalettePtr->operator[]( internalPaletteIdx ).c2;
3938  }
3939  if( internalPalettePtr->operator[]( internalPaletteIdx ).c3 >
3940  internalPaletteMaxValue )
3941  {
3942  internalPaletteMaxValue = internalPalettePtr->operator[]( internalPaletteIdx ).c3;
3943  }
3944  if( internalPalettePtr->operator[]( internalPaletteIdx ).c4 >
3945  internalPaletteMaxValue )
3946  {
3947  internalPaletteMaxValue = internalPalettePtr->operator[]( internalPaletteIdx ).c4;
3948  }
3949  }
3950 
3951  int outDataType = 0;
3952  if( internalPaletteMaxValue < ( std::numeric_limits< unsigned char >::max() ) )
3953  {
3954  outDataType = te::dt::UCHAR_TYPE;
3955  }
3956  else if( internalPaletteMaxValue < ( std::numeric_limits< unsigned short int >::max() ) )
3957  {
3958  outDataType = te::dt::UINT16_TYPE;
3959  }
3960  else
3961  {
3962  outDataType = te::dt::UINT32_TYPE;
3963  }
3964 
3965  std::vector< te::rst::BandProperty* > bandsProperties;
3966  bandsProperties.push_back( new te::rst::BandProperty(
3967  * internalInRasterPtr->getBand( internalInRasterBand )->getProperty() ) );
3968  bandsProperties.push_back( new te::rst::BandProperty(
3969  * internalInRasterPtr->getBand( internalInRasterBand )->getProperty() ) );
3970  bandsProperties.push_back( new te::rst::BandProperty(
3971  * internalInRasterPtr->getBand( internalInRasterBand )->getProperty() ) );
3972  bandsProperties[ 0 ]->m_colorInterp = te::rst::RedCInt;
3973  bandsProperties[ 1 ]->m_colorInterp = te::rst::GreenCInt;
3974  bandsProperties[ 2 ]->m_colorInterp = te::rst::BlueCInt;
3975  bandsProperties[ 0 ]->m_type = outDataType;
3976  bandsProperties[ 1 ]->m_type = outDataType;
3977  bandsProperties[ 2 ]->m_type = outDataType;
3978  bandsProperties[ 0 ]->m_noDataValue = internalPaletteMaxValue + 1;
3979  bandsProperties[ 1 ]->m_noDataValue = internalPaletteMaxValue + 1;
3980  bandsProperties[ 2 ]->m_noDataValue = internalPaletteMaxValue + 1;
3981 
3982  TERP_TRUE_OR_RETURN_FALSE( CreateNewRaster( *internalInRasterPtr->getGrid(),
3983  bandsProperties, rasterInfo, rasterType, outRasterPtr ),
3984  "Output raster creation error" );
3985 
3986  const te::rst::Band& inBand = (*internalInRasterPtr->getBand( internalInRasterBand ));
3987  te::rst::Band& outBandRed = (*outRasterPtr->getBand( 0 ));
3988  te::rst::Band& outBandGreen = (*outRasterPtr->getBand( 1 ));
3989  te::rst::Band& outBandBlue = (*outRasterPtr->getBand( 2 ));
3990  const double rBandNoDataValue = outBandRed.getProperty()->m_noDataValue;
3991  const double gBandNoDataValue = outBandGreen.getProperty()->m_noDataValue;
3992  const double bBandNoDataValue = outBandBlue.getProperty()->m_noDataValue;
3993  const unsigned int nRows = internalInRasterPtr->getNumberOfRows();
3994  const unsigned int nCols = internalInRasterPtr->getNumberOfColumns();
3995  const double inputNoDataValue = inBand.getProperty()->m_noDataValue;
3996  unsigned int col = 0;
3997  double value = 0;
3998 
3999  for( unsigned int row = 0 ; row < nRows ; ++row )
4000  {
4001  for( col = 0 ; col < nCols ; ++col )
4002  {
4003  inBand.getValue( col, row, value );
4004 
4005  if( value == inputNoDataValue )
4006  {
4007  outBandRed.setValue( col, row, rBandNoDataValue );
4008  outBandGreen.setValue( col, row, gBandNoDataValue );
4009  outBandBlue.setValue( col, row, bBandNoDataValue );
4010  }
4011  else
4012  {
4013  size_t index = 0;
4014  for (std::map<double, double>::iterator it = limits.begin(); it != limits.end(); ++it, ++index)
4015  {
4016  if (value >= (*it).first && value < (*it).second)
4017  break;
4018  }
4019 
4020  if (index < limits.size())
4021  {
4022  const te::rst::BandProperty::ColorEntry& cEntry = internalPalettePtr->operator[](static_cast<unsigned int>(index));
4023 
4024  outBandRed.setValue(col, row, cEntry.c1);
4025  outBandGreen.setValue(col, row, cEntry.c2);
4026  outBandBlue.setValue(col, row, cEntry.c3);
4027  }
4028  else
4029  {
4030  outBandRed.setValue(col, row, rBandNoDataValue);
4031  outBandGreen.setValue(col, row, gBandNoDataValue);
4032  outBandBlue.setValue(col, row, bBandNoDataValue);
4033  }
4034  }
4035  }
4036 
4037  if( enableProgress )
4038  {
4039  progressPtr->pulse();
4040  if( ! progressPtr->isActive() ) return false;
4041  }
4042  }
4043  }
4044 
4045  return true;
4046  }
4047 
4048  } // end namespace rp
4049 } // end namespace te
te::rp::GetMeanValue
bool GetMeanValue(const te::rst::Band &band, const unsigned int maxThreads, double &meanValue)
Get the mean of band pixel values.
Definition: src/terralib/rp/Functions.cpp:2104
te::rp::GetDigitalNumberBandMax
double GetDigitalNumberBandMax(std::string bandName)
Returns the maximum digital number of a given sensor/band.
Definition: src/terralib/rp/Functions.cpp:857
te::rst::Grid::getNumberOfRows
unsigned int getNumberOfRows() const
Returns the grid number of rows.
Definition: raster/Grid.cpp:209
te::rst::Raster::setValue
virtual void setValue(unsigned int c, unsigned int r, const double value, std::size_t b=0)
Sets the attribute value in a band of a cell.
Definition: src/terralib/raster/Raster.cpp:235
te::rst::PaletteIdxCInt
Palette indexes color interpretation.
Definition: src/terralib/raster/Enums.h:58
te::gm::CurvePolygon::push_back
void push_back(Curve *ring)
It adds the curve to the curve polygon.
Definition: CurvePolygon.cpp:137
band
unsigned int band
Definition: TsInterpolator.cpp:50
te::da::DataSourceFactory::make
static std::unique_ptr< DataSource > make(const std::string &driver, const te::core::URI &connInfo)
Definition: src/terralib/dataaccess/datasource/DataSourceFactory.cpp:38
te::rp::RasterHandler::release
void release(std::unique_ptr< te::da::DataSource > &dataSourcePtr, std::unique_ptr< te::da::DataSourceTransactor > &transactorPtr, std::unique_ptr< te::da::DataSet > &dataSetPtr, std::unique_ptr< te::rst::Raster > &rasterPtr)
Relase the internal state and give the ownership to the given pointers.
Definition: RasterHandler.cpp:112
te::rp::RasterHandler::reset
void reset()
Reset the internal state (all internal pointed objects are deleted).
Definition: RasterHandler.cpp:66
te::da::DataSource::getTransactor
virtual std::unique_ptr< DataSourceTransactor > getTransactor()=0
It returns the set of parameters used to set up the access channel to the underlying repository...
te::gm::Coord2D::y
double y
y-coordinate.
Definition: Coord2D.h:114
te::rst::GetDataTypeRanges
TERASTEREXPORT void GetDataTypeRanges(const int &dataType, double &min, double &max)
Return the values range of a given data type.
Definition: src/terralib/raster/Utils.cpp:327
Exception.h
Exception class.
te::rst::RasterIterator::getRow
unsigned int getRow() const
Returns the current row in iterator.
Definition: RasterIterator.h:175
te::rp::Contrast::InputParameters::m_type
ContrastType m_type
The contrast type to be applied.
Definition: Contrast.h:87
te::rst::BandProperty
A raster band description.
Definition: BandProperty.h:61
te::gm::Coord2D::x
double x
x-coordinate.
Definition: Coord2D.h:113
te::rst::RGBPalInt
RGB indexed palette interpretation.
Definition: src/terralib/raster/Enums.h:84
te::rst::Grid::getSRID
int getSRID() const
Returns the grid spatial reference system identifier.
Definition: raster/Grid.cpp:265
te::da::DataSetType
A class that models the description of a dataset.
Definition: DataSetType.h:72
te::rst::Raster::getNumberOfColumns
unsigned int getNumberOfColumns() const
Returns the raster number of columns.
Definition: src/terralib/raster/Raster.cpp:215
te::rst::BandProperty::m_nblocksx
int m_nblocksx
The number of blocks in x.
Definition: BandProperty.h:145
te::rp::B3SplineWAFilter
Spline filter type.
Definition: src/terralib/rp/Functions.h:76
te::rp::SaveSensorParams
void SaveSensorParams(std::map< std::string, SpectralSensorParams > &params)
Saves in SpectralSensor.json file the spectral sensors parameters.
Definition: src/terralib/rp/Functions.cpp:721
te::rst::Interpolator
It interpolates one pixel based on a selected algorithm. Methods currently available are Nearest Neig...
Definition: Interpolator.h:55
te::rst::BandProperty::m_nblocksy
int m_nblocksy
The number of blocks in y.
Definition: BandProperty.h:146
te::rst::RasterIterator
This class implements and iterator to "navigate" over a raster, with a predefined number of bands...
Definition: RasterIterator.h:52
te::rp::RasterSlicing
bool RasterSlicing(const te::rst::Raster &inputRaster, const unsigned int inputRasterBand, const bool createPaletteRaster, const unsigned int slicesNumber, const bool eqHistogram, const std::map< std::string, std::string > &rasterInfo, const std::string &rasterType, const bool enableProgress, std::vector< te::rst::BandProperty::ColorEntry > const *const palettePtr, std::unique_ptr< te::rst::Raster > &outRasterPtr, std::map< double, double > limits)
Generate all wavelet planes from the given input raster.
Definition: src/terralib/rp/Functions.cpp:3787
te::rp::TriangleWAFilter
Triangle filter type.
Definition: src/terralib/rp/Functions.h:77
te::common::Version::asString
static std::string asString()
Definition: Version.cpp:60
te::rp::getSensorParams
std::map< std::string, SpectralSensorParams > getSensorParams()
Returns a map with spectral sensors parameters defined in SpectralSensor.json file.
Definition: src/terralib/rp/Functions.cpp:688
te::gm::Envelope::getUpperRightX
const double & getUpperRightX() const
It returns a constant refernce to the x coordinate of the upper right corner.
Definition: geometry/Envelope.h:414
te::rp::ConvertIHS2RGB
bool ConvertIHS2RGB(const te::rst::Raster &inputIHSRaster, const unsigned int intensityBandIdx, const unsigned int hueBandIdx, const unsigned int saturationBandIdx, const double rgbRangeMin, const double rgbRangeMax, te::rst::Raster &outputRGBRaster)
IHS to RGB conversion.
Definition: src/terralib/rp/Functions.cpp:1282
te::rp::GetDataTypeRange
void TERPEXPORT GetDataTypeRange(const int dataType, double &min, double &max)
Returns the real data type range (all values that can be represented by the given data type)...
Definition: src/terralib/rp/Functions.h:215
te::rp::GetBandNames
std::vector< std::string > GetBandNames()
Returns a vector os with band&#39;s names.
Definition: src/terralib/rp/Functions.cpp:594
te::rp::CreateNewMemRaster
bool CreateNewMemRaster(const te::rst::Grid &rasterGrid, std::vector< te::rst::BandProperty * > bandsProperties, RasterHandler &outRasterHandler)
Create a new raster into a new memory datasource.
Definition: src/terralib/rp/Functions.cpp:260
te::common::TaskProgress
This class can be used to inform the progress of a task.
Definition: TaskProgress.h:53
MIN
#define MIN(a, b)
Macro that returns min between two values.
Definition: GroupingAlgorithms.h:46
te::rst::RedCInt
Red channel color interpretation.
Definition: src/terralib/raster/Enums.h:59
te::gm::Envelope::getLowerLeftY
const double & getLowerLeftY() const
It returns a constant refernce to the y coordinate of the lower left corner.
Definition: geometry/Envelope.h:404
TERPEXPORT
#define TERPEXPORT
You can use this macro in order to export/import classes and functions from this module.
Definition: src/terralib/rp/Config.h:139
te::gm::Coord2D
An utility struct for representing 2D coordinates.
Definition: Coord2D.h:40
M_PI
#define M_PI
Definition: src/terralib/rp/Functions.cpp:78
te::rst::Interpolator::getValue
void getValue(const double &c, const double &r, std::complex< double > &v, const std::size_t &b)
Get the interpolated value at specific band.
Definition: Interpolator.h:93
te::rst::BandProperty::m_type
int m_type
The data type of the elements in the band ( See te::dt namespace basic data types for reference )...
Definition: BandProperty.h:133
te::rp::ComposeBands
bool ComposeBands(te::rp::FeederConstRaster &feeder, const std::vector< unsigned int > &inputRasterBands, const te::rst::Interpolator::Method &interpMethod, const std::map< std::string, std::string > &outputRasterInfo, const std::string &outputDataSourceType, std::unique_ptr< te::rst::Raster > &outputRasterPtr)
Compose a set of bands into one multi-band raster.
Definition: src/terralib/rp/Functions.cpp:3155
te::da::DataSource
An abstract class for data providers like a DBMS, Web Services or a regular file. ...
Definition: src/terralib/dataaccess/datasource/DataSource.h:119
te::rst::BandProperty::m_noDataValue
double m_noDataValue
Value to indicate elements where there is no data, default is std::numeric_limits<double>::max().
Definition: BandProperty.h:136
te::gm::Envelope::getUpperRightY
const double & getUpperRightY() const
It returns a constant refernce to the x coordinate of the upper right corner.
Definition: geometry/Envelope.h:424
te::rp::FeederConstRaster::getCurrentObj
virtual te::rst::Raster const * getCurrentObj() const =0
Return the current sequence object.
te::rp::DirectPrincipalComponents
bool DirectPrincipalComponents(const te::rst::Raster &inputRaster, const std::vector< unsigned int > &inputRasterBands, boost::numeric::ublas::matrix< double > &pcaMatrix, te::rst::Raster &pcaRaster, const std::vector< unsigned int > &pcaRasterBands, const unsigned int maxThreads)
Generate all principal components from the given input raster.
Definition: src/terralib/rp/Functions.cpp:2522
te::rp::DecomposeBands
bool DecomposeBands(const te::rst::Raster &inputRaster, const std::vector< unsigned int > &inputRasterBands, const std::vector< std::map< std::string, std::string > > &outputRastersInfos, const std::string &outputDataSourceType, std::vector< boost::shared_ptr< te::rst::Raster > > &outputRastersPtrs)
Decompose a multi-band raster into a set of one-band rasters.
Definition: src/terralib/rp/Functions.cpp:3089
te::rp::GetDigitalNumberBandMin
double GetDigitalNumberBandMin(std::string bandName)
Returns the minimum digital number of a given sensor/band.
Definition: src/terralib/rp/Functions.cpp:862
te::rst::InterpolationMethod
InterpolationMethod
Allowed interpolation methods.
Definition: src/terralib/raster/Enums.h:92
MAX
#define MAX(a, b)
Macro that returns max between two values.
Definition: GroupingAlgorithms.h:50
te::rst::RasterProperty
Raster property.
Definition: RasterProperty.h:58
te::rp::RemapValues
bool RemapValues(const te::rst::Raster &inputRaster, const std::vector< unsigned int > &inputRasterBands, const boost::numeric::ublas::matrix< double > &remapMatrix, te::rst::Raster &outputRaster, const std::vector< unsigned int > &outputRasterBands, const unsigned int maxThreads)
Remap pixel values using a remap function matrix.
Definition: src/terralib/rp/Functions.cpp:2815
te::rst::Grid::geoToGrid
void geoToGrid(const double &x, const double &y, double &col, double &row) const
Get the grid point associated to a spatial location.
Definition: raster/Grid.cpp:311
te::rst::Grid::getResolutionY
double getResolutionY() const
Returns the grid vertical (y-axis) resolution.
Definition: raster/Grid.cpp:259
te::rp::ConvertRGB2IHS
bool ConvertRGB2IHS(const te::rst::Raster &inputRGBRaster, const unsigned int redBandIdx, const unsigned int greenBandIdx, const unsigned int blueBandIdx, const double rgbRangeMin, const double rgbRangeMax, te::rst::Raster &outputIHSRaster)
RGB to IHS conversion.
Definition: src/terralib/rp/Functions.cpp:1005
TERP_TRUE_OR_RETURN_FALSE
#define TERP_TRUE_OR_RETURN_FALSE(value, message)
Checks if value is true. For false values a warning message will be logged and a return of context wi...
Definition: Macros.h:185
te::rp::ConvertHLS2RGB
bool ConvertHLS2RGB(const te::rst::Raster &inputHLSRaster, const unsigned int hueBandIdx, const unsigned int lightBandIdx, const unsigned int saturationBandIdx, const double, const double rgbRangeMax, te::rst::Raster &outputRGBRaster)
HLS to RGB conversion.
Definition: src/terralib/rp/Functions.cpp:1739
te::rp::RemapValuesThreadParams::m_inputRasterBandsPtr
std::vector< unsigned int > const * m_inputRasterBandsPtr
Definition: src/terralib/rp/Functions.cpp:122
te::rst::Raster::getAccessPolicy
te::common::AccessPolicy getAccessPolicy() const
Returns the raster access policy.
Definition: src/terralib/raster/Raster.cpp:91
Contrast.h
Contrast enhancement.
nCols
unsigned int unsigned int nCols
Definition: TsInterpolator.cpp:37
b
int b
Definition: TsRtree.cpp:32
te::rp::Copy2DiskRaster
bool TERPEXPORT Copy2DiskRaster(const te::rst::Raster &inputRaster, const std::string &fileName)
Create a new raster into a GDAL datasource.
Definition: src/terralib/rp/Functions.cpp:317
te::gm::LinearRing
A LinearRing is a LineString that is both closed and simple.
Definition: LinearRing.h:53
te::common::GetPhysProcNumber
TECOMMONEXPORT unsigned int GetPhysProcNumber()
Returns the number of physical processors.
Definition: PlatformUtils.cpp:249
te::rp::Contrast::InputParameters
Contrast input parameters.
Definition: Contrast.h:65
te::common::Singleton::getInstance
static T & getInstance()
It returns a reference to the singleton instance.
Definition: Singleton.h:126
te::rp::NormalizeRaster
bool NormalizeRaster(te::rst::Raster &inraster, double nmin, double nmax)
Normalizes one raster in a given interval.
Definition: src/terralib/rp/Functions.cpp:867
te::gm::LineString::setPoint
void setPoint(std::size_t i, const double &x, const double &y)
It sets the value of the specified point.
Definition: geometry/LineString.cpp:369
te::rp::FeederConstRaster::reset
virtual void reset()=0
Reset the feeder to the first position (subsequent accesses will start from the first sequence obejct...
te::gm::Envelope
An Envelope defines a 2D rectangular region.
Definition: geometry/Envelope.h:51
te::rst::Raster
An abstract class for raster data strucutures.
Definition: src/terralib/raster/Raster.h:71
te::rst::RasterIterator::begin
static RasterIterator begin(Raster *r, const std::vector< unsigned int > &bands)
Returns an iterator referring to the first value.
Definition: RasterIterator.h:236
te::rst::Raster::getNumberOfRows
unsigned int getNumberOfRows() const
Returns the raster number of rows.
Definition: src/terralib/raster/Raster.cpp:210
te::rst::Raster::getNumberOfBands
virtual std::size_t getNumberOfBands() const =0
Returns the number of bands (dimension of cells attribute values) in the raster.
te::rp::Matrix::getColumnsNumber
unsigned int getColumnsNumber() const
The number of current matrix columns.
Definition: Matrix.h:798
te::rp::InverseWaveletAtrous
bool InverseWaveletAtrous(const te::rst::Raster &waveletRaster, const unsigned int levelsNumber, te::rst::Raster &outputRaster, const std::vector< unsigned int > &outputRasterBands)
Regenerate the original raster from its wavelets planes.
Definition: src/terralib/rp/Functions.cpp:3546
te::rp::WaveletAtrousFilterType
WaveletAtrousFilterType
Wavelet Atrous Filter types.
Definition: src/terralib/rp/Functions.h:73
te::rst::Band::getProperty
BandProperty * getProperty()
Returns the band property.
Definition: src/terralib/raster/Band.cpp:575
te
URI C++ Library.
Definition: Attributes.h:37
te::rst::BandProperty::m_blkw
int m_blkw
Block width (pixels).
Definition: BandProperty.h:143
te::rp::GetDigitalNumberBandInfo
std::pair< double, double > GetDigitalNumberBandInfo(std::string bandName)
Returns the maximun and minimum digital numbers of a given sensor/band.
Definition: src/terralib/rp/Functions.cpp:847
te::rst::Grid::getNumberOfColumns
unsigned int getNumberOfColumns() const
Returns the grid number of columns.
Definition: raster/Grid.cpp:193
te::rp::DirectWaveletAtrous
bool DirectWaveletAtrous(const te::rst::Raster &inputRaster, const std::vector< unsigned int > &inputRasterBands, te::rst::Raster &waveletRaster, const unsigned int levelsNumber, const boost::numeric::ublas::matrix< double > &filter)
Generate all wavelet planes from the given input raster.
Definition: src/terralib/rp/Functions.cpp:3418
te::rst::Grid::getResolutionX
double getResolutionX() const
Returns the grid horizontal (x-axis) resolution.
Definition: raster/Grid.cpp:253
compose.bands
list bands
Definition: compose.py:2
RasterHandler.h
Raster tuple.
te::rst::Band::write
virtual void write(int x, int y, void *buffer)=0
It writes a data block from the specified buffer.
te::rst::RasterSummary
boost::ptr_vector< BandSummary > RasterSummary
RasterSummary is just a typedef of a boost::ptr_vector.
Definition: RasterSummary.h:44
p
te::gm::Polygon * p
Definition: examples/geometry/main.cpp:52
te::rp::Convert2DoublesVector
void Convert2DoublesVector(void *inputVector, const int inputVectorDataType, unsigned int inputVectorSize, double *outputVector)
Convert vector elements.
Definition: src/terralib/rp/Functions.cpp:367
te::rp::RasterHandler
RasterHandler.
Definition: RasterHandler.h:47
te::common::EigenVectors
bool EigenVectors(const boost::numeric::ublas::matrix< T > &inputMatrix, boost::numeric::ublas::matrix< T > &eigenVectorsMatrix, boost::numeric::ublas::matrix< T > &eigenValuesMatrix)
Computes the eigenvectors of a given matrix.
Definition: MatrixUtils.h:250
te::rp::GetMeanValueThread
void GetMeanValueThread(GetMeanValueThreadParams *paramsPtr)
Definition: src/terralib/rp/Functions.cpp:2011
te::rst::Band
A raster band description.
Definition: src/terralib/raster/Band.h:63
te::rst::Raster::getBand
virtual const Band * getBand(std::size_t i) const =0
Returns the raster i-th band.
te::rp::GetDetailedExtent
bool GetDetailedExtent(const te::rst::Grid &grid, te::gm::LinearRing &detailedExtent)
Create a datailed extent from the given grid.
Definition: src/terralib/rp/Functions.cpp:3271
te::rst::Raster::getGrid
Grid * getGrid()
It returns the raster grid.
Definition: src/terralib/raster/Raster.cpp:96
te::common::Exception
This class is designed to declare objects to be thrown as exceptions by TerraLib. ...
Definition: src/terralib/common/Exception.h:58
te::rp::ConvertRGB2HLS
bool ConvertRGB2HLS(const te::rst::Raster &inputRGBRaster, const unsigned int redBandIdx, const unsigned int greenBandIdx, const unsigned int blueBandIdx, const double, const double rgbRangeMax, te::rst::Raster &outputHLSRaster)
RGB to HLS conversion.
Definition: src/terralib/rp/Functions.cpp:1542
te::rp::Matrix::reset
void reset()
Reset (clear) the active instance data.
Definition: Matrix.h:502
te::rst::Raster::getValue
virtual void getValue(unsigned int c, unsigned int r, double &value, std::size_t b=0) const
Returns the attribute value of a band of a cell.
Definition: src/terralib/raster/Raster.cpp:230
te::rp::RemapValuesThread
void RemapValuesThread(RemapValuesThreadParams *paramsPtr)
Definition: src/terralib/rp/Functions.cpp:2637
te::rst::Band::setValue
virtual void setValue(unsigned int c, unsigned int r, const double value)=0
Sets the cell attribute value.
te::srs::Converter
A Converter is responsible for the conversion of coordinates between different Coordinate Systems (CS...
Definition: Converter.h:53
te::rp::CreateNewRaster
bool CreateNewRaster(const te::rst::Grid &rasterGrid, const std::vector< te::rst::BandProperty * > &bandsProperties, const std::string &outDataSetName, const std::string &dataSourceType, RasterHandler &outRasterHandler)
Create a new raster into the givem data source.
Definition: src/terralib/rp/Functions.cpp:131
te::rp::CreateWaveletAtrousFilter
boost::numeric::ublas::matrix< double > CreateWaveletAtrousFilter(const WaveletAtrousFilterType &filterType)
Create a Wavele Atrous Filter.
Definition: src/terralib/rp/Functions.cpp:3378
te::rp::GetCovarianceValue
bool GetCovarianceValue(const te::rst::Band &band1, const te::rst::Band &band2, const unsigned int maxThreads, double const *const mean1ValuePtr, double const *const mean2ValuePtr, double &covarianceValue)
Get the covariance of band pixel values.
Definition: src/terralib/rp/Functions.cpp:2375
te::rst::Band::read
virtual void read(int x, int y, void *buffer) const =0
It reads a data block to the specified buffer.
te::rst::Raster::getSRID
int getSRID() const
Returns the raster spatial reference system identifier.
Definition: src/terralib/raster/Raster.cpp:205
te::rst::BandProperty::ColorEntry::c1
short c1
gray, red, cyan or hue.
Definition: BandProperty.h:72
te::rp::RasterResample
bool RasterResample(const te::rst::Raster &inputRaster, const std::vector< unsigned int > &inputRasterBands, const te::rst::Interpolator::Method interpMethod, const unsigned int firstRow, const unsigned int firstColumn, const unsigned int height, const unsigned int width, const unsigned int newheight, const unsigned int newwidth, const std::map< std::string, std::string > &rinfo, const std::string &dataSourceType, std::unique_ptr< te::rst::Raster > &resampledRasterPtr)
Resample a subset of the raster, given a box.
Definition: src/terralib/rp/Functions.cpp:3637
TERP_LOG_AND_RETURN_FALSE
#define TERP_LOG_AND_RETURN_FALSE(message)
Logs a warning message will and return false.
Definition: Macros.h:269
te::rp::Contrast::initialize
bool initialize(const AlgorithmInputParameters &inputParams)
Initialize the algorithm instance making it ready for execution.
Definition: src/terralib/rp/Contrast.cpp:302
te::rp::GetSpectralBandInfo
SpectralSensorParams GetSpectralBandInfo(std::string bandName)
Returns the maximun and minimum reflectance values of a given sensor/band.
Definition: src/terralib/rp/Functions.cpp:748
te::rp::FeederConstRaster
Feeder from a input rasters.
Definition: FeedersRaster.h:46
te::rst::RasterFactory::make
static Raster * make()
It creates and returns an empty raster with default raster driver.
Definition: src/terralib/raster/RasterFactory.cpp:35
te::gm::Envelope::getLowerLeftX
const double & getLowerLeftX() const
It returns a constant reference to the x coordinate of the lower left corner.
Definition: geometry/Envelope.h:394
te::common::GetInverseMatrix
bool GetInverseMatrix(const boost::numeric::ublas::matrix< T > &inputMatrix, boost::numeric::ublas::matrix< T > &outputMatrix)
Matrix inversion.
Definition: MatrixUtils.h:143
te::rp::Contrast::execute
bool execute(AlgorithmOutputParameters &outputParams)
Executes the algorithm using the supplied parameters.
Definition: src/terralib/rp/Contrast.cpp:179
te::rp::GetIndexedDetailedExtent
bool GetIndexedDetailedExtent(const te::rst::Grid &grid, te::gm::LinearRing &indexedDetailedExtent)
Create a indexed (lines,columns) datailed extent from the given grid.
Definition: src/terralib/rp/Functions.cpp:3326
te::rst::Grid::getExtent
te::gm::Envelope * getExtent()
Returns the geographic extension of the grid.
Definition: raster/Grid.cpp:275
te::rst::BandProperty::ColorEntry::c3
short c3
blue, yellow, or saturation.
Definition: BandProperty.h:74
te::rp::InversePrincipalComponents
bool InversePrincipalComponents(const te::rst::Raster &pcaRaster, const boost::numeric::ublas::matrix< double > &pcaMatrix, te::rst::Raster &outputRaster, const std::vector< unsigned int > &outputRasterBands, const unsigned int maxThreads)
Regenerate the original raster from its principal components.
Definition: src/terralib/rp/Functions.cpp:2612
te::rst::Grid::gridToGeo
void gridToGeo(const double &col, const double &row, double &x, double &y) const
Get the spatial location of a grid point.
Definition: raster/Grid.cpp:305
te::rp::RemapValuesThreadParams::m_outputRasterBandsPtr
std::vector< unsigned int > const * m_outputRasterBandsPtr
Definition: src/terralib/rp/Functions.cpp:124
te::core::GetAppLocalDataLocation
TECOREEXPORT std::string GetAppLocalDataLocation()
It returns the writable folder location to store per user data.
Definition: src/terralib/core/utils/Platform.cpp:151
te::rp::CreateFixedStepPalette
void CreateFixedStepPalette(const unsigned int paletteSize, const bool randomize, std::vector< te::rst::BandProperty::ColorEntry > &palette)
Create a fixed step sequential color palette.
Definition: src/terralib/rp/Functions.cpp:3745
te::rst::RasterIterator::end
static RasterIterator end(Raster *r, const std::vector< unsigned int > &bands)
Returns an iterator referring to after the end of the iterator.
Definition: RasterIterator.h:241
te::rp::RasterHandler::getRasterPtr
te::rst::Raster * getRasterPtr()
Returns a pointer the the handled raster instance or NULL if no instance is handled.
Definition: RasterHandler.cpp:107
te::rst::RasterIterator::getColumn
unsigned int getColumn() const
Returns the current column in iterator.
Definition: RasterIterator.h:180
te::rp::Contrast::OutputParameters::m_createdOutRasterDSType
std::string m_createdOutRasterDSType
Output raster data source type (as described in te::raster::RasterFactory ), leave empty if the resul...
Definition: Contrast.h:198
te::rst::BandProperty::ColorEntry
An structure to represent a color tuple.
Definition: BandProperty.h:70
Functions.h
Raster Processing functions.
te::rst::BandProperty::getType
int getType() const
It returns the data type of the elements in the band.
Definition: BandProperty.h:113
te::rp::FeederConstRaster::moveNext
virtual bool moveNext()=0
Advances to the next sequence obeject.
te::rp::CreateNewGdalRaster
bool CreateNewGdalRaster(const te::rst::Grid &rasterGrid, std::vector< te::rst::BandProperty * > bandsProperties, const std::string &fileName, RasterHandler &outRasterHandler)
Create a new raster into a GDAL datasource.
Definition: src/terralib/rp/Functions.cpp:271
te::rst::BandProperty::m_blkh
int m_blkh
Block height (pixels).
Definition: BandProperty.h:144
te::rst::BlueCInt
Blue channel color interpretation.
Definition: src/terralib/raster/Enums.h:61
te::rst::SUMMARY_MIN
Calculate the min value.
Definition: src/terralib/raster/Enums.h:40
te::rst::Grid
A rectified grid is the spatial support for raster data.
Definition: raster/Grid.h:68
te::rp::RemapValuesThreadParams::m_remapMatrixPtr
boost::numeric::ublas::matrix< double > const * m_remapMatrixPtr
Definition: src/terralib/rp/Functions.cpp:125
te::rst::GreenCInt
Green channel color interpretation.
Definition: src/terralib/raster/Enums.h:60
te::rp::Contrast::OutputParameters
Contrast output parameters.
Definition: Contrast.h:188
te::rst::SUMMARY_MAX
Calculate the max value.
Definition: src/terralib/raster/Enums.h:41
te::rp::GetCovarianceValueThread
void GetCovarianceValueThread(GetCovarianceValueThreadParams *paramsPtr)
Definition: src/terralib/rp/Functions.cpp:2302
te::rst::Band::getBlockSize
virtual int getBlockSize() const
It returns the number of bytes ocuppied by a data block.
Definition: src/terralib/raster/Band.cpp:777
nLines
unsigned int nLines
Definition: TsInterpolator.cpp:37
te::rp::GetSensorFilename
std::string GetSensorFilename()
Returns a json filename with spectral sensors parameters.
Definition: src/terralib/rp/Functions.cpp:582
te::rp::GetStdDevValueThread
void GetStdDevValueThread(GetStdDevValueThreadParams *paramsPtr)
Definition: src/terralib/rp/Functions.cpp:2161
te::rp::Matrix::getLinesNumber
unsigned int getLinesNumber() const
The number of current matrix lines.
Definition: Matrix.h:791
te::rst::BandProperty::ColorEntry::c2
short c2
green, magenta, or lightness.
Definition: BandProperty.h:73
col
unsigned int col
Definition: TsInterpolator.cpp:52
te::rp::GetStdDevValue
bool GetStdDevValue(const te::rst::Band &band, const unsigned int maxThreads, double const *const meanValuePtr, double &stdDevValue)
Get the standard deviation of band pixel values.
Definition: src/terralib/rp/Functions.cpp:2229
te::rp::Contrast
Contrast enhancement.
Definition: Contrast.h:57
te::rst::Band::getValue
virtual void getValue(unsigned int c, unsigned int r, double &value) const =0
Returns the cell attribute value.
te::rp::FeederConstRaster::getCurrentOffset
virtual unsigned int getCurrentOffset() const =0
Return the index of the current object.
te::rp::FeederConstRaster::getObjsCount
virtual unsigned int getObjsCount() const =0
Return the total number of feeder objects.
te::rp::ConvertDoublesVector
void ConvertDoublesVector(double *inputVector, unsigned int inputVectorSize, const int outputVectorDataType, void *outputVector)
Convert a doubles vector.
Definition: src/terralib/rp/Functions.cpp:474