Commit f3adf4a9 authored by Manuel Grizonnet's avatar Manuel Grizonnet

WPI:add an example to perform unmixing

parent d754b8ca
......@@ -26,6 +26,7 @@ ADD_SUBDIRECTORY(Tutorials)
ADD_SUBDIRECTORY(Markov)
ADD_SUBDIRECTORY(OBIA)
ADD_SUBDIRECTORY(Simulation)
ADD_SUBDIRECTORY(Hyperspectral)
IF(OTB_USE_VISU_GUI)
ADD_SUBDIRECTORY(Visualization)
......@@ -87,13 +88,14 @@ ELSE(OTB_BINARY_DIR)
ADD_SUBDIRECTORY(Markov)
ADD_SUBDIRECTORY(OBIA)
ADD_SUBDIRECTORY(Simulation)
ADD_SUBDIRECTORY(Hyperspectral)
IF(OTB_USE_VISU_GUI)
ADD_SUBDIRECTORY(Visualization)
ENDIF(OTB_USE_VISU_GUI)
IF(OTB_USE_PATENTED)
ADD_SUBDIRECTORY( Patented )
ADD_SUBDIRECTORY(Patented)
ENDIF(OTB_USE_PATENTED)
IF(OTB_USE_PQXX)
......
#
# Examples on the use of hyperspectral algorithms
#
PROJECT(HyperspectralExamples)
INCLUDE_REGULAR_EXPRESSION("^.*$")
ADD_EXECUTABLE(Unmixing HyperspectralUnmixing.cxx )
TARGET_LINK_LIBRARIES(Unmixing OTBCommon OTBIO OTBHyperspectral)
IF( NOT OTB_DISABLE_CXX_TESTING AND BUILD_TESTING )
SET(BASELINE ${OTB_DATA_ROOT}/Baseline/Examples/Hyperspctral)
SET(BASELINE_IMAGES_OTB ${OTB_DATA_ROOT}/Baseline/OTB/Images)
SET(INPUTDATA_EXAMPLES ${OTB_DATA_ROOT}/Examples)
SET(INPUTDATA ${OTB_DATA_ROOT}/Input)
SET(TEMP ${OTB_BINARY_DIR}/Testing/Temporary)
SET(EXE_TESTS ${CXX_TEST_PATH}/otbHyperspectralExamplesTests)
SET(TOL 0.0)
# tests#
# ------- ProsailModelExampleTest ----------
ADD_TEST(siTvHyperspectralExampleTest ${EXE_TESTS}
#--compare-image ${TOL} ${BASELINE}/siTvUnmixingExampleTest.txt
# ${TEMP}/siTvUnmixingExampleTest.txt
Unmixing
#${INPUTDATA}/Hyperspectral/synthetic/hsi_cube.tif
${INPUTDATA}/cupriteSubHsi.tif
${TEMP}/siTvHyperspectralUnmixingExampleTest.tif
)
INCLUDE_DIRECTORIES(${OTB_SOURCE_DIR}/Testing/Code)
ADD_EXECUTABLE(otbHyperspectralExamplesTests otbHyperspectralExamplesTests.cxx)
TARGET_LINK_LIBRARIES(otbHyperspectralExamplesTests OTBCommon OTBIO OTBHyperspectral OTBCommon OTBTesting)
ENDIF( NOT OTB_DISABLE_CXX_TESTING AND BUILD_TESTING )
#
# Examples on the use of hyperspectral algorithms
#
PROJECT(HyperspectralExamples)
INCLUDE_REGULAR_EXPRESSION("^.*$")
ADD_EXECUTABLE(Unmixing HyperspectralUnmixing.cxx )
TARGET_LINK_LIBRARIES(Unmixing OTBCommon OTBIO OTBHyperspectral)
IF( NOT OTB_DISABLE_CXX_TESTING AND BUILD_TESTING )
SET(BASELINE ${OTB_DATA_ROOT}/Baseline/Examples/Hyperspctral)
SET(BASELINE_IMAGES_OTB ${OTB_DATA_ROOT}/Baseline/OTB/Images)
SET(INPUTDATA_EXAMPLES ${OTB_DATA_ROOT}/Examples)
SET(INPUTDATA ${OTB_DATA_ROOT}/Input)
SET(TEMP ${OTB_BINARY_DIR}/Testing/Temporary)
SET(EXE_TESTS ${CXX_TEST_PATH}/otbSimulationExamplesTests)
SET(TOL 0.0)
# tests#
# ------- ProsailModelExampleTest ----------
ADD_TEST(siTvHyperspectralExampleTest ${EXE_TESTS}
--compare-image ${TOL} ${BASELINE}/siTvUnmixingExampleTest.txt
${TEMP}/siTvUnmixingExampleTest.txt
Unmixing
${INPUTDATA}/cuprite.tif
${TEMP}/siTvUnmixingExampleTest.txt
)
INCLUDE_DIRECTORIES(${OTB_SOURCE_DIR}/Testing/Code)
ADD_EXECUTABLE(otbHyperspectralExamplesTests otbHyperspectralExamplesTests.cxx)
TARGET_LINK_LIBRARIES(otbHyperspectralExamplesTests OTBCommon OTBIO OTBHyperspectral OTBCommon OTBTesting)
ENDIF( NOT OTB_DISABLE_CXX_TESTING AND BUILD_TESTING )
/*=========================================================================
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.
=========================================================================*/
#include "otbVectorImage.h"
#include "otbImageFileReader.h"
#include "otbImageFileWriter.h"
#include "otbStreamingStatisticsVectorImageFilter.h"
#include "otbEigenvalueLikelihoodMaximisation.h"
#include "otbVcaImageFilter.h"
#include "otbUnConstrainedLeastSquareImageFilter.h"
#include "otbVectorImageToMatrixImageFilter.h"
int main(int argc, char * argv[])
{
if (argc != 3)
{
std::cerr << "Usage: " << argv[0];
std::cerr << "inputFileName outputFileName " << std::endl;
return EXIT_FAILURE;
}
const char * infname = argv[1];
const char * outfname = argv[2];
const unsigned int Dimension = 2;
typedef double PixelType;
typedef otb::VectorImage<PixelType, Dimension> ImageType;
typedef otb::ImageFileReader<ImageType> ReaderType;
typedef otb::ImageFileWriter<FloatImageType> WriterType;
typedef otb::StreamingStatisticsVectorImageFilter<ImageType> StreamingStatisticsVectorImageFilterType;
typedef otb::EigenvalueLikelihoodMaximisation<double> ELMType;
typedef otb::VCAImageFilter<ImageType> VCAFilterType;
typedef otb::UnConstrainedLeastSquareImageFilter<ImageType,ImageType,double> UCLSUnmixingFilterType;
typedef otb::VectorImageToMatrixImageFilter<ImageType> VectorImageToMatrixImageFilterType;
typedef vnl_vector<double> VectorType;
typedef vnl_matrix<double> MatrixType;
/// / Noise filtering
// typedef otb::LocalActivityVectorImageFilter< ImageType, ImageType > NoiseFilterType;
// // Image filtering
// typedef otb::MNFImageFilter< ImageType, ImageType,
// NoiseFilterType, otb::Transform::FORWARD > FilterType;
//read image
ReaderType::Pointer reader = ReaderType::New();
reader->SetFileName(infname);
reader->UpdateOutputInformation();
std::cout << "Computing stats" << std::endl;
StreamingStatisticsVectorImageFilterType::Pointer stats = StreamingStatisticsVectorImageFilterType::New();
stats->SetInput(reader->GetOutput());
otb::StandardWriterWatcher watcher(stats->GetStreamer(), stats->GetFilter(), "Computing image stats");
stats->Update();
VectorType mean (stats->GetMean().GetDataPointer(), stats->GetMean().GetSize());
MatrixType covariance = stats->GetCovariance().GetVnlMatrix();
MatrixType correlation = stats->GetCorrelation().GetVnlMatrix();
/*
* Estimate Endmembers Numbers
*/
std::cout << "Estimate Endmembers Numbers by ELM" << std::endl;
ELMType::Pointer elm = ELMType::New();
elm->SetCovariance(covariance);
elm->SetCorrelation(correlation);
elm->SetNumberOfPixels(reader->GetOutput()->GetLargestPossibleRegion().GetNumberOfPixels());
elm->Compute();
const unsigned int nbEndmembers = elm->GetNumberOfEndmembers();
std::cout << "ELM : " << nbEndmembers << " estimated endmembers" << std::endl;
std::cout << "ELM likelihood values: " << elm->GetLikelihood() << std::endl;
/*
* Estimate Endmembers
*/
std::cout << "Estimate Endmembers by VCA" << std::endl;
VCAFilterType::Pointer vca = VCAFilterType::New();
vca->SetNumberOfEndmembers(nbEndmembers);
vca->SetInput(reader->GetOutput());
ImageType::Pointer endmembersImage;
endmembersImage = vca->GetOutput();
/*
* Transform Endmembers image to matrix representation
*/
std::cout << "Endmembers extracted" << std::endl;
std::cout << "Converting endmembers to matrix" << std::endl;
VectorImageToMatrixImageFilterType::Pointer endMember2Matrix = VectorImageToMatrixImageFilterType::New();
endMember2Matrix->SetInput(endmembersImage);
endMember2Matrix->Update();
MatrixType endMembersMatrix = endMember2Matrix->GetMatrix();
std::cout << "Endmembers matrix : " << endMembersMatrix << std::endl;
/*
* Unmix
*/
UCLSUnmixingFilterType::Pointer unmixer =
UCLSUnmixingFilterType::New();
unmixer->SetInput(reader->GetOutput());
unmixer->SetMatrix(endMembersMatrix);
unmixer->SetNumberOfThreads(1); // FIXME : currently buggy
ImageType::Pointer abundanceMap;
abundanceMap = unmixer->GetOutput();
/*
* Write abundance map
*/
std::cout << "Write abundance map " << outfname << std::endl;
WriterType::Pointer writer = WriterType::New();
writer->SetInput(abundanceMap);
writer->SetFileName(outfname);
writer->Update();
return EXIT_SUCCESS;
} // end main
/*=========================================================================
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.
=========================================================================*/
#include "otbVectorImage.h"
#include "otbImageFileReader.h"
#include "otbImageFileWriter.h"
#include "otbStreamingStatisticsVectorImageFilter.h"
#include "otbEigenvalueLikelihoodMaximisation.h"
#include "otbVcaImageFilter.h"
#include "otbUnConstrainedLeastSquareImageFilter.h"
#include "otbVectorImageToMatrixImageFilter.h"
#include "itkCastImageFilter.h"
int main(int argc, char * argv[])
{
if (argc != 3)
{
std::cerr << "Usage: " << argv[0];
std::cerr << "inputFileName outputFileName " << std::endl;
return EXIT_FAILURE;
}
const char * infname = argv[1];
const char * outfname = argv[2];
const unsigned int Dimension = 2;
typedef double PixelType;
typedef otb::VectorImage<PixelType, Dimension> ImageType;
typedef otb::VectorImage<float, Dimension> FloatImageType;
typedef itk::CastImageFilter<ImageType, FloatImageType> DoubleToFloatFilterType;
typedef otb::ImageFileReader<ImageType> ReaderType;
typedef otb::ImageFileWriter<FloatImageType> WriterType;
typedef otb::StreamingStatisticsVectorImageFilter<ImageType> StreamingStatisticsVectorImageFilterType;
typedef otb::EigenvalueLikelihoodMaximisation<double> ELMType;
typedef otb::VCAImageFilter<ImageType> VCAFilterType;
typedef otb::UnConstrainedLeastSquareImageFilter<ImageType,ImageType,double> UCLSUnmixingFilterType;
typedef otb::VectorImageToMatrixImageFilter<ImageType> VectorImageToMatrixImageFilterType;
typedef vnl_vector<double> VectorType;
typedef vnl_matrix<double> MatrixType;
/// / Noise filtering
// typedef otb::LocalActivityVectorImageFilter< ImageType, ImageType > NoiseFilterType;
// // Image filtering
// typedef otb::MNFImageFilter< ImageType, ImageType,
// NoiseFilterType, otb::Transform::FORWARD > FilterType;
//read image
ReaderType::Pointer reader = ReaderType::New();
reader->SetFileName(infname);
reader->UpdateOutputInformation();
std::cout << "Computing stats" << std::endl;
StreamingStatisticsVectorImageFilterType::Pointer stats = StreamingStatisticsVectorImageFilterType::New();
stats->SetInput(reader->GetOutput());
otb::StandardWriterWatcher watcher(stats->GetStreamer(), stats->GetFilter(), "Computing image stats");
stats->Update();
VectorType mean (stats->GetMean().GetDataPointer(), stats->GetMean().GetSize());
MatrixType covariance = stats->GetCovariance().GetVnlMatrix();
MatrixType correlation = stats->GetCorrelation().GetVnlMatrix();
/*
* Estimate Endmembers Numbers
*/
std::cout << "Estimate Endmembers Numbers by ELM" << std::endl;
ELMType::Pointer elm = ELMType::New();
elm->SetCovariance(covariance);
elm->SetCorrelation(correlation);
elm->SetNumberOfPixels(reader->GetOutput()->GetLargestPossibleRegion().GetNumberOfPixels());
elm->Compute();
const unsigned int nbEndmembers = elm->GetNumberOfEndmembers();
std::cout << "ELM : " << nbEndmembers << " estimated endmembers" << std::endl;
std::cout << "ELM likelihood values: " << elm->GetLikelihood() << std::endl;
/*
* Estimate Endmembers
*/
std::cout << "Estimate Endmembers by VCA" << std::endl;
VCAFilterType::Pointer vca = VCAFilterType::New();
vca->SetNumberOfEndmembers(nbEndmembers);
vca->SetInput(reader->GetOutput());
ImageType::Pointer endmembersImage;
endmembersImage = vca->GetOutput();
/*
* Transform Endmembers image to matrix representation
*/
std::cout << "Endmembers extracted" << std::endl;
std::cout << "Converting endmembers to matrix" << std::endl;
VectorImageToMatrixImageFilterType::Pointer endMember2Matrix = VectorImageToMatrixImageFilterType::New();
endMember2Matrix->SetInput(endmembersImage);
endMember2Matrix->Update();
MatrixType endMembersMatrix = endMember2Matrix->GetMatrix();
std::cout << "Endmembers matrix : " << endMembersMatrix << std::endl;
/*
* Unmix
*/
UCLSUnmixingFilterType::Pointer unmixer =
UCLSUnmixingFilterType::New();
unmixer->SetInput(reader->GetOutput());
unmixer->SetMatrix(endMembersMatrix);
unmixer->SetNumberOfThreads(1); // FIXME : currently buggy
ImageType::Pointer abundanceMap;
abundanceMap = unmixer->GetOutput();
/*
* Write abundance map
*/
//std::cout << "Write abundance map" << outputImageName << std::endl;
DoubleToFloatFilterType::Pointer cast = DoubleToFloatFilterType::New();
cast->SetInput( abundanceMap );
FloatImageType::Pointer abundanceMapFloat = cast->GetOutput();
WriterType::Pointer writer = WriterType::New();
writer->SetInput(abundanceMapFloat);
writer->SetFileName(outfname);
writer->Update();
return EXIT_SUCCESS;
} // end main
/*=========================================================================
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 otbMultiScaleTest for the test driver
// and all it expects is that you have a function called RegisterTests
#include <iostream>
#include "otbTestMain.h"
void RegisterTests()
{
REGISTER_TEST(ProsailModel);
REGISTER_TEST(LAIFromNDVIImageTransform);
REGISTER_TEST(LAIAndPROSAILToSensorResponse);
}
#undef main
#define main ProsailModel
#include "ProsailModel.cxx"
#undef main
#define main LAIFromNDVIImageTransform
#include "LAIFromNDVIImageTransform.cxx"
#undef main
#define main LAIAndPROSAILToSensorResponse
#include "LAIAndPROSAILToSensorResponse.cxx"
/*=========================================================================
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 otbMultiScaleTest for the test driver
// and all it expects is that you have a function called RegisterTests
#include <iostream>
#include "otbTestMain.h"
void RegisterTests()
{
REGISTER_TEST(Unmixing);
}
#undef main
#define main Unmixing
#include "HyperspectralUnmixing.cxx"
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment