Skip to content
Snippets Groups Projects
Commit 011961d7 authored by Julien Malik's avatar Julien Malik
Browse files

ENH: remove MatricMultiplyImageFilter, already in OTB

parent 2b383f0f
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Code/BasicFilters/otbMatrixMultiplyImageFilter.h deleted 100644 → 0
+ 0
− 168
View file @ 2b383f0f
/*=========================================================================
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 __otbMatrixMultiplyImageFilter_h
#define __otbMatrixMultiplyImageFilter_h
#include "itkMacro.h"
#include "otbUnaryFunctorImageFilter.h"
namespace otb
{
namespace Functor {
/** \class MatrixMultiplyFunctor
*
* \brief Apply a matrix multiplication to a VariableLengthVector
*
*/
template<class TInput, class TOutput, class TPrecision>
class MatrixMultiplyFunctor
{
public:
typedef TInput InputType;
typedef TOutput OutputType;
typedef TPrecision PrecisionType;
typedef vnl_vector<PrecisionType> VectorType;
typedef vnl_matrix<PrecisionType> MatrixType;
MatrixMultiplyFunctor() {}
~MatrixMultiplyFunctor() {}
unsigned int GetOutputSize()
{
return m_Matrix.rows();
}
bool operator !=(const MatrixMultiplyFunctor& other) const
{
return true;
}
bool operator ==(const MatrixMultiplyFunctor& other) const
{
return !(*this != other);
}
const MatrixType& GetMatrix()
{
return m_Matrix;
}
void SetMatrix(const MatrixType& m)
{
m_Matrix = m;
}
OutputType operator ()(const InputType& in)
{
// TODO : support different types & Precision correctly
// TODO : do the multiplication by hand ?
VectorType inVector(in.GetDataPointer(), in.Size());
VectorType outVector = m_Matrix * inVector;
// make a wrapper around outVector memory (data are not copied here)
OutputType outWrapper(outVector.data_block(), outVector.size());
// do the real copy here
return OutputType(outWrapper);
}
private:
MatrixType m_Matrix;
};
}
/** \class MatrixMultiplyImageFilter
*
* \brief Applies pixel-wise matrix multiplication to a VectorImage
*
* \ingroup Streamed
* \ingroup Threaded
*/
template <class TInputImage, class TOutputImage, class TPrecision>
class ITK_EXPORT MatrixMultiplyImageFilter :
public otb::UnaryFunctorImageFilter<TInputImage, TOutputImage,
Functor::MatrixMultiplyFunctor<typename TInputImage::PixelType,
typename TOutputImage::PixelType, TPrecision> >
{
public:
/** Standard class typedefs. */
typedef MatrixMultiplyImageFilter Self;
typedef otb::UnaryFunctorImageFilter
<TInputImage,
TOutputImage,
Functor::MatrixMultiplyFunctor<
typename TInputImage::PixelType,
typename TOutputImage::PixelType,
TPrecision>
> Superclass;
typedef itk::SmartPointer<Self> Pointer;
typedef itk::SmartPointer<const Self> ConstPointer;
typedef Functor::MatrixMultiplyFunctor<
typename TInputImage::PixelType,
typename TOutputImage::PixelType,
TPrecision> FunctorType;
typedef typename FunctorType::MatrixType MatrixType;
/** Method for creation through the object factory. */
itkNewMacro(Self);
/** Run-time type information (and related methods). */
itkTypeMacro(MatrixMultiplyImageFilter, otb::UnaryFunctorImageFilter);
/** Pixel types. */
typedef typename TInputImage::PixelType InputPixelType;
typedef typename TOutputImage::PixelType OutputPixelType;
const MatrixType& GetMatrix()
{
return this->GetFunctor().GetMatrix();
}
void SetMatrix(const MatrixType& m)
{
FunctorType f;
f.SetMatrix(m);
this->SetFunctor(f);
}
protected:
MatrixMultiplyImageFilter();
virtual ~MatrixMultiplyImageFilter() {}
void PrintSelf(std::ostream& os, itk::Indent indent) const;
private:
MatrixMultiplyImageFilter(const Self &); //purposely not implemented
void operator =(const Self&); //purposely not implemented
};
} // end namespace otb
#ifndef OTB_MANUAL_INSTANTIATION
#include "otbMatrixMultiplyImageFilter.txx"
#endif
#endif
Code/BasicFilters/otbMatrixMultiplyImageFilter.txx deleted 100644 → 0
+ 0
− 45
View file @ 2b383f0f
/*=========================================================================
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 __otbMatrixMultiplyImageFilter_txx
#define __otbMatrixMultiplyImageFilter_txx
#include "otbMatrixMultiplyImageFilter.h"
namespace otb
{
/**
*
*/
template <class TInputImage, class TOutputImage, class TPrecision>
MatrixMultiplyImageFilter<TInputImage, TOutputImage, TPrecision>
::MatrixMultiplyImageFilter()
{
}
template <class TInputImage, class TOutputImage, class TPrecision>
void
MatrixMultiplyImageFilter<TInputImage, TOutputImage, TPrecision>
::PrintSelf(std::ostream& os, itk::Indent indent) const
{
Superclass::PrintSelf(os, indent);
}
} // end namespace
#endif
Testing/otbMatrixMultiplyImageFilter.cxx deleted 100644 → 0
+ 0
− 155
View file @ 2b383f0f
/*=========================================================================
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 "otbMatrixMultiplyImageFilter.h"
#include "otbVectorImage.h"
#include "otbImageFileReader.h"
#include "otbImageFileWriter.h"
#include "otbVectorImageToMatrixImageFilter.h"
const unsigned int Dimension = 2;
typedef double PrecisionType;
typedef otb::VectorImage<PrecisionType, Dimension> ImageType;
typedef otb::ImageFileReader<ImageType> ReaderType;
typedef otb::MatrixMultiplyImageFilter<ImageType,ImageType,double> MatrixMultiplyImageFilterType;
typedef otb::VectorImageToMatrixImageFilter<ImageType> VectorImageToMatrixImageFilterType;
typedef otb::ImageFileWriter<ImageType> WriterType;
typedef VectorImageToMatrixImageFilterType::MatrixType MatrixType;
int otbMatrixMultiplyImageFilterNewTest(int argc, char * argv[])
{
MatrixMultiplyImageFilterType::Pointer filter = MatrixMultiplyImageFilterType::New();
std::cout << filter << std::endl;
return EXIT_SUCCESS;
}
int otbMatrixMultiplyImageFilterTest(int argc, char * argv[])
{
ImageType::Pointer image = ImageType::New();
const unsigned int Size = 10;
const unsigned int NbComponentIn = 4;
const unsigned int NbComponentOut = 2;
ImageType::RegionType region;
region.SetIndex(0, 0);
region.SetIndex(1, 0);
region.SetSize(0, Size);
region.SetSize(1, Size);
image->SetRegions(region);
image->SetNumberOfComponentsPerPixel(NbComponentIn);
image->Allocate();
itk::ImageRegionIteratorWithIndex<ImageType> it(image, image->GetLargestPossibleRegion());
for (it.GoToBegin(); !it.IsAtEnd(); ++it)
{
ImageType::IndexType idx = it.GetIndex();
ImageType::PixelType value;
value.SetSize(NbComponentIn);
for (unsigned int k = 0; k < NbComponentIn; ++k)
{
value[k] = idx[0] + idx[1] * Size + k;
}
it.Set(value);
}
MatrixType matrix(NbComponentOut,NbComponentIn);
for (unsigned int i = 0; i < matrix.rows(); ++i)
{
for (unsigned int j = 0; j < matrix.cols(); ++j)
{
matrix(i, j) = (i + j) * j;
}
}
MatrixMultiplyImageFilterType::Pointer mul = MatrixMultiplyImageFilterType::New();
mul->SetInput(image);
mul->SetMatrix(matrix);
mul->Update();
ImageType::Pointer outImage = mul->GetOutput();
itk::ImageRegionIteratorWithIndex<ImageType> outIt(outImage, outImage->GetLargestPossibleRegion());
for (outIt.GoToBegin(); !outIt.IsAtEnd(); ++outIt)
{
ImageType::IndexType idx = outIt.GetIndex();
ImageType::PixelType value = outIt.Get();
vnl_vector<PrecisionType> vectorValue(value.GetDataPointer(), value.GetSize());
vnl_vector<PrecisionType> input(NbComponentIn);
for (unsigned int k = 0; k < NbComponentIn; ++k)
{
input[k] = idx[0] + idx[1] * Size + k;
}
vnl_vector<PrecisionType> expected = matrix * input;
vnl_vector<PrecisionType> diff = expected - vectorValue;
if ( diff.magnitude() > 1E-10 )
{
std::cerr << "At index " << idx << ". Expected " << expected << " , got " << vectorValue << std::endl;
return EXIT_FAILURE;
}
}
return EXIT_SUCCESS;
}
int otbMatrixMultiplyImageFilterRealImageTest(int argc, char * argv[])
{
const char * inputImage = argv[1];
const char * inputEndmembers = argv[2];
const char * outputImage = argv[3];
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());
std::cout << "computing endMember2Matrix" << std::endl;
endMember2Matrix->Update();
std::cout << "computing endMember2Matrix done" << std::endl;
MatrixType endMembers = endMember2Matrix->GetMatrix();
std::cout << "endMembers : " << endMembers.rows() << " " << endMembers.cols() << std::endl;
MatrixType pinv = vnl_matrix_inverse<double>(endMembers);
std::cout << "pinv : " << pinv.rows() << " " << pinv.cols() << std::endl;
MatrixMultiplyImageFilterType::Pointer mul = MatrixMultiplyImageFilterType::New();
mul->SetInput(readerImage->GetOutput());
mul->SetMatrix(pinv);
std::cout << "computing MatrixMultiplyImageFilterType" << std::endl;
mul->Update();
WriterType::Pointer writer = WriterType::New();
writer->SetFileName(outputImage);
writer->SetInput(mul->GetOutput());
writer->Update();
return EXIT_SUCCESS;
}
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