Commit 68d13c95 authored by Julien Malik's avatar Julien Malik

ADD: first import

parents
ADD_EXECUTABLE(otbAvirisBandSelection otbAvirisBandSelection.cxx)
TARGET_LINK_LIBRARIES(otbAvirisBandSelection OTBCommon OTBIO OTBBasicFilters )
/*=========================================================================
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 "otbMultiChannelExtractROI.h"
#include "otbVectorRescaleIntensityImageFilter.h"
const unsigned int Dimension = 2;
typedef double PixelType;
typedef otb::VectorImage<PixelType, Dimension> ImageType;
typedef otb::ImageFileReader<ImageType> ReaderType;
typedef otb::MultiChannelExtractROI<PixelType,PixelType> MultiChannelExtractROIType;
typedef otb::VectorRescaleIntensityImageFilter<ImageType,ImageType> RescaleImageFilterType;
typedef otb::ImageFileWriter<ImageType> WriterType;
int main(int argc, char * argv[])
{
const char * infname = argv[1];
const char * outfname = argv[2];
ReaderType::Pointer reader = ReaderType::New();
reader->SetFileName(infname);
MultiChannelExtractROIType::Pointer extractROI = MultiChannelExtractROIType::New();
extractROI->SetInput(reader->GetOutput());
unsigned int i = 0;
for (i = 5; i < 107; i++)
{
extractROI->SetChannel(i);
}
for (i = 122; i < 152; i++)
{
extractROI->SetChannel(i);
}
for (i = 182; i < 202; i++)
{
extractROI->SetChannel(i);
}
//extractROI->GenerateOutputInformation();
RescaleImageFilterType::Pointer scale = RescaleImageFilterType::New();
scale->SetInput(extractROI->GetOutput());
ImageType::PixelType inputMin;
ImageType::PixelType inputMax;
ImageType::PixelType outputMin;
ImageType::PixelType outputMax;
//const unsigned int nbBands = extractROI->GetOutput()->GetNumberOfComponentsPerPixel();
const unsigned int nbBands = 152;
scale->SetAutomaticInputMinMaxComputation(false);
inputMin.SetSize(nbBands);
inputMin.Fill(0.0);
scale->SetInputMinimum(inputMin);
inputMax.SetSize(nbBands);
inputMax.Fill(100000.0);
scale->SetInputMaximum(inputMax);
outputMin.SetSize(nbBands);
outputMin.Fill(0.0);
scale->SetOutputMinimum(outputMin);
outputMax.SetSize(nbBands);
outputMax.Fill(10.0);
scale->SetOutputMaximum(outputMax);
WriterType::Pointer writer = WriterType::New();
writer->SetFileName(outfname);
writer->SetInput(scale->GetOutput());
writer->Update();
return EXIT_SUCCESS;
}
PROJECT(Hyperspectral)
cmake_minimum_required(VERSION 2.6)
IF(NOT LIBRARY_OUTPUT_PATH)
SET (LIBRARY_OUTPUT_PATH ${Hyperspectral_BINARY_DIR}/bin CACHE INTERNAL "Single output directory for building all libraries." FORCE)
ENDIF(NOT LIBRARY_OUTPUT_PATH)
IF(NOT EXECUTABLE_OUTPUT_PATH)
SET (EXECUTABLE_OUTPUT_PATH ${Hyperspectral_BINARY_DIR}/bin CACHE INTERNAL "Single output directory for building all executables." FORCE)
ENDIF(NOT EXECUTABLE_OUTPUT_PATH)
INCLUDE_REGULAR_EXPRESSION("^.*$")
# Find OTB and load its settings.
FIND_PACKAGE(OTB)
IF(OTB_FOUND)
INCLUDE(${OTB_USE_FILE})
ELSE(OTB_FOUND)
MESSAGE(FATAL_ERROR "OTB not found. Please set OTB_DIR")
ENDIF(OTB_FOUND)
INCLUDE_DIRECTORIES(
${Hyperspectral_SOURCE_DIR}/Code/BasicFilters
${Hyperspectral_SOURCE_DIR}/Code/Hyperspectral
${Hyperspectral_SOURCE_DIR}/Code/Vahine
)
ADD_SUBDIRECTORY(Code)
ADD_SUBDIRECTORY(Application)
ENABLE_TESTING()
ADD_SUBDIRECTORY(Testing)
## This file should be placed in the root directory of your project.
## Then modify the CMakeLists.txt file in the root directory of your
## project to incorporate the testing dashboard.
## # The following are required to uses Dart and the Cdash dashboard
## ENABLE_TESTING()
## INCLUDE(Dart)
set(CTEST_PROJECT_NAME "Hyper")
#set(CTEST_NIGHTLY_START_TIME "20:00:00 CEST")
#set(CTEST_NIGHTLY_IDENT_LOCATION "http://host2.orfeo-toolbox.org/nightly/libNightlyNumber")
#set(CTEST_DROP_METHOD "http")
#set(CTEST_DROP_SITE "dash.orfeo-toolbox.org")
#set(CTEST_DROP_LOCATION "/submit.php?project=OTB")
#set(CTEST_DROP_SITE_CDASH TRUE)
This diff is collapsed.
/*=========================================================================
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 __otbStreamingStatisticsVectorImageFilter_txx
#define __otbStreamingStatisticsVectorImageFilter_txx
#include "otbStreamingStatisticsVectorImageFilter2.h"
#include "itkImageRegionIterator.h"
#include "itkImageRegionConstIteratorWithIndex.h"
#include "itkNumericTraits.h"
#include "itkProgressReporter.h"
#include "otbMacro.h"
namespace otb
{
template<class TInputImage>
PersistentStreamingStatisticsVectorImageFilter2<TInputImage>
::PersistentStreamingStatisticsVectorImageFilter2()
{
// first output is a copy of the image, DataObject created by
// superclass
//
// allocate the data objects for the outputs which are
// just decorators around vector/matric types
this->itk::ProcessObject::SetNthOutput(1, this->MakeOutput(1).GetPointer());
this->itk::ProcessObject::SetNthOutput(2, this->MakeOutput(2).GetPointer());
this->itk::ProcessObject::SetNthOutput(3, this->MakeOutput(3).GetPointer());
}
template<class TInputImage>
itk::DataObject::Pointer
PersistentStreamingStatisticsVectorImageFilter2<TInputImage>
::MakeOutput(unsigned int output)
{
switch (output)
{
case 0:
return static_cast<itk::DataObject*>(TInputImage::New().GetPointer());
break;
case 1:
return static_cast<itk::DataObject*>(RealPixelObjectType::New().GetPointer());
break;
case 2:
case 3:
return static_cast<itk::DataObject*>(MatrixObjectType::New().GetPointer());
break;
default:
// might as well make an image
return static_cast<itk::DataObject*>(TInputImage::New().GetPointer());
break;
}
}
template<class TInputImage>
typename PersistentStreamingStatisticsVectorImageFilter2<TInputImage>::RealPixelObjectType*
PersistentStreamingStatisticsVectorImageFilter2<TInputImage>
::GetMeanOutput()
{
return static_cast<RealPixelObjectType*>(this->itk::ProcessObject::GetOutput(1));
}
template<class TInputImage>
const typename PersistentStreamingStatisticsVectorImageFilter2<TInputImage>::RealPixelObjectType*
PersistentStreamingStatisticsVectorImageFilter2<TInputImage>
::GetMeanOutput() const
{
return static_cast<const RealPixelObjectType*>(this->itk::ProcessObject::GetOutput(1));
}
template<class TInputImage>
typename PersistentStreamingStatisticsVectorImageFilter2<TInputImage>::MatrixObjectType*
PersistentStreamingStatisticsVectorImageFilter2<TInputImage>
::GetCorrelationOutput()
{
return static_cast<MatrixObjectType*>(this->itk::ProcessObject::GetOutput(2));
}
template<class TInputImage>
const typename PersistentStreamingStatisticsVectorImageFilter2<TInputImage>::MatrixObjectType*
PersistentStreamingStatisticsVectorImageFilter2<TInputImage>
::GetCorrelationOutput() const
{
return static_cast<const MatrixObjectType*>(this->itk::ProcessObject::GetOutput(2));
}
template<class TInputImage>
typename PersistentStreamingStatisticsVectorImageFilter2<TInputImage>::MatrixObjectType*
PersistentStreamingStatisticsVectorImageFilter2<TInputImage>
::GetCovarianceOutput()
{
return static_cast<MatrixObjectType*>(this->itk::ProcessObject::GetOutput(3));
}
template<class TInputImage>
const typename PersistentStreamingStatisticsVectorImageFilter2<TInputImage>::MatrixObjectType*
PersistentStreamingStatisticsVectorImageFilter2<TInputImage>
::GetCovarianceOutput() const
{
return static_cast<const MatrixObjectType*>(this->itk::ProcessObject::GetOutput(3));
}
template<class TInputImage>
void
PersistentStreamingStatisticsVectorImageFilter2<TInputImage>
::GenerateOutputInformation()
{
Superclass::GenerateOutputInformation();
if (this->GetInput())
{
this->GetOutput()->CopyInformation(this->GetInput());
this->GetOutput()->SetLargestPossibleRegion(this->GetInput()->GetLargestPossibleRegion());
if (this->GetOutput()->GetRequestedRegion().GetNumberOfPixels() == 0)
{
this->GetOutput()->SetRequestedRegion(this->GetOutput()->GetLargestPossibleRegion());
}
}
}
template<class TInputImage>
void
PersistentStreamingStatisticsVectorImageFilter2<TInputImage>
::AllocateOutputs()
{
// This is commented to prevent the streaming of the whole image for the first stream strip
// It shall not cause any problem because the output image of this filter is not intended to be used.
//InputImagePointer image = const_cast< TInputImage * >( this->GetInput() );
//this->GraftOutput( image );
// Nothing that needs to be allocated for the remaining outputs
}
template<class TInputImage>
void
PersistentStreamingStatisticsVectorImageFilter2<TInputImage>
::Reset()
{
TInputImage * inputPtr = const_cast<TInputImage *>(this->GetInput());
inputPtr->UpdateOutputInformation();
unsigned int numberOfThreads = this->GetNumberOfThreads();
unsigned int numberOfComponent = inputPtr->GetNumberOfComponentsPerPixel();
RealPixelType zeroRealPixel;
zeroRealPixel.SetSize(numberOfComponent);
zeroRealPixel.Fill(itk::NumericTraits<RealType>::Zero);
this->GetMeanOutput()->Set(zeroRealPixel);
MatrixType zeroMatrix;
zeroMatrix.SetSize(numberOfComponent, numberOfComponent);
zeroMatrix.Fill(itk::NumericTraits<RealType>::Zero);
this->GetCovarianceOutput()->Set(zeroMatrix);
this->GetCorrelationOutput()->Set(zeroMatrix);
m_FirstOrderAccumulators.resize(numberOfThreads);
std::fill(m_FirstOrderAccumulators.begin(), m_FirstOrderAccumulators.end(), zeroRealPixel);
m_SecondOrderAccumulators.resize(numberOfThreads);
std::fill(m_SecondOrderAccumulators.begin(), m_SecondOrderAccumulators.end(), zeroMatrix);
}
template<class TInputImage>
void
PersistentStreamingStatisticsVectorImageFilter2<TInputImage>
::Synthetize()
{
TInputImage * inputPtr = const_cast<TInputImage *>(this->GetInput());
const unsigned int nbPixels = inputPtr->GetLargestPossibleRegion().GetNumberOfPixels();
const unsigned int numberOfComponent = inputPtr->GetNumberOfComponentsPerPixel();
RealPixelType streamFirstOrderAccumulator(numberOfComponent);
streamFirstOrderAccumulator.Fill(itk::NumericTraits<RealType>::Zero);
MatrixType streamSecondOrderAccumulator(numberOfComponent, numberOfComponent);
streamSecondOrderAccumulator.Fill(itk::NumericTraits<RealType>::Zero);
const unsigned int numberOfThreads = this->GetNumberOfThreads();
for (unsigned int i = 0; i < numberOfThreads; ++i)
{
streamFirstOrderAccumulator += m_FirstOrderAccumulators [i];
streamSecondOrderAccumulator += m_SecondOrderAccumulators[i];
}
const double regul = static_cast<double>(nbPixels) / (nbPixels - 1);
this->GetMeanOutput()->Set(streamFirstOrderAccumulator / nbPixels);
const RealPixelType& mean = this->GetMeanOutput()->Get();
MatrixType cor = streamSecondOrderAccumulator / nbPixels;
this->GetCorrelationOutput()->Set(cor);
MatrixType cov = cor;
for (unsigned int r = 0; r < numberOfComponent; ++r)
{
for (unsigned int c = 0; c < numberOfComponent; ++c)
{
cov(r, c) = regul * (cov(r, c) - mean[r] * mean[c]);
}
}
this->GetCovarianceOutput()->Set(cov);
}
template<class TInputImage>
void
PersistentStreamingStatisticsVectorImageFilter2<TInputImage>
::ThreadedGenerateData(const RegionType& outputRegionForThread, int threadId)
{
// Support progress methods/callbacks
itk::ProgressReporter progress(this, threadId, outputRegionForThread.GetNumberOfPixels());
// Grab the input
InputImagePointer inputPtr = const_cast<TInputImage *>(this->GetInput());
RealPixelType& threadFirstOrder = m_FirstOrderAccumulators [threadId];
MatrixType& threadSecondOrder = m_SecondOrderAccumulators[threadId];
itk::ImageRegionConstIteratorWithIndex<TInputImage> it(inputPtr, outputRegionForThread);
MatrixType pixelCorr;
for (it.GoToBegin(); !it.IsAtEnd(); ++it, progress.CompletedPixel())
{
PixelType vectorValue = it.Get();
threadFirstOrder += vectorValue;
for (unsigned int r = 0; r < threadSecondOrder.Rows(); ++r)
{
for (unsigned int c = 0; c < threadSecondOrder.Cols(); ++c)
{
threadSecondOrder(r,c) += vectorValue[r] * vectorValue[c];
}
}
}
}
template <class TImage>
void
PersistentStreamingStatisticsVectorImageFilter2<TImage>
::PrintSelf(std::ostream& os, itk::Indent indent) const
{
Superclass::PrintSelf(os, indent);
os << indent << "Mean: " << this->GetMeanOutput()->Get() << std::endl;
os << indent << "Covariance: " << this->GetCovarianceOutput()->Get() << std::endl;
os << indent << "Correlation: " << this->GetCorrelationOutput()->Get() << std::endl;
}
} // end namespace otb
#endif
/*=========================================================================
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 MERCHANT2ABILITY or FITNESS FOR A PARTICULAR
PURPOSE. See the above copyright notices for more information.
=========================================================================*/
#ifndef __otbVectorImageToMatrixFilter_h
#define __otbVectorImageToMatrixFilter_h
#include "otbPersistentImageFilter.h"
#include "otbPersistentFilterStreamingDecorator.h"
#include "itkNumericTraits.h"
#include "itkArray.h"
#include "itkSimpleDataObjectDecorator.h"
#include "itkImageRegionSplitter.h"
#include "itkVariableSizeMatrix.h"
#include "itkVariableLengthVector.h"
namespace otb
{
/** \class PersistentVectorImageToMatrixFilter
* \brief Compute covariance & correlation of a large image using streaming
*
* This filter persists its temporary data. It means that if you Update it n times on n different
* requested regions, the output statistics will be the statitics of the whole set of n regions.
*
* To reset the temporary data, one should call the Reset() function.
*
* To get the statistics once the regions have been processed via the pipeline, use the Synthetize() method.
*
* \sa PersistentImageFilter
* \ingroup Streamed
* \ingroup Multithreaded
* \ingroup MathematicalStatisticsImageFilters
*
*/
template<class TInputImage>
class ITK_EXPORT PersistentVectorImageToMatrixFilter :
public PersistentImageFilter<TInputImage, TInputImage>
{
public:
/** Standard Self typedef */
typedef PersistentVectorImageToMatrixFilter Self;
typedef PersistentImageFilter<TInputImage, TInputImage> Superclass;
typedef itk::SmartPointer<Self> Pointer;
typedef itk::SmartPointer<const Self> ConstPointer;
/** Method for creation through the object factory. */
itkNewMacro(Self);
/** Runtime information support. */
itkTypeMacro(PersistentVectorImageToMatrixFilter, PersistentImageFilter);
/** Image related typedefs. */
typedef TInputImage ImageType;
typedef typename TInputImage::Pointer InputImagePointer;
typedef typename TInputImage::RegionType RegionType;
typedef typename TInputImage::SizeType SizeType;
typedef typename TInputImage::IndexType IndexType;
typedef typename TInputImage::PixelType PixelType;
typedef typename TInputImage::InternalPixelType InternalPixelType;
itkStaticConstMacro(InputImageDimension, unsigned int, TInputImage::ImageDimension);
/** Image related typedefs. */
itkStaticConstMacro(ImageDimension, unsigned int, TInputImage::ImageDimension);
/** Type to use for computations. */
typedef typename itk::NumericTraits<InternalPixelType>::RealType RealType;
typedef itk::VariableLengthVector<RealType> RealPixelType;
/** Smart Pointer type to a DataObject. */
typedef typename itk::DataObject::Pointer DataObjectPointer;
/** Type of DataObjects used for scalar outputs */
typedef itk::VariableSizeMatrix<RealType> MatrixType;
typedef itk::SimpleDataObjectDecorator<MatrixType> MatrixObjectType;
/** Return the computed matrix as a DataObject */
MatrixObjectType* GetMatrixOutput();
const MatrixObjectType* GetMatrixOutput() const;
/** Make a DataObject of the correct type to be used as the specified
* output.
*/
virtual DataObjectPointer MakeOutput(unsigned int idx);
/** Pass the input through unmodified. Do this by Grafting in the
* AllocateOutputs method.
*/
virtual void AllocateOutputs();
virtual void GenerateOutputInformation();
virtual void Synthetize(void);
virtual void Reset(void);
protected:
PersistentVectorImageToMatrixFilter();
virtual ~PersistentVectorImageToMatrixFilter() {}
virtual void PrintSelf(std::ostream& os, itk::Indent indent) const;
void ThreadedGenerateData(const RegionType& outputRegionForThread, int threadId);
private:
PersistentVectorImageToMatrixFilter(const Self &); //purposely not implemented
void operator =(const Self&); //purposely not implemented
}; // end of class PersistentVectorImageToMatrixFilter
/**===========================================================================*/
/** \class VectorImageToMatrixImageFilter
* \brief This class streams the whole input image through the PersistentStatisticsImageFilter.
*
* This way, it allows to compute the first order global statistics of this image. It calls the
* Reset() method of the PersistentStatisticsImageFilter before streaming the image and the
* Synthetize() method of the PersistentStatisticsImageFilter after having streamed the image
* to compute the statistics. The accessor on the results are wrapping the accessors of the
* internal PersistentStatisticsImageFilter.
*
* \sa PersistentVectorImageToMatrixFilter
* \sa PersistentImageFilter
* \sa PersistentFilterStreamingDecorator
* \sa StreamingImageVirtualWriter
* \ingroup Streamed
* \ingroup Multithreaded
* \ingroup MathematicalStatisticsImageFilters
*/
template<class TInputImage>
class ITK_EXPORT VectorImageToMatrixImageFilter :
public PersistentFilterStreamingDecorator<PersistentVectorImageToMatrixFilter<TInputImage> >
{
public:
/** Standard Self typedef */
typedef VectorImageToMatrixImageFilter Self;
typedef PersistentFilterStreamingDecorator
<PersistentVectorImageToMatrixFilter<TInputImage> > Superclass;
typedef itk::SmartPointer<Self> Pointer;
typedef itk::SmartPointer<const Self> ConstPointer;
/** Type macro */
itkNewMacro(Self);
/** Creation through object factory macro */
itkTypeMacro(VectorImageToMatrixImageFilter, PersistentFilterStreamingDecorator);
typedef TInputImage InputImageType;
typedef typename Superclass::FilterType PersistentFilterType;
typedef typename PersistentFilterType::MatrixType MatrixType;
typedef typename PersistentFilterType::MatrixObjectType MatrixObjectType;
void SetInput(InputImageType * input)
{
this->GetFilter()->SetInput(input);
}
const InputImageType * GetInput()
{
return this->GetFilter()->GetInput();
}
/** Return the computed Mean. */
const MatrixType& GetMatrix() const
{
return this->GetMatrixOutput()->Get();
}
MatrixObjectType* GetMatrixOutput()
{
return this->GetFilter()->GetMatrixOutput();
}
const MatrixObjectType* GetMatrixOutput() const
{
return this->GetFilter()->GetMatrixOutput();
}
protected:
/** Constructor */
VectorImageToMatrixImageFilter() {};
/** Destructor */
virtual ~VectorImageToMatrixImageFilter() {}
private:
VectorImageToMatrixImageFilter(const Self &); //purposely not implemented
void operator =(const Self&); //purposely not implemented
};
} // end namespace otb
#ifndef OTB_MANUAL_INSTANTIATION
#include "otbVectorImageToMatrixImageFilter.txx"
#endif
#endif
/*=========================================================================
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 __otbVectorImageToMatrixFilter_txx
#define __otbVectorImageToMatrixFilter_txx
#include "otbVectorImageToMatrixImageFilter.h"
#include "itkImageRegionIterator.h"
#include "itkImageRegionConstIteratorWithIndex.h"
#include "itkNumericTraits.h"
#include "itkProgressReporter.h"
#include "otbMacro.h"
namespace otb
{
template<class TInputImage>
PersistentVectorImageToMatrixFilter<TInputImage>
::PersistentVectorImageToMatrixFilter()
{
// first output is a copy of the image, DataObject created by
// superclass
//
// allocate the data objects for the outputs which are
// just decorators around vector/matric types
std::cout << "PersistentVectorImageToMatrixFilter" << std::endl;
this->itk::ProcessObject::SetNthOutput(1, this->MakeOutput(1).GetPointer());
}
template<class TInputImage>
itk::DataObject::Pointer
PersistentVectorImageToMatrixFilter<TInputImage>
::MakeOutput(unsigned int output)
{
std::cout << "MakeOutput" << output << std::endl;
switch (output)
{
case 0:
return static_cast<itk::DataObject*>(TInputImage::New().GetPointer());
break;
case 1:
return static_cast<itk::DataObject*>(MatrixObjectType::New().GetPointer());
break;
default:
// might as well make an image
return static_cast<itk::DataObject*>(TInputImage::New().GetPointer());
break;
}
}
template<class TInputImage>
typename PersistentVectorImageToMatrixFilter<TInputImage>::MatrixObjectType*
PersistentVectorImageToMatrixFilter<TInputImage>
::GetMatrixOutput()
{
return static_cast<MatrixObjectType*>(this->itk::ProcessObject::GetOutput(1));
}
template<class TInputImage>
const typename PersistentVectorImageToMatrixFilter<TInputImage>::MatrixObjectType*
PersistentVectorImageToMatrixFilter<TInputImage>
::GetMatrixOutput() const
{