Commit 6cd21bd0 authored by Julien Michel's avatar Julien Michel
Browse files

Ajout du filtre de conversion VectorImage -> ImageList.

parent 13bf108d
/*=========================================================================
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.
=========================================================================*/
#ifndef _otbVectorImageToImageListFilter_h
#define _otbVectorImageToImageListFilter_h
#include "otbImageToImageListFilter.h"
namespace otb
{
/** \class VectorImageToImageListFilter
* \brief This class aims at converting a multi-band image to a list of scalar images.
*
* This class takes a multi-band image represented as an otb::VectorImage and produces a list
* of scalar images corresponding to each band of the input image.
*
* Casting is done through standard cast operation.
*
* \ingroup Streamed
*/
template <class TVectorImageType, class TImageList>
class ITK_EXPORT VectorImageToImageListFilter
: public ImageToImageListFilter<TVectorImageType,typename TImageList::ImageType>
{
public:
/** Standard typedefs */
typedef VectorImageToImageListFilter Self;
typedef ImageToImageListFilter<TVectorImageType,
typename TImageList::ImageType> Superclass;
typedef itk::SmartPointer<Self> Pointer;
typedef itk::SmartPointer<const Self> ConstPointer;
/** Type macro */
itkNewMacro(Self);
/** Creation through object factory macro */
itkTypeMacro(VectorImageToImageListFilter,ImageToImageListFilter);
/** Template parameters typedefs */
typedef TVectorImageType InputVectorImageType;
typedef typename InputVectorImageType::ConstPointer InputVectorImageConstPointerType;
typedef TImageList OutputImageListType;
typedef typename OutputImageListType::Pointer OutputImageListPointerType;
typedef typename OutputImageListType::ImageType OutputImageType;
typedef typename OutputImageType::Pointer OutputImagePointerType;
protected:
/** Main computation method */
virtual void GenerateData(void);
/** Constructor */
VectorImageToImageListFilter() {};
/** Destructor */
virtual ~VectorImageToImageListFilter() {};
/**PrintSelf method */
virtual void PrintSelf(std::ostream& os, itk::Indent indent) const;
private:
VectorImageToImageListFilter(const Self&); //purposely not implemented
void operator=(const Self&); //purposely not implemented
};
}// End namespace otb
#ifndef OTB_MANUAL_INSTANTIATION
#include "otbVectorImageToImageListFilter.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.
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 _otbVectorImageToImageListFilter_txx
#define _otbVectorImageToImageListFilter_txx
#include "otbVectorImageToImageListFilter.h"
#include "itkImageRegionIterator.h"
#include "itkImageRegionConstIterator.h"
#include <vector>
#include "otbMacro.h"
namespace otb
{
/**
* Main computation method
*/
template <class TVectorImageType, class TImageList>
void
VectorImageToImageListFilter<TVectorImageType,TImageList>
::GenerateData(void)
{
OutputImageListPointerType outputPtr = this->GetOutput();
InputVectorImageConstPointerType inputPtr = this->GetInput();
typedef itk::ImageRegionConstIterator<InputVectorImageType> InputIteratorType;
typedef itk::ImageRegionIterator<OutputImageType> OutputIteratorType;
InputIteratorType inputIt(inputPtr,inputPtr->GetRequestedRegion());
std::vector<OutputIteratorType> outputIteratorList;
for(unsigned int i = 0;i<inputPtr->GetNumberOfComponentsPerPixel();++i)
{
typename OutputImageType::Pointer tmpImagePtr = OutputImageType::New();
tmpImagePtr->SetRegions(inputPtr->GetRequestedRegion());
tmpImagePtr->Allocate();
outputPtr->PushBack(tmpImagePtr);
outputIteratorList.push_back(OutputIteratorType(outputPtr->Back(),inputPtr->GetRequestedRegion()));
outputIteratorList.back().GoToBegin();
}
inputIt.GoToBegin();
while(!inputIt.IsAtEnd())
{
unsigned int counter = 0;
for(typename std::vector<OutputIteratorType>::iterator it = outputIteratorList.begin();
it!=outputIteratorList.end();++it)
{
(*it).Set(static_cast<typename OutputImageType::PixelType>(inputIt.Get()[counter]));
++(*it);
++counter;
}
++inputIt;
}
}
/**
* PrintSelf Method
*/
template <class TVectorImageType, class TImageList>
void
VectorImageToImageListFilter<TVectorImageType,TImageList>
::PrintSelf(std::ostream& os, itk::Indent indent) const
{
Superclass::PrintSelf(os, indent);
}
} // End namespace otb
#endif
......@@ -231,6 +231,69 @@ ADD_TEST(bfTVectorRescaleIntensityImageFilter ${BASICFILTERS_TESTS}
0 255
)
# ------- otb::VectorImageToImageListFilter ----------------------------
ADD_TEST(bfTuVectorImageToImageListFilterNew ${BASICFILTERS_TESTS}
otbVectorImageToImageListFilterNew)
ADD_TEST(bfTvVectorImageToImageListFilter ${BASICFILTERS_TESTS}
--compare-n-images ${EPSILON} 3
${BASELINE}/bfTvVectorImageToImageListFilterBand1.png
${TEMP}/bfTvVectorImageToImageListFilterBand1.png
${BASELINE}/bfTvVectorImageToImageListFilterBand2.png
${TEMP}/bfTvVectorImageToImageListFilterBand2.png
${BASELINE}/bfTvVectorImageToImageListFilterBand3.png
${TEMP}/bfTvVectorImageToImageListFilterBand3.png
otbVectorImageToImageListFilter
${INPUTDATA}/poupees.png
${TEMP}/bfTvVectorImageToImageListFilterBand1.png
${TEMP}/bfTvVectorImageToImageListFilterBand2.png
${TEMP}/bfTvVectorImageToImageListFilterBand3.png
)
# ------- otb::ImageListToVectorImageFilter ----------------------------
#ADD_TEST(bfTuImageListToVectorImageFilterNew ${BASICFILTERS_TESTS}
# otbImageListToVectorImageFilterNew)
#ADD_TEST(bfTvImageListToVectorImageFilter ${BASICFILTERS_TESTS}
# --compare-image ${EPSILON}
# ${BASELINE}/bfTvImageListToVectorImageFilter.png
# ${TEMP}/bfTvImageListToVectorImageFilter.png
# otbImageListToVectorImageFilter
# ${INPUTDATA}/poupees.c1
# ${INPUTDATA}/poupees.c2
# ${INPUTDATA}/poupees.c3
# ${TEMP}/bfTvImageListToVectorImageFilter.png
#)
# ------- otb::ImageListToImageListApplyFilter ----------------------------
#ADD_TEST(bfTuImageListToImageListApplyFilterNew ${BASICFILTERS_TESTS}
# otbImageListToImageListApplyFilterNew)
#ADD_TEST(bfTvImageListToImageListApplyFilter ${BASICFILTERS_TESTS}
# --compare-n-images ${EPSILON} 3
# ${BASELINE}/bfTvImageListToImageListApplyFilterBand1.png
# ${TEMP}/bfTvImageListToImageListApplyFilterBand1.png
# ${BASELINE}/bfTvImageListToImageListApplyFilterBand2.png
# ${TEMP}/bfTvImageListToImageListApplyFilterBand2.png
# ${BASELINE}/bfTvImageListToImageListApplyFilterBand3.png
# ${TEMP}/bfTvImageListToImageListApplyFilterBand3.png
# otbImageListToImageListApplyFilter
# ${INPUTDATA}/poupees.c1
# ${INPUTDATA}/poupees.c2
# ${INPUTDATA}/poupees.c3
# ${TEMP}/bfTvImageListToImageListApplyFilterBand1.png
# ${TEMP}/bfTvImageListToImageListApplyFilterBand2.png
# ${TEMP}/bfTvImageListToImageListApplyFilterBand3.png
#)
# A enrichir
SET(BasicFilters_SRCS
otbLeeFilter.cxx
......@@ -260,6 +323,12 @@ otbSpectralAngleDistanceImageFilterNew.cxx
otbSpectralAngleDistanceImageFilter.cxx
otbVectorRescaleIntensityImageFilterNew.cxx
otbVectorRescaleIntensityImageFilter.cxx
otbVectorImageToImageListFilterNew.cxx
otbVectorImageToImageListFilter.cxx
#otbImageListToVectorImageFilterNew.cxx
#otbImageListToVectorImageFilter.cxx
#otbImageListToImageListApplyFilterNew.cxx
#otbImageListToImageListApplyFilter.cxx
)
......
......@@ -54,4 +54,10 @@ REGISTER_TEST(otbSpectralAngleDistanceImageFilterNew);
REGISTER_TEST(otbSpectralAngleDistanceImageFilter);
REGISTER_TEST(otbVectorRescaleIntensityImageFilterNew);
REGISTER_TEST(otbVectorRescaleIntensityImageFilter);
REGISTER_TEST(otbVectorImageToImageListFilterNew);
REGISTER_TEST(otbVectorImageToImageListFilter);
//REGISTER_TEST(otbImageListToVectorImageFilterNew);
//REGISTER_TEST(otbImageListToVectorImageFilter);
//REGISTER_TEST(otbImageListToVectorImageApplyFilterNew);
//REGISTER_TEST(otbImageListToVectorImageApplyFilter);
}
/*=========================================================================
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 "itkExceptionObject.h"
#include "otbVectorImageToImageListFilter.h"
#include "otbVectorImage.h"
#include "otbImageList.h"
#include "otbImage.h"
#include "otbImageFileReader.h"
#include "otbImageFileWriter.h"
int otbVectorImageToImageListFilter(int argc, char * argv[])
{
try
{
const unsigned int Dimension = 2;
char * infname = argv[1];
char * outfname1 = argv[2];
char * outfname2 = argv[3];
char * outfname3 = argv[4];
typedef unsigned char PixelType;
typedef otb::Image<PixelType, Dimension> ImageType;
typedef otb::VectorImage<PixelType,Dimension> VectorImageType;
typedef otb::ImageList<ImageType> ImageListType;
// IO
typedef otb::ImageFileReader<VectorImageType> ReaderType;
typedef otb::ImageFileWriter<ImageType> WriterType;
typedef otb::VectorImageToImageListFilter<VectorImageType,ImageListType> VectorImageToImageListFilterType;
// Instantiating object
VectorImageToImageListFilterType::Pointer filter = VectorImageToImageListFilterType::New();
ReaderType::Pointer reader = ReaderType::New();
reader->SetFileName(infname);
reader->Update();
filter->SetInput(reader->GetOutput());
filter->Update();
WriterType::Pointer writer = WriterType::New();
writer->SetFileName(outfname1);
writer->SetInput(filter->GetOutput()->GetNthElement(0));
writer->Update();
writer = WriterType::New();
writer->SetFileName(outfname2);
writer->SetInput(filter->GetOutput()->GetNthElement(1));
writer->Update();
writer = WriterType::New();
writer->SetFileName(outfname3);
writer->SetInput(filter->GetOutput()->GetNthElement(2));
writer->Update();
}
catch( itk::ExceptionObject & err )
{
std::cout << "Exception itk::ExceptionObject thrown !" << std::endl;
std::cout << err << std::endl;
return EXIT_FAILURE;
}
catch( ... )
{
std::cout << "Unknown exception thrown !" << std::endl;
return EXIT_FAILURE;
}
return EXIT_SUCCESS;
}
/*=========================================================================
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 "itkExceptionObject.h"
#include "otbVectorImageToImageListFilter.h"
#include "otbVectorImage.h"
#include "otbImageList.h"
#include "otbImage.h"
int otbVectorImageToImageListFilterNew(int argc, char * argv[])
{
try
{
const unsigned int Dimension = 2;
typedef unsigned char PixelType;
typedef otb::Image<PixelType, Dimension> ImageType;
typedef otb::VectorImage<PixelType,Dimension> VectorImageType;
typedef otb::ImageList<ImageType> ImageListType;
typedef otb::VectorImageToImageListFilter<VectorImageType,ImageListType> VectorImageToImageListFilterType;
// Instantiating object
VectorImageToImageListFilterType::Pointer object = VectorImageToImageListFilterType::New();
}
catch( itk::ExceptionObject & err )
{
std::cout << "Exception itk::ExceptionObject thrown !" << std::endl;
std::cout << err << std::endl;
return EXIT_FAILURE;
}
catch( ... )
{
std::cout << "Unknown exception thrown !" << std::endl;
return EXIT_FAILURE;
}
return EXIT_SUCCESS;
}
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