From 38a0b18ccaa76670da9c4ff754d9df0118097ca3 Mon Sep 17 00:00:00 2001
From: Aurelien Bricier <aurelien.bricier@c-s.fr>
Date: Mon, 11 Oct 2010 16:32:19 +0200
Subject: [PATCH] ENH: one-pass computation of the 11 flusser moments
---
.../otbFlusserImageFunction.h | 280 ++++++-------
.../otbFlusserImageFunction.txx | 381 +++++++-----------
.../FlusserMomentImageExample.cxx | 19 +-
.../FeatureExtraction/otbFlusserImage.cxx | 37 +-
4 files changed, 312 insertions(+), 405 deletions(-)
diff --git a/Code/FeatureExtraction/otbFlusserImageFunction.h b/Code/FeatureExtraction/otbFlusserImageFunction.h
index d5eac225cd..9bacb02c61 100644
--- a/Code/FeatureExtraction/otbFlusserImageFunction.h
+++ b/Code/FeatureExtraction/otbFlusserImageFunction.h
@@ -1,138 +1,142 @@
-/*=========================================================================
-
- 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 __otbFlusserImageFunction_h
-#define __otbFlusserImageFunction_h
-
-#include "otbRealMomentImageFunction.h"
-
-#include <complex>
-
-namespace otb
-{
-
-/**
- * \class FlusserImageFunction
- * \brief Calculate the Flusser's invariant parameters.
- *
- * Calculate the Flusser's invariant over an image defined as:
- *
- * - \f$ \psi_{1} = c_{11} \f$
- * - \f$ \psi_{2} = c_{21} c_{12} \f$
- * - \f$ \psi_{3} = Re (c_{20} c_{12}^{2} )\f$
- * - \f$ \psi_{4} = Im (c_{20} c_{12}^{2} )\f$
- * - \f$ \psi_{5} = Re (c_{30} c_{12}^{3} )\f$
- * - \f$ \psi_{6} = Im (c_{30} c_{12}^{3} )\f$
- * - \f$ \psi_{7} = c_{22} \f$
- * - \f$ \psi_{8} = Re (c_{31} c_{12}^{2} )\f$
- * - \f$ \psi_{9} = Im (c_{31} c_{12}^{2} )\f$
- * - \f$ \psi_{10} = Re (c_{40} c_{12}^{4} )\f$
- * - \f$ \psi_{11} = Im (c_{40} c_{12}^{4} )\f$
- *
- * With :
- *
- * \f[ c_{p,q}=\int_{-\infty}^{\infty} \int_{-\infty}^{\infty} (x+iy)^{p} \cdot (x-iy)^{q} \cdot f(x,y) \cdot
- dx \cdot dy \f]
- *
- * And:
- * - \f$(x,y)\f$ pixel localization;
- * - \f$ f(x,y)\f$ the pixel value over the \f$(x,y)\f$ coordinate.
- *
- * This class is templated over the input image type and the
- * coordinate representation type (e.g. float or double).
- *
- * \ingroup ImageFunctions
- */
-
-template <class TInput,
- class TOutput = double,
- class TPrecision = double,
- class TCoordRep = float>
-class ITK_EXPORT FlusserImageFunction :
- public RealMomentImageFunction<TInput, TOutput, TPrecision, TCoordRep>
-{
-public:
- /** Standard class typedefs. */
- typedef FlusserImageFunction Self;
- typedef RealMomentImageFunction<TInput, TOutput, TPrecision, TCoordRep> Superclass;
- typedef itk::SmartPointer<Self> Pointer;
- typedef itk::SmartPointer<const Self> ConstPointer;
-
- /** Run-time type information (and related methods). */
- itkTypeMacro(FlusserImageFunction, RealMomentImageFunction);
-
- /** Method for creation through the object factory. */
- itkNewMacro(Self);
-
- /** InputImageType typedef support. */
- typedef TInput InputType;
- typedef typename Superclass::IndexType IndexType;
- typedef typename Superclass::ContinuousIndexType ContinuousIndexType;
- typedef typename Superclass::PointType PointType;
-
- typedef typename Superclass::RealType RealType;
- typedef typename std::complex<RealType> ComplexType;
-
- /** Type for calculation precision */
- typedef typename Superclass::PrecisionType PrecisionType;
-
- /** Dimension of the underlying image. */
- itkStaticConstMacro(ImageDimension, unsigned int,
- InputType::ImageDimension);
-
- /** Evalulate the function at specified index */
- virtual RealType EvaluateAtIndex(const IndexType& index) const;
-
- /** Evaluate the function at non-integer positions */
- virtual RealType Evaluate(const PointType& point) const
- {
- IndexType index;
- this->ConvertPointToNearestIndex(point, index);
- return this->EvaluateAtIndex(index);
- }
- virtual RealType EvaluateAtContinuousIndex(
- const ContinuousIndexType& cindex) const
- {
- IndexType index;
- this->ConvertContinuousIndexToNearestIndex(cindex, index);
- return this->EvaluateAtIndex(index);
- }
-
- /** Get/Set the radius of the neighborhood over which the
- statistics are evaluated */
- itkSetClampMacro(MomentNumber, short, 1, 11);
- itkGetConstReferenceMacro(MomentNumber, short);
-
-protected:
- FlusserImageFunction();
- virtual ~FlusserImageFunction() {}
- void PrintSelf(std::ostream& os, itk::Indent indent) const;
-
-private:
- FlusserImageFunction(const Self &); //purposely not implemented
- void operator =(const Self&); //purposely not implemented
-
- short m_MomentNumber;
-};
-
-} // namespace otb
-
-#ifndef OTB_MANUAL_INSTANTIATION
-#include "otbFlusserImageFunction.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 __otbFlusserImageFunction_h
+#define __otbFlusserImageFunction_h
+
+#include "itkImageFunction.h"
+#include "itkFixedArray.h"
+
+namespace otb
+{
+
+/**
+ * \class FlusserImageFunction
+ * \brief Calculate the Flusser's invariant parameters.
+ *
+ * Calculate the Flusser's invariant over an image defined as:
+ *
+ * - \f$ \psi_{1} = c_{11} \f$
+ * - \f$ \psi_{2} = c_{21} c_{12} \f$
+ * - \f$ \psi_{3} = Re (c_{20} c_{12}^{2} )\f$
+ * - \f$ \psi_{4} = Im (c_{20} c_{12}^{2} )\f$
+ * - \f$ \psi_{5} = Re (c_{30} c_{12}^{3} )\f$
+ * - \f$ \psi_{6} = Im (c_{30} c_{12}^{3} )\f$
+ * - \f$ \psi_{7} = c_{22} \f$
+ * - \f$ \psi_{8} = Re (c_{31} c_{12}^{2} )\f$
+ * - \f$ \psi_{9} = Im (c_{31} c_{12}^{2} )\f$
+ * - \f$ \psi_{10} = Re (c_{40} c_{12}^{4} )\f$
+ * - \f$ \psi_{11} = Im (c_{40} c_{12}^{4} )\f$
+ *
+ * With :
+ *
+ * \f[ c_{p,q}=\int_{-\infty}^{\infty} \int_{-\infty}^{\infty} (x+iy)^{p} \cdot (x-iy)^{q} \cdot f(x,y) \cdot
+ dx \cdot dy \f]
+ *
+ * And:
+ * - \f$(x,y)\f$ pixel localization;
+ * - \f$ f(x,y)\f$ the pixel value over the \f$(x,y)\f$ coordinate.
+ *
+ * This class is templated over the input image type and the
+ * coordinate representation type (e.g. float or double).
+ *
+ * \ingroup ImageFunctions
+ */
+
+template <class TInputImage, class TCoordRep = float >
+class ITK_EXPORT FlusserImageFunction :
+ public itk::ImageFunction< TInputImage,
+ itk::FixedArray<
+ ITK_TYPENAME itk::NumericTraits<typename TInputImage::PixelType>::RealType,
+ 11 >,
+ TCoordRep >
+{
+public:
+ /** Standard class typedefs. */
+ typedef FlusserImageFunction Self;
+ typedef itk::ImageFunction< TInputImage,
+ itk::FixedArray<
+ ITK_TYPENAME itk::NumericTraits<
+ typename TInputImage::PixelType>::RealType,
+ 11 >,
+ TCoordRep > Superclass;
+ typedef itk::SmartPointer<Self> Pointer;
+ typedef itk::SmartPointer<const Self> ConstPointer;
+
+ /** Run-time type information (and related methods). */
+ itkTypeMacro(FlusserImageFunction, ImageFunction);
+
+ /** Method for creation through the object factory. */
+ itkNewMacro(Self);
+
+ /** InputImageType typedef support. */
+ typedef TInputImage InputImageType;
+ typedef typename Superclass::IndexType IndexType;
+ typedef typename Superclass::ContinuousIndexType ContinuousIndexType;
+ typedef typename Superclass::PointType PointType;
+
+ typedef typename Superclass::OutputType RealType;
+ typedef typename RealType::ValueType ScalarRealType;
+
+ /** Dimension of the underlying image. */
+ itkStaticConstMacro(ImageDimension, unsigned int,
+ InputImageType::ImageDimension);
+
+ /** Evalulate the function at specified index */
+ virtual RealType EvaluateAtIndex(const IndexType& index) const;
+
+ /** Evaluate the function at non-integer positions */
+ virtual RealType Evaluate(const PointType& point) const
+ {
+ IndexType index;
+ this->ConvertPointToNearestIndex(point, index);
+ return this->EvaluateAtIndex(index);
+ }
+ virtual RealType EvaluateAtContinuousIndex(
+ const ContinuousIndexType& cindex) const
+ {
+ IndexType index;
+ this->ConvertContinuousIndexToNearestIndex(cindex, index);
+ return this->EvaluateAtIndex(index);
+ }
+
+ /** Get/Set the radius of the neighborhood over which the
+ * statistics are evaluated
+ */
+ itkSetMacro( NeighborhoodRadius, unsigned int );
+ itkGetConstReferenceMacro( NeighborhoodRadius, unsigned int );
+
+protected:
+ FlusserImageFunction();
+ virtual ~FlusserImageFunction() {}
+ void PrintSelf(std::ostream& os, itk::Indent indent) const;
+
+private:
+ FlusserImageFunction(const Self &); //purposely not implemented
+ void operator =(const Self&); //purposely not implemented
+
+ unsigned int m_NeighborhoodRadius;
+};
+
+} // namespace otb
+
+#ifndef OTB_MANUAL_INSTANTIATION
+#include "otbFlusserImageFunction.txx"
+#endif
+
+#endif
+
diff --git a/Code/FeatureExtraction/otbFlusserImageFunction.txx b/Code/FeatureExtraction/otbFlusserImageFunction.txx
index e8274de0b8..a0f5bf8d2e 100644
--- a/Code/FeatureExtraction/otbFlusserImageFunction.txx
+++ b/Code/FeatureExtraction/otbFlusserImageFunction.txx
@@ -1,230 +1,151 @@
-/*=========================================================================
-
- 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 __otbFlusserImageFunction_txx
-#define __otbFlusserImageFunction_txx
-
-#include "otbFlusserImageFunction.h"
-#include "otbComplexMomentImageFunction.h"
-#include "itkNumericTraits.h"
-#include "itkMacro.h"
-#include <complex>
-
-namespace otb
-{
-
-/**
- * Constructor
- */
-template <class TInput, class TOutput, class TPrecision, class TCoordRep>
-FlusserImageFunction<TInput, TOutput, TPrecision, TCoordRep>
-::FlusserImageFunction()
-{
- m_MomentNumber = -1;
-}
-
-/**
- *
- */
-template <class TInput, class TOutput, class TPrecision, class TCoordRep>
-void
-FlusserImageFunction<TInput, TOutput, TPrecision, TCoordRep>
-::PrintSelf(std::ostream& os, itk::Indent indent) const
-{
- this->Superclass::PrintSelf(os, indent);
- os << indent << " m_MomentNumber : " << m_MomentNumber << std::endl;
-}
-
-template <class TInput, class TOutput, class TPrecision, class TCoordRep>
-typename FlusserImageFunction<TInput, TOutput, TPrecision, TCoordRep>::RealType
-FlusserImageFunction<TInput, TOutput, TPrecision, TCoordRep>
-::EvaluateAtIndex(const IndexType& index) const
-{
-
- RealType FlusserValue(itk::NumericTraits<RealType>::Zero);
- ComplexType FlusserValueComplex(itk::NumericTraits<ComplexType>::Zero);
-
- typedef otb::ComplexMomentImageFunction<InputType, ComplexType> CMType;
- typename CMType::Pointer function = CMType::New();
-
- if (!this->GetInputImage())
- {
- return (itk::NumericTraits<RealType>::max());
- }
-
- if (!this->IsInsideBuffer(index))
- {
- return (itk::NumericTraits<RealType>::max());
- }
-
- assert(m_MomentNumber > 0);
- assert(m_MomentNumber < 12);
-
- function->SetInputImage(this->GetInputImage());
- function->SetNeighborhoodRadius(this->GetNeighborhoodRadius());
-
- switch (m_MomentNumber)
- {
- case 1:
- {
- ComplexType C11;
- function->SetP(1);
- function->SetQ(1);
- C11 = function->EvaluateAtIndex(index);
- FlusserValue = C11.real();
- }
- break;
- case 2:
- {
- ComplexType C21, C12;
- function->SetP(2);
- function->SetQ(1);
- C21 = function->EvaluateAtIndex(index);
- function->SetP(1);
- function->SetQ(2);
- C12 = function->EvaluateAtIndex(index);
-
- FlusserValue = vcl_abs(C21 * C12);
- }
- break;
- case 3:
- {
- ComplexType C20, C12;
- function->SetP(2);
- function->SetQ(0);
- C20 = function->EvaluateAtIndex(index);
- function->SetP(1);
- function->SetQ(2);
- C12 = function->EvaluateAtIndex(index);
- FlusserValueComplex = C20 * vcl_pow(C12, 2);
- FlusserValue = FlusserValueComplex.real();
- }
- break;
- case 4:
- {
- ComplexType C20, C12;
- function->SetP(2);
- function->SetQ(0);
- C20 = function->EvaluateAtIndex(index);
- function->SetP(1);
- function->SetQ(2);
- C12 = function->EvaluateAtIndex(index);
- FlusserValueComplex = C20 * vcl_pow(C12, 2);
- FlusserValue = FlusserValueComplex.imag();
- }
- break;
- case 5:
- {
- ComplexType C30, C12;
- function->SetP(3);
- function->SetQ(0);
- C30 = function->EvaluateAtIndex(index);
- function->SetP(1);
- function->SetQ(2);
- C12 = function->EvaluateAtIndex(index);
-
- FlusserValueComplex = C30 * vcl_pow(C12, 3);
- FlusserValue = FlusserValueComplex.real();
- }
- break;
- case 6:
- {
- ComplexType C30, C12;
- function->SetP(3);
- function->SetQ(0);
- C30 = function->EvaluateAtIndex(index);
- function->SetP(1);
- function->SetQ(2);
- C12 = function->EvaluateAtIndex(index);
-
- FlusserValueComplex = C30 * vcl_pow(C12, 3);
- FlusserValue = FlusserValueComplex.imag();
- }
- break;
- case 7:
- {
- ComplexType C22;
- function->SetP(2);
- function->SetQ(2);
- C22 = function->EvaluateAtIndex(index);
- FlusserValue = C22.real();
- }
- break;
- case 8:
- {
- ComplexType C31, C12;
- function->SetP(3);
- function->SetQ(1);
- C31 = function->EvaluateAtIndex(index);
- function->SetP(1);
- function->SetQ(2);
- C12 = function->EvaluateAtIndex(index);
- FlusserValueComplex = C31 * vcl_pow(C12, 2);
- FlusserValue = FlusserValueComplex.real();
- }
- break;
- case 9:
- {
- ComplexType C31, C12;
- function->SetP(3);
- function->SetQ(1);
- C31 = function->EvaluateAtIndex(index);
- function->SetP(1);
- function->SetQ(2);
- C12 = function->EvaluateAtIndex(index);
- FlusserValueComplex = C31 * vcl_pow(C12, 2);
- FlusserValue = FlusserValueComplex.imag();
- }
- break;
- case 10:
- {
- ComplexType C40, C12;
- function->SetP(4);
- function->SetQ(0);
- C40 = function->EvaluateAtIndex(index);
- function->SetP(1);
- function->SetQ(2);
- C12 = function->EvaluateAtIndex(index);
- FlusserValueComplex = C40 * vcl_pow(C12, 4);
- FlusserValue = FlusserValueComplex.real();
- }
- break;
- case 11:
- {
- ComplexType C40, C12;
- function->SetP(4);
- function->SetQ(0);
- C40 = function->EvaluateAtIndex(index);
- function->SetP(1);
- function->SetQ(2);
- C12 = function->EvaluateAtIndex(index);
- FlusserValueComplex = C40 * vcl_pow(C12, 4);
- FlusserValue = FlusserValueComplex.imag();
- }
- break;
-
- default:
- itkWarningMacro("Hu's invariant parameters are between 1 and 7");
- }
-
- return (static_cast<RealType>(FlusserValue));
-
-}
-
-} // 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.
+
+
+ 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 __otbFlusserImageFunction_txx
+#define __otbFlusserImageFunction_txx
+
+#include "otbFlusserImageFunction.h"
+#include "itkConstNeighborhoodIterator.h"
+#include "itkNumericTraits.h"
+#include "itkMacro.h"
+#include <complex>
+
+namespace otb
+{
+
+/**
+ * Constructor
+ */
+template <class TInputImage, class TCoordRep>
+FlusserImageFunction<TInputImage, TCoordRep>
+::FlusserImageFunction()
+{
+ m_NeighborhoodRadius = 1;
+}
+
+template <class TInputImage, class TCoordRep>
+void
+FlusserImageFunction<TInputImage, TCoordRep>
+::PrintSelf(std::ostream& os, itk::Indent indent) const
+{
+ this->Superclass::PrintSelf(os, indent);
+}
+
+template <class TInputImage, class TCoordRep>
+typename FlusserImageFunction<TInputImage,TCoordRep>::RealType
+FlusserImageFunction<TInputImage,TCoordRep>
+::EvaluateAtIndex(const IndexType& index) const
+{
+ // Build moments vector
+ RealType moments;
+
+ // Initialize moments
+ moments.Fill( itk::NumericTraits< ScalarRealType >::Zero );
+
+ // Check for input image
+ if( !this->GetInputImage() )
+ {
+ return moments;
+ }
+
+ // Check for out of buffer
+ if ( !this->IsInsideBuffer( index ) )
+ {
+ return moments;
+ }
+
+ // Define complex type
+ typedef std::complex<ScalarRealType> ComplexType;
+
+ // Define and intialize cumulants for complex moments
+ ComplexType c11, c12, c21,c20, c30, c22, c31, c40;
+ c11 = itk::NumericTraits<ComplexType>::Zero;
+ c12 = itk::NumericTraits<ComplexType>::Zero;
+ c21 = itk::NumericTraits<ComplexType>::Zero;
+ c20 = itk::NumericTraits<ComplexType>::Zero;
+ c30 = itk::NumericTraits<ComplexType>::Zero;
+ c22 = itk::NumericTraits<ComplexType>::Zero;
+ c31 = itk::NumericTraits<ComplexType>::Zero;
+ c40 = itk::NumericTraits<ComplexType>::Zero;
+
+ ScalarRealType c00 = itk::NumericTraits<ScalarRealType>::Zero;
+
+ // Create an N-d neighborhood kernel, using a zeroflux boundary condition
+ typename InputImageType::SizeType kernelSize;
+ kernelSize.Fill( m_NeighborhoodRadius );
+
+ itk::ConstNeighborhoodIterator<InputImageType>
+ it(kernelSize, this->GetInputImage(), this->GetInputImage()->GetBufferedRegion());
+
+ // Set the iterator at the desired location
+ it.SetLocation(index);
+
+ // Walk the neighborhood
+ const unsigned int size = it.Size();
+ for (unsigned int i = 0; i < size; ++i)
+ {
+ // Retrieve value, and centered-reduced position
+ ScalarRealType value = static_cast<ScalarRealType>(it.GetPixel(i));
+ ScalarRealType x = static_cast<ScalarRealType>(it.GetOffset(i)[0]);
+ ScalarRealType y = static_cast<ScalarRealType>(it.GetOffset(i)[1]);
+
+ // Build complex value
+ ComplexType xpy(x,y),xmy(x,-y);
+
+ // Update cumulants
+ c00+=value;
+ c11+=xpy*xmy*value;
+ c12+=xpy*xmy*xmy*value;
+ c21+=xpy*xpy*xmy*value;
+ c20+=xpy*xpy*value;
+ c30+=xpy*xpy*xpy*value;
+ c22+=xpy*xpy*xmy*xmy*value;
+ c31+=xpy*xpy*xpy*xmy*value;
+ c40+=xpy*xpy*xpy*xpy*value;
+ }
+
+ // Nomalisation
+ c11/=vcl_pow(c00, (1+1)/2);
+ c12/=vcl_pow(c00, (1+2)/2);
+ c21/=vcl_pow(c00, (2+1)/2);
+ c20/=vcl_pow(c00, (2+0)/2);
+ c30/=vcl_pow(c00, (3+0)/2);
+ c22/=vcl_pow(c00, (2+2)/2);
+ c31/=vcl_pow(c00, (3+1)/2);
+ c40/=vcl_pow(c00, (4+0)/2);
+
+ // Compute moments combinations
+ moments[0] = static_cast<ScalarRealType>(c11.real());
+ moments[1] = static_cast<ScalarRealType>((c21*c12).real());
+ moments[2] = static_cast<ScalarRealType>((c20*c12*c12).real());
+ moments[3] = static_cast<ScalarRealType>((c20*c12*c12).imag());
+ moments[4] = static_cast<ScalarRealType>((c30*c12*c12*c12).real());
+ moments[5] = static_cast<ScalarRealType>((c30*c12*c12*c12).imag());
+ moments[6] = static_cast<ScalarRealType>(c22.real());
+ moments[7] = static_cast<ScalarRealType>((c31*c12*c12).real());
+ moments[8] = static_cast<ScalarRealType>((c31*c12*c12).imag());
+ moments[9] = static_cast<ScalarRealType>((c40*c12*c12*c12*c12).real());
+ moments[10] = static_cast<ScalarRealType>((c40*c12*c12*c12*c12).imag());
+
+ // Return result
+ return moments;
+}
+
+} // namespace otb
+
+#endif
+
diff --git a/Examples/FeatureExtraction/FlusserMomentImageExample.cxx b/Examples/FeatureExtraction/FlusserMomentImageExample.cxx
index 0dd0f47b45..97a9b9be0f 100644
--- a/Examples/FeatureExtraction/FlusserMomentImageExample.cxx
+++ b/Examples/FeatureExtraction/FlusserMomentImageExample.cxx
@@ -47,13 +47,12 @@ int main(int argc, char * argv[])
if (argc != 3)
{
std::cerr << "Usage: " << argv[0] << " inputImageFile ";
- std::cerr << " moment_number" << std::endl;
+ std::cerr << " neighborhood_radius" << std::endl;
return EXIT_FAILURE;
}
const char * inputFilename = argv[1];
-
- unsigned int mMomentNumber((unsigned char) ::atoi(argv[2]));
+ const unsigned int radius = atoi(argv[2]);
typedef unsigned char InputPixelType;
const unsigned int Dimension = 2;
@@ -75,9 +74,8 @@ int main(int argc, char * argv[])
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
- typedef float MomentType;
- typedef otb::FlusserImageFunction<InputImageType,
- MomentType> FlusserType;
+ typedef otb::FlusserImageFunction<InputImageType> FlusserType;
+ typedef FlusserType::RealType MomentType;
FlusserType::Pointer fmFunction = FlusserType::New();
// Software Guide : EndCodeSnippet
@@ -122,7 +120,7 @@ int main(int argc, char * argv[])
// Software Guide : BeginCodeSnippet
fmFunction->SetInputImage(image);
- fmFunction->SetMomentNumber(mMomentNumber);
+ fmFunction->SetNeighborhoodRadius(radius);
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
@@ -134,8 +132,11 @@ int main(int argc, char * argv[])
// Software Guide : BeginCodeSnippet
MomentType Result = fmFunction->EvaluateAtIndex(center);
- std::cout << "The moment of order " << mMomentNumber <<
- " is equal to " << Result << std::endl;
+ for (unsigned int j=0; j<11; j++)
+ {
+ std::cout << "The moment of order " << j+1 <<
+ " is equal to " << Result[j] << std::endl;
+ }
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
diff --git a/Testing/Code/FeatureExtraction/otbFlusserImage.cxx b/Testing/Code/FeatureExtraction/otbFlusserImage.cxx
index 5a49e32aa6..7ba36af02f 100644
--- a/Testing/Code/FeatureExtraction/otbFlusserImage.cxx
+++ b/Testing/Code/FeatureExtraction/otbFlusserImage.cxx
@@ -33,7 +33,6 @@ int otbFlusserImage(int argc, char * argv[])
{
const char * inputFilename = argv[1];
const char * outputFilename = argv[2];
- unsigned int Number = 1;
typedef unsigned char InputPixelType;
const unsigned int Dimension = 2;
@@ -41,48 +40,30 @@ int otbFlusserImage(int argc, char * argv[])
typedef itk::Image<InputPixelType, Dimension> InputImageType;
typedef otb::ImageFileReader<InputImageType> ReaderType;
typedef std::complex<float> ComplexType;
- typedef float RealType;
- typedef otb::FlusserImageFunction<InputImageType, float, float> FunctionType;
-
- InputImageType::RegionType region;
- InputImageType::SizeType size;
- InputImageType::IndexType start;
-
- start.Fill(0);
- size[0] = 50;
- size[1] = 50;
+ typedef otb::FlusserImageFunction<InputImageType> FunctionType;
+ typedef FunctionType::RealType RealType;
ReaderType::Pointer reader = ReaderType::New();
- FunctionType::Pointer function = FunctionType::New();
+ FunctionType::Pointer function = FunctionType::New();
reader->SetFileName(inputFilename);
-
- InputImageType::Pointer image = reader->GetOutput();
-
- region.SetIndex(start);
- region.SetSize(size);
-
- image->SetRegions(region);
- image->Update();
- function->SetInputImage(image);
+ reader->Update();
+ function->SetInputImage(reader->GetOutput());
InputImageType::IndexType index;
index[0] = 10;
index[1] = 10;
+ function->SetNeighborhoodRadius(3);
RealType Result;
+ Result = function->EvaluateAtIndex(index);
std::ofstream outputStream(outputFilename);
outputStream << std::setprecision(10) << "Flusser Image moments: [12]" << std::endl;
- for (Number = 1; Number < 12; Number++)
+ for (unsigned int j = 1; j < 12; j++)
{
- //OTB-FA-00024-CS
- function->SetMomentNumber(Number);
- //OTB-FA-00025-CS
- function->SetNeighborhoodRadius(3);
- Result = function->EvaluateAtIndex(index);
- outputStream << "Flusser(" << Number << ") = " << Result << std::endl;
+ outputStream << "Flusser(" << j << ") = " << Result[j-1] << std::endl;
}
outputStream.close();
--
GitLab