diff --git a/CMakeLists.txt b/CMakeLists.txt index 9e8fda6cccadda3aef21bb48abc0da88eed44b4c..33d751f741a7d0c0f4c836d412dc4522f5cd3354 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -146,7 +146,7 @@ OPTION(OTB_USE_JPEG2000 "Use to support jpeg2000 image file format." ON) #----------------------------------------------------------------------------- # Use the muparser library. -OPTION(OTB_COMPILE_MUPARSER "Use to enable muParser fonctionalities." OFF) +OPTION(OTB_COMPILE_MUPARSER "Use to enable muParser fonctionalities." ON) #----------------------------------------------------------------------------- diff --git a/Code/BasicFilters/CMakeLists.txt b/Code/BasicFilters/CMakeLists.txt index 0b50e4f2bb62387e1f11868e124e1ed294d1579b..730ef790e694761615f524177df3188d6549d355 100644 --- a/Code/BasicFilters/CMakeLists.txt +++ b/Code/BasicFilters/CMakeLists.txt @@ -4,6 +4,11 @@ FILE(GLOB OTBBasicFilters_SRCS "*.cxx" ) ADD_LIBRARY(OTBBasicFilters ${OTBBasicFilters_SRCS}) TARGET_LINK_LIBRARIES (OTBBasicFilters OTBCommon ITKBasicFilters otbedison) + +IF(OTB_COMPILE_MUPARSER) + TARGET_LINK_LIBRARIES(OTBBasicFilters otbmuparser) +ENDIF(OTB_COMPILE_MUPARSER) + IF(OTB_LIBRARY_PROPERTIES) SET_TARGET_PROPERTIES(OTBBasicFilters PROPERTIES ${OTB_LIBRARY_PROPERTIES}) ENDIF(OTB_LIBRARY_PROPERTIES) diff --git a/Code/BasicFilters/otbBinaryParserImageFilter.h b/Code/BasicFilters/otbBinaryParserImageFilter.h new file mode 100644 index 0000000000000000000000000000000000000000..8549634d056ccf27f0e0a6ea66380c75894ad4f3 --- /dev/null +++ b/Code/BasicFilters/otbBinaryParserImageFilter.h @@ -0,0 +1,96 @@ +/*========================================================================= + + Program: ORFEO Toolbox + Language: C++ + Date: $Date$ + Version: $Revision$ + + + Copyright (c) Centre National d'Etudes Spatiales. All rights reserved. + See OTBCopyright.txt for details. + + Some parts of this code are derived from ITK. See ITKCopyright.txt + for details. + + + This software is distributed WITHOUT ANY WARRANTY; without even + the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR + PURPOSE. See the above copyright notices for more information. + +=========================================================================*/ + +#ifndef __otbBinaryParserImageFilter_h +#define __otbBinaryParserImageFilter_h + +#include "itkInPlaceImageFilter.h" +#include "itkImageRegionIteratorWithIndex.h" +#include "itkArray.h" + +#include "otbParser.h" + +namespace otb +{ + +template< class TImage > +class ITK_EXPORT BinaryParserImageFilter + : public itk::InPlaceImageFilter< TImage > +{ +public: + /** Standard class typedefs. */ + typedef BinaryParserImageFilter< TImage > Self; + typedef itk::InPlaceImageFilter< TImage > Superclass; + typedef itk::SmartPointer< Self > Pointer; + typedef itk::SmartPointer< const Self > ConstPointer; + + /** Method for creation through the object factory. */ + itkNewMacro(Self); + + /** Run-time type information (and related methods). */ + itkTypeMacro(BinaryParserImageFilter, InPlaceImageFilter); + + /** Some convenient typedefs. */ + typedef TImage ImageType; + typedef typename ImageType::ConstPointer ImagePointer; + typedef typename ImageType::RegionType ImageRegionType; + typedef typename ImageType::PixelType PixelType; + typedef Parser ParserType; + + void SetInput1( const ImageType * image1); + ImageType * GetInput1(); + void SetInput2( const ImageType * image2); + ImageType * GetInput2(); + void SetExpression(const std::string& expression); + + std::string GetExpression() const; + +protected : + BinaryParserImageFilter(); + virtual ~BinaryParserImageFilter(); + virtual void PrintSelf(std::ostream& os, itk::Indent indent) const; + + void BeforeThreadedGenerateData(); + void ThreadedGenerateData(const ImageRegionType& outputRegionForThread, int threadId ); + void AfterThreadedGenerateData(); + +private : + BinaryParserImageFilter(const Self&); //purposely not implemented + void operator=(const Self&); //purposely not implemented + + std::string m_Expression; + std::vector<ParserType::Pointer> m_VParser; + std::vector<PixelType> m_VImage1; + std::vector<PixelType> m_VImage2; + + long m_UnderflowCount; + long m_OverflowCount; + itk::Array<long> m_ThreadUnderflow; + itk::Array<long> m_ThreadOverflow; +}; + +}//end namespace otb + +#ifndef ITK_MANUAL_INSTANTIATION +#include "otbBinaryParserImageFilter.txx" +#endif + +#endif diff --git a/Code/BasicFilters/otbBinaryParserImageFilter.txx b/Code/BasicFilters/otbBinaryParserImageFilter.txx new file mode 100644 index 0000000000000000000000000000000000000000..c82fca8318e80b09f0b990bd4d29552120837517 --- /dev/null +++ b/Code/BasicFilters/otbBinaryParserImageFilter.txx @@ -0,0 +1,214 @@ +/*========================================================================= + + Program: ORFEO Toolbox + Language: C++ + Date: $Date$ + Version: $Revision$ + + + Copyright (c) Centre National d'Etudes Spatiales. All rights reserved. + See OTBCopyright.txt for details. + + Some parts of this code are derived from ITK. See ITKCopyright.txt + for details. + + + This software is distributed WITHOUT ANY WARRANTY; without even + the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR + PURPOSE. See the above copyright notices for more information. + +=========================================================================*/ +#ifndef __otbBinaryParserImageFilter_txx +#define __otbBinaryParserImageFilter_txx +#include "otbBinaryParserImageFilter.h" + +#include "itkImageRegionIterator.h" +#include "itkImageRegionConstIterator.h" +#include "itkNumericTraits.h" +#include "itkProgressReporter.h" +#include "otbMacro.h" + +namespace otb +{ + +/** Constructor */ +template <class TImage> +BinaryParserImageFilter<TImage> +::BinaryParserImageFilter() +{ + this->SetNumberOfRequiredInputs( 2 ); + this->InPlaceOff(); + + m_UnderflowCount = 0; + m_OverflowCount = 0; + m_ThreadUnderflow.SetSize(1); + m_ThreadOverflow.SetSize(1); +} + +/** Destructor */ +template <class TImage> +BinaryParserImageFilter<TImage> +::~BinaryParserImageFilter() +{ +} + +template <class TImage> +void BinaryParserImageFilter<TImage> +::PrintSelf(std::ostream& os, itk::Indent indent) const +{ + Superclass::PrintSelf(os, indent); + + os << indent << "Expression: " << m_Expression << std::endl; + os << indent << "Computed values follow:" << std::endl; + os << indent << "UnderflowCount: " << m_UnderflowCount << std::endl; + os << indent << "OverflowCount: " << m_OverflowCount << std::endl; + os << indent << "itk::NumericTraits<PixelType>::NonpositiveMin() : " + << itk::NumericTraits<PixelType>::NonpositiveMin() << std::endl; + os << indent << "itk::NumericTraits<PixelType>::max() : " + << itk::NumericTraits<PixelType>::max() << std::endl; +} + +template <class TImage> +void BinaryParserImageFilter<TImage> +::SetInput1( const ImageType * image1) +{ + this->SetNthInput(0, const_cast<TImage *>( image1 )); +} + +template <typename TImage> +TImage * BinaryParserImageFilter<TImage> +::GetInput1() +{ + return const_cast<TImage *>(this->GetInput(0)); +} + +template <class TImage> +void BinaryParserImageFilter<TImage> +::SetInput2( const ImageType * image2) +{ + this->SetNthInput(1, const_cast<TImage *>( image2 )); +} + +template <typename TImage> +TImage * BinaryParserImageFilter<TImage> +::GetInput2() +{ + return const_cast<TImage *>(this->GetInput(1)); +} + +template< typename TImage > +void BinaryParserImageFilter<TImage> +::SetExpression(const std::string& expression) +{ + if (m_Expression != expression) + m_Expression = expression; + this->Modified(); +} + +template< typename TImage > +std::string BinaryParserImageFilter<TImage> +::GetExpression() const +{ + return m_Expression; +} + +template< typename TImage > +void BinaryParserImageFilter<TImage> +::BeforeThreadedGenerateData() +{ + typename std::vector<ParserType::Pointer>::iterator itParser; + int nbThreads = this->GetNumberOfThreads(); + unsigned int i; + + // Allocate and initialize the thread temporaries + m_ThreadUnderflow.SetSize(nbThreads); + m_ThreadUnderflow.Fill(0); + m_ThreadOverflow.SetSize(nbThreads); + m_ThreadOverflow.Fill(0); + m_VParser.resize(nbThreads); + m_VImage1.resize(nbThreads); + m_VImage2.resize(nbThreads); + for(itParser = m_VParser.begin(); itParser < m_VParser.end(); itParser++) + *itParser = ParserType::New(); + + for(i = 0; i < nbThreads; i++) + { + m_VParser.at(i)->SetExpr(m_Expression); + m_VParser.at(i)->DefineVar("Image1", &m_VImage1.at(i)); + m_VParser.at(i)->DefineVar("Image2", &m_VImage2.at(i)); + } +} + +template< typename TImage > +void BinaryParserImageFilter<TImage> +::AfterThreadedGenerateData() +{ + int nbThreads = this->GetNumberOfThreads(); + + m_UnderflowCount = 0; + m_OverflowCount = 0; + + // Accumulate counts for each thread + for( int i = 0; i < nbThreads; i++) + { + m_UnderflowCount += m_ThreadUnderflow[i]; + m_OverflowCount += m_ThreadOverflow[i]; + } + + if((m_UnderflowCount != 0) || (m_OverflowCount!=0)) + otbWarningMacro(<< std::endl + << "The Following Parsed Expression : " + << this->GetExpression() << std::endl + << "Generated " << m_UnderflowCount << " Underflow(s) " + << "And " << m_OverflowCount << " Overflow(s) " << std::endl + << "The Parsed Expression, The Inputs And The Output " + << "Type May Be Incompatible !"); + +} + +template< typename TImage > +void BinaryParserImageFilter<TImage> +::ThreadedGenerateData(const ImageRegionType& outputRegionForThread, + int threadId) +{ + PixelType value; + itk::ImageRegionConstIterator<TImage> it1 (this->GetInput1(), outputRegionForThread); + itk::ImageRegionConstIterator<TImage> it2 (this->GetInput2(), outputRegionForThread); + itk::ImageRegionIterator<TImage> ot (this->GetOutput(), outputRegionForThread); + + // support progress methods/callbacks + itk::ProgressReporter progress(this, threadId, outputRegionForThread.GetNumberOfPixels()); + + while (!it1.IsAtEnd()) + { + m_VImage1.at(threadId) = it1.Get(); + m_VImage2.at(threadId) = it2.Get(); + + value = static_cast<double>(m_VParser.at(threadId)->Eval()); + + if (value < itk::NumericTraits<PixelType>::NonpositiveMin()) + { + ot.Set(itk::NumericTraits<PixelType>::NonpositiveMin()); + m_ThreadUnderflow[threadId]++; + } + else if (value > itk::NumericTraits<PixelType>::max()) + { + ot.Set(itk::NumericTraits<PixelType>::max()); + m_ThreadOverflow[threadId]++; + } + else + { + ot.Set(static_cast<PixelType>(value)); + } + + ++it1; + ++it2; + ++ot; + + progress.CompletedPixel(); + } +} + +}// end namespace otb + +#endif diff --git a/Code/Common/otbParser.cxx b/Code/Common/otbParser.cxx index ec169247f797eaf5abccff89c5ed0bad0695e38f..29569543e7cc319407727b4c4576ae08e602154c 100644 --- a/Code/Common/otbParser.cxx +++ b/Code/Common/otbParser.cxx @@ -22,9 +22,7 @@ namespace otb { -void -Parser -::PrintSelf(std::ostream& os, itk::Indent indent) const +void Parser::PrintSelf(std::ostream& os, itk::Indent indent) const { Superclass::PrintSelf(os, indent); } @@ -50,25 +48,30 @@ Parser::Parser() InitConst(); } +Parser::~Parser() +{ + +} + void Parser::InitConst() { - m_Parser.DefineConst( "e", CONST_E ); - m_Parser.DefineConst( "log2e", CONST_LOG2E ); - m_Parser.DefineConst( "log10e", CONST_LOG10E ); - m_Parser.DefineConst( "ln2", CONST_LN2 ); - m_Parser.DefineConst( "ln10", CONST_LN10 ); - m_Parser.DefineConst( "pi", CONST_PI ); - m_Parser.DefineConst( "euler", CONST_EULER ); + m_InternalParser.DefineConst( "e", CONST_E ); + m_InternalParser.DefineConst( "log2e", CONST_LOG2E ); + m_InternalParser.DefineConst( "log10e", CONST_LOG10E ); + m_InternalParser.DefineConst( "ln2", CONST_LN2 ); + m_InternalParser.DefineConst( "ln10", CONST_LN10 ); + m_InternalParser.DefineConst( "pi", CONST_PI ); + m_InternalParser.DefineConst( "euler", CONST_EULER ); } void Parser::InitFun() { - m_Parser.DefineFun("ndvi", NDVI); + m_InternalParser.DefineFun("ndvi", NDVI); } void Parser::SetExpr(const std::string & Expression) { - m_Parser.SetExpr(Expression); + m_InternalParser.SetExpr(Expression); } Parser::ValueType Parser::Eval() @@ -76,7 +79,7 @@ Parser::ValueType Parser::Eval() Parser::ValueType result; try { - result = m_Parser.Eval(); + result = m_InternalParser.Eval(); } catch(ParserType::ExceptionType &e) { @@ -89,7 +92,7 @@ void Parser::DefineVar(const std::string &sName, Parser::ValueType *fVar) { try { - m_Parser.DefineVar(sName, fVar); + m_InternalParser.DefineVar(sName, fVar); } catch(ParserType::ExceptionType &e) { @@ -99,7 +102,12 @@ void Parser::DefineVar(const std::string &sName, Parser::ValueType *fVar) void Parser::ClearVar() { - m_Parser.ClearVar(); + m_InternalParser.ClearVar(); +} + +const std::string& Parser::GetExpr() const +{ + return m_InternalParser.GetExpr(); } void Parser::ExceptionHandler(ParserType::ExceptionType &e) diff --git a/Code/Common/otbParser.h b/Code/Common/otbParser.h index c845117fd8002fa4160d556c9e0d423497d58793..26a77a1b260dc5cc76c5191751565f3c381bac1f 100644 --- a/Code/Common/otbParser.h +++ b/Code/Common/otbParser.h @@ -44,12 +44,12 @@ public: virtual void InitConst(); virtual void InitFun(); + virtual void SetExpr(const std::string & Expression); ValueType Eval(); void DefineVar(const std::string &sName, ValueType *fVar); void ClearVar(); - - + const std::string& GetExpr() const; typedef mu::Parser::exception_type ExceptionType; virtual void ExceptionHandler(ExceptionType &e); @@ -57,7 +57,7 @@ public: protected: Parser(); - virtual ~Parser() {} + virtual ~Parser(); virtual void PrintSelf(std::ostream& os, itk::Indent indent) const; @@ -65,7 +65,7 @@ private: Parser(const Self &); //purposely not implemented void operator =(const Self &); //purposely not implemented - mu::Parser m_Parser; + mu::Parser m_InternalParser; //---------- User Defined Functions ----------//BEGIN diff --git a/Testing/Code/BasicFilters/CMakeLists.txt b/Testing/Code/BasicFilters/CMakeLists.txt index 0504f2b4cd796cf4f5c5d074f7cd38e7bb532e35..780154ac2b205949939e00bbed67f2108bd25e7f 100644 --- a/Testing/Code/BasicFilters/CMakeLists.txt +++ b/Testing/Code/BasicFilters/CMakeLists.txt @@ -28,6 +28,9 @@ SET(BASICFILTERS_TESTS9 ${CXX_TEST_PATH}/otbBasicFiltersTests9) SET(BASICFILTERS_TESTS10 ${CXX_TEST_PATH}/otbBasicFiltersTests10) SET(BASICFILTERS_TESTS11 ${CXX_TEST_PATH}/otbBasicFiltersTests11) SET(BASICFILTERS_TESTS12 ${CXX_TEST_PATH}/otbBasicFiltersTests12) +IF(OTB_COMPILE_MUPARSER) + SET(BASICFILTERS_TESTS13 ${CXX_TEST_PATH}/otbBasicFiltersTests13) +ENDIF(OTB_COMPILE_MUPARSER) # ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ # ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ otbBasicFiltersTests1 ~~~~~~~~~~~~~~~~~~~~~~~~~~ @@ -1621,6 +1624,22 @@ ADD_TEST(bfTvImageAndVectorImageOperationFilterTest ${BASICFILTERS_TESTS12} ) +# ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ +# ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ otbBasicFiltersTests13 ~~~~~~~~~~~~~~~~~~~~~~~~~ +# ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + +# ---------------------------- otbNaryParserFunctorImageFilter ---------------------------- +IF(OTB_COMPILE_MUPARSER) + +ADD_TEST(bfTuBinaryParserImageFilterNew ${BASICFILTERS_TESTS13} + otbBinaryParserImageFilterNew) + +ADD_TEST(bfTvBinaryParserImageFilter ${BASICFILTERS_TESTS13} + otbBinaryParserImageFilter) + +ENDIF(OTB_COMPILE_MUPARSER) + + # A enrichir SET(BasicFilters_SRCS1 otbBasicFiltersTests1.cxx @@ -1849,6 +1868,15 @@ otbImageAndVectorImageOperationFilterNew.cxx otbImageAndVectorImageOperationFilterTest.cxx ) +IF(OTB_COMPILE_MUPARSER) +# A enrichir +SET(BasicFilters_SRCS13 +otbBasicFiltersTests13.cxx +otbBinaryParserImageFilterNew.cxx +otbBinaryParserImageFilter.cxx +) +ENDIF(OTB_COMPILE_MUPARSER) + OTB_ADD_EXECUTABLE(otbBasicFiltersTests1 "${BasicFilters_SRCS1}" "OTBBasicFilters;OTBIO;OTBTesting") OTB_ADD_EXECUTABLE(otbBasicFiltersTests2 "${BasicFilters_SRCS2}" "OTBBasicFilters;OTBIO;OTBTesting") OTB_ADD_EXECUTABLE(otbBasicFiltersTests3 "${BasicFilters_SRCS3}" "OTBBasicFilters;OTBIO;OTBTesting") @@ -1861,6 +1889,8 @@ OTB_ADD_EXECUTABLE(otbBasicFiltersTests9 "${BasicFilters_SRCS9}" "OTBBasicFilter OTB_ADD_EXECUTABLE(otbBasicFiltersTests10 "${BasicFilters_SRCS10}" "OTBBasicFilters;OTBIO;OTBTesting") OTB_ADD_EXECUTABLE(otbBasicFiltersTests11 "${BasicFilters_SRCS11}" "OTBBasicFilters;ITKAlgorithms;OTBIO;OTBTesting") OTB_ADD_EXECUTABLE(otbBasicFiltersTests12 "${BasicFilters_SRCS12}" "OTBBasicFilters;ITKAlgorithms;OTBIO;OTBTesting") - +IF(OTB_COMPILE_MUPARSER) + OTB_ADD_EXECUTABLE(otbBasicFiltersTests13 "${BasicFilters_SRCS13}" "OTBBasicFilters;OTBIO;OTBTesting") +ENDIF(OTB_COMPILE_MUPARSER) ENDIF( NOT OTB_DISABLE_CXX_TESTING AND BUILD_TESTING ) diff --git a/Testing/Code/BasicFilters/otbBasicFiltersTests13.cxx b/Testing/Code/BasicFilters/otbBasicFiltersTests13.cxx new file mode 100644 index 0000000000000000000000000000000000000000..9de9ecf8bd5e10cbac9a0a489a71326122500e7a --- /dev/null +++ b/Testing/Code/BasicFilters/otbBasicFiltersTests13.cxx @@ -0,0 +1,31 @@ +/*========================================================================= + + Program: ORFEO Toolbox + Language: C++ + Date: $Date$ + Version: $Revision$ + + + Copyright (c) Centre National d'Etudes Spatiales. All rights reserved. + See OTBCopyright.txt for details. + + + This software is distributed WITHOUT ANY WARRANTY; without even + the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR + PURPOSE. See the above copyright notices for more information. + +=========================================================================*/ + +// this file defines the otbCommonTest for the test driver +// and all it expects is that you have a function called RegisterTests +#if defined(_MSC_VER) +#pragma warning ( disable : 4786 ) +#endif + +#include "otbTestMain.h" + +void RegisterTests() +{ + REGISTER_TEST(otbBinaryParserImageFilterNew); + REGISTER_TEST(otbBinaryParserImageFilter); +} diff --git a/Testing/Code/BasicFilters/otbBinaryParserImageFilter.cxx b/Testing/Code/BasicFilters/otbBinaryParserImageFilter.cxx new file mode 100644 index 0000000000000000000000000000000000000000..39a9308c65ee5325e09edefdad421a0748653dc3 --- /dev/null +++ b/Testing/Code/BasicFilters/otbBinaryParserImageFilter.cxx @@ -0,0 +1,139 @@ +/*========================================================================= + + Program: ORFEO Toolbox + Language: C++ + Date: $Date$ + Version: $Revision$ + + + Copyright (c) Centre National d'Etudes Spatiales. All rights reserved. + See OTBCopyright.txt for details. + + + This software is distributed WITHOUT ANY WARRANTY; without even + the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR + PURPOSE. See the above copyright notices for more information. + +=========================================================================*/ +#if defined(_MSC_VER) +#pragma warning ( disable : 4786 ) +#endif + +#include "itkExceptionObject.h" +#include <iostream> + +#include "otbMath.h" +#include "otbImage.h" +#include "otbBinaryParserImageFilter.h" + +int otbBinaryParserImageFilter( int argc, char* argv[]) +{ + typedef double PixelType; + typedef otb::Image<PixelType, 2> ImageType; + typedef otb::BinaryParserImageFilter<ImageType> FilterType; + + unsigned int i; + const unsigned int N = 100; + + ImageType::SizeType size; + size.Fill(N); + ImageType::IndexType index; + index.Fill(0); + ImageType::RegionType region; + region.SetSize(size); + region.SetIndex(index); + + ImageType::Pointer image1 = ImageType::New(); + ImageType::Pointer image2 = ImageType::New(); + + image1->SetLargestPossibleRegion( region ); + image1->SetBufferedRegion( region ); + image1->SetRequestedRegion( region ); + image1->Allocate(); + + image2->SetLargestPossibleRegion( region ); + image2->SetBufferedRegion( region ); + image2->SetRequestedRegion( region ); + image2->Allocate(); + + typedef itk::ImageRegionIteratorWithIndex<ImageType> IteratorType; + IteratorType it1(image1, region); + IteratorType it2(image2, region); + + for (it1.GoToBegin(), it2.GoToBegin(); !it1.IsAtEnd(); ++it1, ++it2) + { + ImageType::IndexType i1 = it1.GetIndex(); + ImageType::IndexType i2 = it2.GetIndex(); + + it1.Set( i1[0] + i1[1] -50 ); + it2.Set( i2[0] * i2[1] ); + } + + + FilterType::Pointer filter = FilterType::New(); + std::cout << "Number Of Threads : " << filter->GetNumberOfThreads() << std::endl; + + filter->SetInput1(image1); + filter->SetInput2(image2); + + filter->SetExpression("cos(2 * pi * Image1)/(2 * pi * Image2 + 1E-3) + ndvi(Image1,Image2)*sqrt(2)"); + filter->Update(); + + std::cout << "--- Standard Use" << std::endl; + std::cout << "Parsed Expression : " << filter->GetExpression() << std::endl; + std::cout << "About the Filter : " << filter << std::endl; + + ImageType::Pointer output = filter->GetOutput(); + IteratorType it(output, region); + + for (it1.GoToBegin(), it2.GoToBegin(), it.GoToBegin(); !it1.IsAtEnd(); ++it1, ++it2, ++it) + { + ImageType::IndexType i1 = it1.GetIndex(); + ImageType::IndexType i2 = it2.GetIndex(); + double px1 = ( i1[0] + i1[1] -50 ); + double px2 = ( i2[0] * i2[1] ); + double result = it.Get(); + double ndvi_expected; + + if ( vcl_abs( px1 + px2) < 1E-6 ) + ndvi_expected = 0.0; + else + ndvi_expected = (px2-px1)/(px2+px1); + + double expected = vcl_cos( 2 * otb::CONST_PI * px1 ) / ( 2 * otb::CONST_PI * px2 + 1E-3 ) + + ndvi_expected * vcl_sqrt(2); + + /* + std::cout << "Pixel_1 = " << it1.Get() << " Pixel_2 = " << it2.Get() + << " Result = " << it.Get() << " Expected = " << expected << std::endl; + */ + + if ( vcl_abs( result - expected ) > 1E-10 ) + { + itkGenericExceptionMacro( << "Error > 1E-10 -> Test Failled" << std::endl + << "Pixel_1 = " << it1.Get() << " Pixel_2 = " << it2.Get() + << " Result = " << it.Get() << " Expected = " << expected + << std::endl ); + } + } + std::cout << "Process Completed" << std::endl; + + /** Edge Effect Handling */ + + std::cout << "--- Edge Effect Handling" << std::endl; + std::cout << "- +/-inf section" << std::endl; + filter->SetExpression("Image1 / Image2"); + filter->Update(); + + std::cout << "- nan section" << std::endl; + it1.GoToBegin(); it2.GoToBegin(); it.GoToBegin(); + for(i=1; i<=50; i++ , ++it1, ++it2, ++it){} + if(isnan(it.Get())) + std::cout << "Pixel_1 = " << it1.Get() << " Pixel_2 = " << it2.Get() << " Result = " << it.Get() << " Expected = nan" << std::endl; + else + itkGenericExceptionMacro( << "Error > Bad Edge Effect Handling -> Test Failled" << std::endl + << "Pixel_1 = " << it1.Get() << " Pixel_2 = " << it2.Get() + << " Result = " << it.Get() << " Expected = nan" << std::endl ); + + return EXIT_SUCCESS; +} diff --git a/Testing/Code/BasicFilters/otbBinaryParserImageFilterNew.cxx b/Testing/Code/BasicFilters/otbBinaryParserImageFilterNew.cxx new file mode 100644 index 0000000000000000000000000000000000000000..85a2f1f18b9085773a8c62be49e718ea14c94141 --- /dev/null +++ b/Testing/Code/BasicFilters/otbBinaryParserImageFilterNew.cxx @@ -0,0 +1,39 @@ +/*========================================================================= + + Program: ORFEO Toolbox + Language: C++ + Date: $Date$ + Version: $Revision$ + + + Copyright (c) Centre National d'Etudes Spatiales. All rights reserved. + See OTBCopyright.txt for details. + + + This software is distributed WITHOUT ANY WARRANTY; without even + the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR + PURPOSE. See the above copyright notices for more information. + +=========================================================================*/ +#if defined(_MSC_VER) +#pragma warning ( disable : 4786 ) +#endif + +#include "itkExceptionObject.h" +#include <iostream> + +#include "otbImage.h" +#include "otbBinaryParserImageFilter.h" + +int otbBinaryParserImageFilterNew( int argc, char* argv[]) +{ + typedef double PixelType; + typedef otb::Image<PixelType, 2> ImageType; + typedef otb::BinaryParserImageFilter<ImageType> FilterType; + + FilterType::Pointer filter = FilterType::New(); + + std::cout << "Number Of Threads : " << filter->GetNumberOfThreads() << std::endl; + + return EXIT_SUCCESS; +}