diff --git a/Modules/Applications/AppHyperspectral/app/otbHyperspectralUnmixing.cxx b/Modules/Applications/AppHyperspectral/app/otbHyperspectralUnmixing.cxx index f4a0bd1e5f04dcac2d33eb6229c1c2805f4b9564..371b301c40e4d51548c26dead99926726b7abccf 100644 --- a/Modules/Applications/AppHyperspectral/app/otbHyperspectralUnmixing.cxx +++ b/Modules/Applications/AppHyperspectral/app/otbHyperspectralUnmixing.cxx @@ -23,8 +23,6 @@ #include "otbUnConstrainedLeastSquareImageFilter.h" #include "otbISRAUnmixingImageFilter.h" -#include "otbNCLSUnmixingImageFilter.h" -//#include "otbFCLSUnmixingImageFilter.h" #include "otbMDMDNMFImageFilter.h" @@ -34,8 +32,6 @@ namespace Wrapper { typedef otb::UnConstrainedLeastSquareImageFilter<DoubleVectorImageType, DoubleVectorImageType, double> UCLSUnmixingFilterType; typedef otb::ISRAUnmixingImageFilter<DoubleVectorImageType, DoubleVectorImageType, double> ISRAUnmixingFilterType; -typedef otb::NCLSUnmixingImageFilter<DoubleVectorImageType, DoubleVectorImageType, double> NCLSUnmixingFilterType; -//typedef otb::FCLSUnmixingImageFilter<DoubleVectorImageType, DoubleVectorImageType, double> FCLSUnmixingFilterType; typedef otb::MDMDNMFImageFilter<DoubleVectorImageType, DoubleVectorImageType> MDMDNMFUnmixingFilterType; typedef otb::VectorImageToMatrixImageFilter<DoubleVectorImageType> VectorImageToMatrixImageFilterType; @@ -65,12 +61,12 @@ enum UnMixingMethod { UnMixingMethod_UCLS, //UnMixingMethod_FCLS, - UnMixingMethod_NCLS, + // UnMixingMethod_NCLS, UnMixingMethod_ISRA, UnMixingMethod_MDMDNMF, }; -const char* UnMixingMethodNames [] = { "UCLS", "FCLS", "NCLS", "ISRA", "MDMDNMF", }; +const char* UnMixingMethodNames [] = { "UCLS", "ISRA", "MDMDNMF", }; class HyperspectralUnmixing : public Application @@ -104,9 +100,7 @@ private: "be estimated using the VertexComponentAnalysis application.\n\n" "The application allows estimating the abundance maps with several algorithms:\n\n" "* Unconstrained Least Square (ucls)\n" -// "* Fully Constrained Least Square (fcls)\n" "* Image Space Reconstruction Algorithm (isra)\n" - "* Non-negative constrained\n" "* Least Square (ncls)\n" "* Minimum Dispersion Constrained Non Negative Matrix Factorization (MDMDNMF)." ); @@ -135,12 +129,6 @@ private: AddChoice("ua.ucls", "UCLS"); SetParameterDescription("ua.ucls", "Unconstrained Least Square"); -// AddChoice("ua.fcls", "FCLS"); -// SetParameterDescription("ua.fcls", "Fully constrained Least Square"); - - AddChoice("ua.ncls", "NCLS"); - SetParameterDescription("ua.ncls", "Non-negative constrained Least Square"); - AddChoice("ua.isra", "ISRA"); SetParameterDescription("ua.isra", "Image Space Reconstruction Algorithm"); @@ -217,36 +205,6 @@ private: } break; - case UnMixingMethod_NCLS: - { - otbAppLogINFO("NCLS Unmixing"); - - NCLSUnmixingFilterType::Pointer unmixer = - NCLSUnmixingFilterType::New(); - - unmixer->SetInput(inputImage); - unmixer->GetModifiableFunctor().SetEndmembersMatrix(endMembersMatrix); - abundanceMap = unmixer->GetOutput(); - m_ProcessObjects.push_back(unmixer.GetPointer()); - - } - break; - /* - case UnMixingMethod_FCLS: - { - otbAppLogINFO("FCLS Unmixing"); - - FCLSUnmixingFilterType::Pointer unmixer = - FCLSUnmixingFilterType::New(); - - unmixer->SetInput(inputImage); - unmixer->SetEndmembersMatrix(endMembersMatrix); - abundanceMap = unmixer->GetOutput(); - m_ProcessObjects.push_back(unmixer.GetPointer()); - - } - break; - */ case UnMixingMethod_MDMDNMF: { otbAppLogINFO("MDMD-NMF Unmixing"); diff --git a/Modules/Applications/AppHyperspectral/test/CMakeLists.txt b/Modules/Applications/AppHyperspectral/test/CMakeLists.txt index a430980733388e9732cf69f2270e0396902f904c..9d7a04bf7e2451f4075256ec446f6c3ccedfe36f 100644 --- a/Modules/Applications/AppHyperspectral/test/CMakeLists.txt +++ b/Modules/Applications/AppHyperspectral/test/CMakeLists.txt @@ -30,17 +30,6 @@ otb_test_application(NAME apTvHyHyperspectralUnmixing_UCLS ${BASELINE}/apTvHyHyperspectralUnmixing_UCLS.tif ${TEMP}/apTvHyHyperspectralUnmixing_UCLS.tif) -otb_test_application(NAME apTvHyHyperspectralUnmixing_NCLS - APP HyperspectralUnmixing - OPTIONS -in ${INPUTDATA}/Hyperspectral/synthetic/hsi_cube.tif - -ie ${INPUTDATA}/Hyperspectral/synthetic/endmembers.tif - -out ${TEMP}/apTvHyHyperspectralUnmixing_NCLS.tif double - -ua ncls - VALID --compare-image ${EPSILON_9} - ${BASELINE}/apTvHyHyperspectralUnmixing_NCLS.tif - ${TEMP}/apTvHyHyperspectralUnmixing_NCLS.tif) - - #----------- VertexComponentAnalysis TESTS ---------------- otb_test_application(NAME apTvHyVertexComponentAnalysis APP VertexComponentAnalysis diff --git a/Modules/Core/Transform/include/otbStreamingWarpImageFilter.hxx b/Modules/Core/Transform/include/otbStreamingWarpImageFilter.hxx index dcd43c075fe385163904dc59014297722b57abf6..741f7b617dc6ed9f96f6724b7a14592f23fd7ab4 100644 --- a/Modules/Core/Transform/include/otbStreamingWarpImageFilter.hxx +++ b/Modules/Core/Transform/include/otbStreamingWarpImageFilter.hxx @@ -230,24 +230,9 @@ StreamingWarpImageFilter<TInputImage, TOutputImage, TDisplacementField> { inputPtr->SetRequestedRegion(inputRequestedRegion); } - else - { - - inputFinalSize.Fill(0); - inputRequestedRegion.SetSize(inputFinalSize); - inputFinalIndex.Fill(0); - inputRequestedRegion.SetIndex(inputFinalIndex); - - // store what we tried to request (prior to trying to crop) - inputPtr->SetRequestedRegion(inputRequestedRegion); - -// // build an exception -// itk::InvalidRequestedRegionError e(__FILE__, __LINE__); -// e.SetLocation(ITK_LOCATION); -// e.SetDescription("Requested region is (at least partially) outside the largest possible region."); -// e.SetDataObject(inputPtr); -// throw e; - } + // Here we do not throw an exception, we just do nothing since + // resampling filter can legitimately ask for a region outside of + // input's largest possible region. } diff --git a/Modules/Hyperspectral/Unmixing/include/otbNCLSUnmixingImageFilter.h b/Modules/Hyperspectral/Unmixing/include/otbNCLSUnmixingImageFilter.h deleted file mode 100644 index aa3ef08ce40742422f2a3db3e12875be8ad46fe1..0000000000000000000000000000000000000000 --- a/Modules/Hyperspectral/Unmixing/include/otbNCLSUnmixingImageFilter.h +++ /dev/null @@ -1,131 +0,0 @@ -/* - * Copyright (C) 2005-2019 Centre National d'Etudes Spatiales (CNES) - * - * This file is part of Orfeo Toolbox - * - * https://www.orfeo-toolbox.org/ - * - * Licensed under the Apache License, Version 2.0 (the "License"); - * you may not use this file except in compliance with the License. - * You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ - -#ifndef otbNCLSUnmixingImageFilter_h -#define otbNCLSUnmixingImageFilter_h - -#include "itkMacro.h" -#include "itkNumericTraits.h" -#include "otbFunctorImageFilter.h" -#include "vnl/algo/vnl_svd.h" -#include <boost/shared_ptr.hpp> - -namespace otb -{ - -namespace Functor { - -/** \class NCLSUnmixingFunctor - * - * \brief Performs fully constrained least squares on a pixel - * - * \sa NCLSUnmixingImageFilter - * - * \ingroup OTBUnmixing - */ -template<class TInput, class TOutput, class TPrecision> -class NCLSUnmixingFunctor -{ -public: - typedef NCLSUnmixingFunctor<TInput, TOutput, TPrecision> Self; - - typedef TInput InputType; - typedef TOutput OutputType; - typedef TPrecision PrecisionType; - - typedef vnl_vector<PrecisionType> VectorType; - typedef vnl_matrix<PrecisionType> MatrixType; - - NCLSUnmixingFunctor(); - virtual ~NCLSUnmixingFunctor() = default; - - size_t OutputSize(const std::array<size_t, 1>& nbBands) const; - - void SetEndmembersMatrix(const MatrixType& U); - const MatrixType& GetEndmembersMatrix(void) const; - - void SetMaxIteration(unsigned int val) - { - m_MaxIteration = val; - } - - unsigned int GetMaxIteration() const - { - return m_MaxIteration; - } - - OutputType operator ()(const InputType& in) const; - -private: - - static bool IsNonNegative(PrecisionType val) - { - return val >= 0; - } - - typedef vnl_svd<PrecisionType> SVDType; - typedef boost::shared_ptr<SVDType> SVDPointerType; - - MatrixType m_U; - MatrixType m_Ut; - MatrixType m_UtUinv; - SVDPointerType m_Svd; // SVD of U - unsigned int m_OutputSize; - unsigned int m_MaxIteration; -}; -} - -/** \class NCLSUnmixingImageFilter - * - * \brief Performs fully constrained least squares on each pixel of a VectorImage - * - * This filter takes as input a multiband image and a matrix. - * If the matrix is called \f$A\f$, it solves, for each pixel \f$p\f$, the system - * \f$A \cdot x = p\f$ in the least square sense, with additional constraints on the solution - * \f$\hat{x}\f$ ensuring positivity (each component is positive) and additivity (the sum of - * all components is 1). - * - * The main use of this filter is to unmix an hyperspectral dataset, - * where \f$A\f$ is the mixing matrix, in which each row corresponds to an endmember signature. - * - * The number of rows in \f$A\f$ must match the input image number of bands. - * The number of bands in the output image will be the number of columns of $A$ - * - * References - * "Fully Constrained Least-Squares Based Linear Unmixing." Daniel Heinz, - * Chein-I Chang, and Mark L.G. Althouse. IEEE. 1999. - * - * \ingroup Hyperspectral - * \ingroup Streamed - * \ingroup Threaded - * - * \ingroup OTBUnmixing - */ -template <typename TInputImage, typename TOutputImage, typename TPrecision> -using NCLSUnmixingImageFilter = FunctorImageFilter<Functor::NCLSUnmixingFunctor<typename TInputImage::PixelType, typename TOutputImage::PixelType, TPrecision>>; - -} // end namespace otb - -#ifndef OTB_MANUAL_INSTANTIATION -#include "otbNCLSUnmixingImageFilter.hxx" -#endif - -#endif - diff --git a/Modules/Hyperspectral/Unmixing/include/otbNCLSUnmixingImageFilter.hxx b/Modules/Hyperspectral/Unmixing/include/otbNCLSUnmixingImageFilter.hxx deleted file mode 100644 index fc36006d0716dee9a1605fd5dea32b22f1db1da3..0000000000000000000000000000000000000000 --- a/Modules/Hyperspectral/Unmixing/include/otbNCLSUnmixingImageFilter.hxx +++ /dev/null @@ -1,112 +0,0 @@ -/* - * Copyright (C) 2005-2019 Centre National d'Etudes Spatiales (CNES) - * - * This file is part of Orfeo Toolbox - * - * https://www.orfeo-toolbox.org/ - * - * Licensed under the Apache License, Version 2.0 (the "License"); - * you may not use this file except in compliance with the License. - * You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ - -#ifndef otbNCLSUnmixingImageFilter_hxx -#define otbNCLSUnmixingImageFilter_hxx - -#include "otbNCLSUnmixingImageFilter.h" -#include <algorithm> - -namespace otb -{ - -namespace Functor -{ - -template <class TInput, class TOutput, class TPrecision> -NCLSUnmixingFunctor<TInput, TOutput, TPrecision> -::NCLSUnmixingFunctor() - : m_OutputSize(0), - m_MaxIteration(100) -{ -} - -template <class TInput, class TOutput, class TPrecision> -size_t -NCLSUnmixingFunctor<TInput, TOutput, TPrecision> -::OutputSize(const std::array<size_t,1> &) const -{ - return m_OutputSize; -} - -template <class TInput, class TOutput, class TPrecision> -void -NCLSUnmixingFunctor<TInput, TOutput, TPrecision> -::SetEndmembersMatrix(const MatrixType& U) -{ - m_U = U; - m_Ut = m_U.transpose(); - m_UtUinv = SVDType(m_Ut * m_U).inverse(); - m_OutputSize = m_U.cols(); - m_Svd.reset( new SVDType(m_U) ); -} - - -template <class TInput, class TOutput, class TPrecision> -const typename NCLSUnmixingFunctor<TInput, TOutput, TPrecision>::MatrixType& -NCLSUnmixingFunctor<TInput, TOutput, TPrecision> -::GetEndmembersMatrix() const -{ - return m_U; -} - -template <class TInput, class TOutput, class TPrecision> -typename NCLSUnmixingFunctor<TInput, TOutput, TPrecision>::OutputType -NCLSUnmixingFunctor<TInput, TOutput, TPrecision> -::operator ()(const InputType& in) const -{ - // TODO : support different types between input and output ? - VectorType inVector(in.Size()); - for (unsigned int i = 0; i < in.GetSize(); ++i ) - { - inVector[i] = in[i]; - } - - // Initialize with Unconstrained Least Square solution - VectorType uclsVector = m_Svd->solve(inVector); - - unsigned int nbEndmembers = m_OutputSize; - - // Apply NCLS iterations - VectorType lambda(nbEndmembers); - VectorType nclsVector = uclsVector; - VectorType correction(uclsVector.size()); - - for (unsigned int i = 0; i < m_MaxIteration; ++i) - { - // Error in original paper : divergence - // lambda = m_Ut * (inVector - m_U * nclsVector); - lambda = m_Ut * (m_U * nclsVector - inVector); - correction = m_UtUinv * lambda; - nclsVector -= correction; - } - - OutputType out(nclsVector.size()); - for (unsigned int i = 0; i < out.GetSize(); ++i ) - { - out[i] = nclsVector[i]; - } - return out; -} - -} // end namespace Functor -} // end namespace otb - -#endif diff --git a/Modules/Hyperspectral/Unmixing/test/CMakeLists.txt b/Modules/Hyperspectral/Unmixing/test/CMakeLists.txt index 5430a0a12ebc70f5b61b08e4f1bf3551b5dbc050..5982a7a388390d56198f4b1a30f44d29e48b82a7 100644 --- a/Modules/Hyperspectral/Unmixing/test/CMakeLists.txt +++ b/Modules/Hyperspectral/Unmixing/test/CMakeLists.txt @@ -23,7 +23,6 @@ otb_module_test() set(OTBUnmixingTests otbUnmixingTestDriver.cxx otbMDMDNMFImageFilter.cxx -otbNCLSUnmixingImageFilter.cxx otbISRAUnmixingImageFilter.cxx otbUnConstrainedLeastSquareImageFilter.cxx otbSparseUnmixingImageFilter.cxx @@ -51,16 +50,6 @@ otb_add_test(NAME hyTvMDMDNMFImageFilterTest COMMAND otbUnmixingTestDriver ${TEMP}/hyTvMDMDNMFImageFilterTest.tif 100) -otb_add_test(NAME hyTvNCLSUnmixingImageFilterTest COMMAND otbUnmixingTestDriver - --compare-image ${EPSILON_9} - ${BASELINE}/hyTvNCLSImageFilterTest.tif - ${TEMP}/hyTvNCLSImageFilterTest.tif - otbNCLSUnmixingImageFilterTest - ${INPUTDATA}/Hyperspectral/synthetic/hsi_cube.tif - ${INPUTDATA}/Hyperspectral/synthetic/endmembers.tif - ${TEMP}/hyTvNCLSImageFilterTest.tif - 10) - otb_add_test(NAME hyTvISRAUnmixingImageFilterTest COMMAND otbUnmixingTestDriver --compare-image ${EPSILON_9} ${BASELINE}/hyTvISRAUnmixingImageFilterTest.tif diff --git a/Modules/Hyperspectral/Unmixing/test/otbNCLSUnmixingImageFilter.cxx b/Modules/Hyperspectral/Unmixing/test/otbNCLSUnmixingImageFilter.cxx deleted file mode 100644 index c119b4b987a9bc002456ff513b6e57043b8a5af0..0000000000000000000000000000000000000000 --- a/Modules/Hyperspectral/Unmixing/test/otbNCLSUnmixingImageFilter.cxx +++ /dev/null @@ -1,76 +0,0 @@ -/* - * Copyright (C) 2005-2019 Centre National d'Etudes Spatiales (CNES) - * - * This file is part of Orfeo Toolbox - * - * https://www.orfeo-toolbox.org/ - * - * Licensed under the Apache License, Version 2.0 (the "License"); - * you may not use this file except in compliance with the License. - * You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ - -#include "otbNCLSUnmixingImageFilter.h" - -#include "otbVectorImage.h" -#include "otbImageFileReader.h" -#include "otbImageFileWriter.h" -#include "otbVectorImageToMatrixImageFilter.h" -#include "otbStandardWriterWatcher.h" - -const unsigned int Dimension = 2; -typedef double PixelType; - -typedef otb::VectorImage<PixelType, Dimension> ImageType; -typedef otb::ImageFileReader<ImageType> ReaderType; -typedef otb::NCLSUnmixingImageFilter<ImageType, ImageType, PixelType> UnmixingImageFilterType; -typedef otb::VectorImageToMatrixImageFilter<ImageType> VectorImageToMatrixImageFilterType; -typedef otb::ImageFileWriter<ImageType> WriterType; - -int otbNCLSUnmixingImageFilterTest(int itkNotUsed(argc), char * argv[]) -{ - const char * inputImage = argv[1]; - const char * inputEndmembers = argv[2]; - const char * outputImage = argv[3]; - int maxIter = atoi(argv[4]); - - ReaderType::Pointer readerImage = ReaderType::New(); - readerImage->SetFileName(inputImage); - - ReaderType::Pointer readerEndMembers = ReaderType::New(); - readerEndMembers->SetFileName(inputEndmembers); - VectorImageToMatrixImageFilterType::Pointer endMember2Matrix = VectorImageToMatrixImageFilterType::New(); - endMember2Matrix->SetInput(readerEndMembers->GetOutput()); - - endMember2Matrix->Update(); - - typedef VectorImageToMatrixImageFilterType::MatrixType MatrixType; - MatrixType endMembers = endMember2Matrix->GetMatrix(); - MatrixType pinv = vnl_matrix_inverse<PixelType>(endMembers); - - UnmixingImageFilterType::Pointer unmixer = UnmixingImageFilterType::New(); - - unmixer->SetInput(readerImage->GetOutput()); - unmixer->GetModifiableFunctor().SetMaxIteration(maxIter); - //unmixer->SetNumberOfThreads(1); - unmixer->GetModifiableFunctor().SetEndmembersMatrix(endMember2Matrix->GetMatrix()); - - WriterType::Pointer writer = WriterType::New(); - writer->SetFileName(outputImage); - writer->SetInput(unmixer->GetOutput()); - writer->SetNumberOfDivisionsStrippedStreaming(10); - - otb::StandardWriterWatcher w4(writer, unmixer,"NCLSUnmixingImageFilter"); - - writer->Update(); - - return EXIT_SUCCESS; -} diff --git a/Modules/Hyperspectral/Unmixing/test/otbUnmixingTestDriver.cxx b/Modules/Hyperspectral/Unmixing/test/otbUnmixingTestDriver.cxx index 22727d311a867a90dc53e01c0f650cee999cf0b2..155b26f20f03ae5421a62839375d304f4f3c7978 100644 --- a/Modules/Hyperspectral/Unmixing/test/otbUnmixingTestDriver.cxx +++ b/Modules/Hyperspectral/Unmixing/test/otbUnmixingTestDriver.cxx @@ -24,7 +24,6 @@ void RegisterTests() { REGISTER_TEST(otbMDMDNMFImageFilterTest); REGISTER_TEST(otbMDMDNMFImageFilterTest2); - REGISTER_TEST(otbNCLSUnmixingImageFilterTest); REGISTER_TEST(otbISRAUnmixingImageFilterTest); REGISTER_TEST(otbUnConstrainedLeastSquareImageFilterTest); REGISTER_TEST(otbSparseUnmixingImageFilterTest);