]> Creatis software - FrontAlgorithms.git/commitdiff
...
authorLeonardo Florez-Valencia <florez-l@javeriana.edu.co>
Fri, 1 Apr 2016 19:00:45 +0000 (14:00 -0500)
committerLeonardo Florez-Valencia <florez-l@javeriana.edu.co>
Fri, 1 Apr 2016 19:00:45 +0000 (14:00 -0500)
lib/fpa/Image/Functors/FluxMedialness.h [new file with mode: 0644]
lib/fpa/Image/Functors/FluxMedialness.hxx [new file with mode: 0644]
lib/fpa/Image/Functors/GradientImageFunctionBase.h [new file with mode: 0644]
lib/fpa/Image/Functors/GradientImageFunctionBase.hxx [new file with mode: 0644]
lib/fpa/Image/Functors/GulsunTekMedialness.h [new file with mode: 0644]
lib/fpa/Image/Functors/GulsunTekMedialness.hxx [new file with mode: 0644]
lib/fpa/Image/Functors/MFluxMedialness.h [new file with mode: 0644]
lib/fpa/Image/Functors/MFluxMedialness.hxx [new file with mode: 0644]
plugins/fpa/CMakeLists.txt
plugins/fpa/GradientBaseImageFunctionSource.cxx [new file with mode: 0644]
plugins/fpa/GradientBaseImageFunctionSource.h [new file with mode: 0644]

diff --git a/lib/fpa/Image/Functors/FluxMedialness.h b/lib/fpa/Image/Functors/FluxMedialness.h
new file mode 100644 (file)
index 0000000..1995040
--- /dev/null
@@ -0,0 +1,78 @@
+#ifndef __FPA__IMAGE__FUNCTORS__FLUXMEDIALNESS__H__
+#define __FPA__IMAGE__FUNCTORS__FLUXMEDIALNESS__H__
+
+#include <fpa/Image/Functors/GradientImageFunctionBase.h>
+
+namespace fpa
+{
+  namespace Image
+  {
+    namespace Functors
+    {
+      /**
+       */
+      template< class _TGradient >
+      class FluxMedialness
+        : public GradientImageFunctionBase< _TGradient >
+      {
+      public:
+        typedef FluxMedialness                          Self;
+        typedef GradientImageFunctionBase< _TGradient > Superclass;
+        typedef itk::SmartPointer< Self >               Pointer;
+        typedef itk::SmartPointer< const Self >         ConstPointer;
+
+        itkStaticConstMacro( Dimension, unsigned int, Superclass::Dimension );
+
+        typedef typename Superclass::TOutput TOutput;
+        typedef typename Superclass::TScalar TScalar;
+        typedef typename Superclass::TIndex  TIndex;
+        typedef typename Superclass::TVector TVector;
+        typedef typename Superclass::TPoint  TPoint;
+
+        typedef std::vector< double > TRCandidates;
+
+      public:
+        itkNewMacro( Self );
+        itkTypeMacro( FluxMedialness, GradientImageFunctionBase );
+
+        itkGetConstMacro( RadiusStep, double );
+        itkGetConstMacro( MinRadius, double );
+        itkGetConstMacro( MaxRadius, double );
+        itkGetConstMacro( RadialSampling, unsigned int );
+
+        itkSetMacro( RadiusStep, double );
+        itkSetMacro( MinRadius, double );
+        itkSetMacro( MaxRadius, double );
+        itkSetMacro( RadialSampling, unsigned int );
+
+      protected:
+        FluxMedialness( );
+        virtual ~FluxMedialness( );
+
+        virtual TOutput _Evaluate( const TIndex& i ) const;
+
+      private:
+        // Purposely not implemented.
+        FluxMedialness( const Self& );
+        void operator=( const Self& );
+
+      protected:
+        double       m_MinRadius;
+        double       m_MaxRadius;
+        unsigned int m_RadialSampling;
+        double       m_RadiusStep;
+      };
+
+    } // ecapseman
+
+  } // ecapseman
+
+} // ecapseman
+
+#ifndef ITK_MANUAL_INSTANTIATION
+#include <fpa/Image/Functors/FluxMedialness.hxx>
+#endif // ITK_MANUAL_INSTANTIATION
+
+#endif // __FPA__IMAGE__FUNCTORS__FLUXMEDIALNESS__H__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Image/Functors/FluxMedialness.hxx b/lib/fpa/Image/Functors/FluxMedialness.hxx
new file mode 100644 (file)
index 0000000..8f8651b
--- /dev/null
@@ -0,0 +1,108 @@
+#ifndef __FPA__IMAGE__FUNCTORS__FLUXMEDIALNESS__HXX__
+#define __FPA__IMAGE__FUNCTORS__FLUXMEDIALNESS__HXX__
+
+#include <cmath>
+#include <vnl/vnl_math.h>
+#include <itkLineConstIterator.h>
+
+// -------------------------------------------------------------------------
+template< class _TGradient >
+fpa::Image::Functors::FluxMedialness< _TGradient >::
+FluxMedialness( )
+  : Superclass( ),
+    m_MinRadius( double( 0 ) ),
+    m_MaxRadius( double( 1 ) ),
+    m_RadialSampling( 4 ),
+    m_RadiusStep( double( 1 ) )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TGradient >
+fpa::Image::Functors::FluxMedialness< _TGradient >::
+~FluxMedialness( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TGradient >
+typename fpa::Image::Functors::FluxMedialness< _TGradient >::
+TOutput fpa::Image::Functors::FluxMedialness< _TGradient >::
+_Evaluate( const TIndex& i ) const
+{
+  itk::Object::GlobalWarningDisplayOff( );
+
+  double pi2n = double( 2 ) * double( vnl_math::pi );
+  pi2n /= int( this->m_RadialSampling );
+  const _TGradient* img = this->GetInputImage( );
+  // Gradient in central pixel
+  //TVector grad_idx = img->GetPixel( i );
+
+  double Flux1 = 0;
+  double Flux = 0;
+
+  TRCandidates FluxFinal;
+  TRCandidates radiusGenerated;
+  double dR = double( 0 );
+  double optR = double( 0 );
+  TPoint center;
+  img->TransformIndexToPhysicalPoint( i, center );
+  double radius = double(0);
+
+  for( unsigned int cx = 0; cx < Self::Dimension - 1; cx++ )
+  {
+    for( unsigned int cy = cx + 1; cy < Self::Dimension; cy++ )
+    {
+      dR = double( 0 );
+      FluxFinal.clear();
+      radiusGenerated.clear();
+      radius = this->m_MinRadius;
+      while( radius <= this->m_MaxRadius )
+      {
+        Flux = 0;
+        for( unsigned int I_radial = 0; I_radial < this->m_RadialSampling ; I_radial++ )
+        {
+          Flux1 = 0;
+
+          // Direction of first profile
+          typename TPoint::VectorType dir1;
+          dir1.Fill( double( 0 ) );
+          dir1[ cx ] = std::cos( pi2n * double( I_radial ) );
+          dir1[ cy ] = std::sin( pi2n * double( I_radial ) );
+          dir1 *= (radius);
+          TIndex rIdx;
+          if (img->TransformPhysicalPointToIndex( center + dir1, rIdx ))
+          {
+            TVector grad_rIdx = img->GetPixel( rIdx );
+            TVector u_i1;
+            u_i1.SetVnlVector( ( center - ( center + dir1 ) ).GetVnlVector( ) );
+            u_i1.Normalize( );
+            // dot product
+            Flux1 = grad_rIdx * u_i1;
+          }
+
+          Flux += Flux1;
+
+        } // rof
+        //std::cout<<" radius:"<<radius<<std::endl;
+        //std::cout<<"Center:"<<center[0]<<" "<<center[1]<<std::endl;
+        //std::cout<<"edge:"<<center[0]+radius*std::cos( pi2n * double( 0 ) )<<std::endl;
+        Flux /= this->m_RadialSampling;
+        FluxFinal.push_back(Flux);
+        radiusGenerated.push_back(radius);
+        radius += this->m_RadiusStep;
+
+      }     //elihw
+
+      dR= *( std::max_element( FluxFinal.begin(), FluxFinal.end() ) );
+      optR= (dR>optR)? dR:optR;
+
+    } // rof
+
+  } // rof
+  return( TScalar(optR) );
+}
+
+#endif // __FPA__IMAGE__FUNCTORS__FLUXMEDIALNESS__HXX__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Image/Functors/GradientImageFunctionBase.h b/lib/fpa/Image/Functors/GradientImageFunctionBase.h
new file mode 100644 (file)
index 0000000..0b53f11
--- /dev/null
@@ -0,0 +1,73 @@
+#ifndef __FPA__IMAGE__FUNCTORS__GRADIENTIMAGEFUNCTIONBASE__H__
+#define __FPA__IMAGE__FUNCTORS__GRADIENTIMAGEFUNCTIONBASE__H__
+
+#include <itkImageFunction.h>
+
+namespace fpa
+{
+  namespace Image
+  {
+    namespace Functors
+    {
+      /**
+       * Base class to compute values based on image gradients (vector).
+       * It allows incremental computation of the gradient.
+       */
+      template< class _TGradient >
+      class GradientImageFunctionBase
+        : public itk::ImageFunction< _TGradient, typename _TGradient::PixelType::ValueType, typename _TGradient::PixelType::ValueType >
+      {
+      public:
+        // Types from input arguments
+        typedef _TGradient                     TGradient;
+        typedef typename _TGradient::PixelType TVector;
+        typedef typename TVector::ValueType    TScalar;
+        itkStaticConstMacro( Dimension, unsigned int, _TGradient::ImageDimension );
+
+        // Standard itk types
+        typedef GradientImageFunctionBase                          Self;
+        typedef itk::ImageFunction< _TGradient, TScalar, TScalar > Superclass;
+        typedef itk::SmartPointer< Self >                          Pointer;
+        typedef itk::SmartPointer< const Self >                    ConstPointer;
+
+        // Types from base itk::ImageFunction
+        typedef typename Superclass::InputType           TInput;
+        typedef typename Superclass::OutputType          TOutput;
+        typedef typename Superclass::PointType           TPoint;
+        typedef typename Superclass::ContinuousIndexType TContIndex;
+        typedef typename Superclass::IndexType           TIndex;
+
+      public:
+        itkTypeMacro( GradientImageFunctionBase, itkImageFunction );
+
+      public:
+        virtual void Prepare( ) const;
+        virtual TOutput Evaluate( const TPoint& p ) const;
+        virtual TOutput EvaluateAtIndex( const TIndex& i ) const;
+        virtual TOutput EvaluateAtContinuousIndex( const TContIndex& i ) const;
+
+      protected:
+        GradientImageFunctionBase( );
+        virtual ~GradientImageFunctionBase( );
+
+        virtual TOutput _Evaluate( const TIndex& i ) const = 0;
+
+      private:
+        // Purposely not implemented.
+        GradientImageFunctionBase( const Self& );
+        void operator=( const Self& );
+      };
+
+    } // ecapseman
+
+  } // ecapseman
+
+} // ecapseman
+
+#ifndef ITK_MANUAL_INSTANTIATION
+#include <fpa/Image/Functors/GradientImageFunctionBase.hxx>
+#endif // ITK_MANUAL_INSTANTIATION
+
+#endif // __FPA__IMAGE__FUNCTORS__GRADIENTIMAGEFUNCTIONBASE__H__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Image/Functors/GradientImageFunctionBase.hxx b/lib/fpa/Image/Functors/GradientImageFunctionBase.hxx
new file mode 100644 (file)
index 0000000..34acf1c
--- /dev/null
@@ -0,0 +1,59 @@
+#ifndef __FPA__IMAGE__FUNCTORS__GRADIENTIMAGEFUNCTIONBASE__HXX__
+#define __FPA__IMAGE__FUNCTORS__GRADIENTIMAGEFUNCTIONBASE__HXX__
+
+// -------------------------------------------------------------------------
+template< class _TGradient >
+void fpa::Image::Functors::GradientImageFunctionBase< _TGradient >::
+Prepare( ) const
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TGradient >
+typename fpa::Image::Functors::GradientImageFunctionBase< _TGradient >::
+TOutput fpa::Image::Functors::GradientImageFunctionBase< _TGradient >::
+Evaluate( const TPoint& p ) const
+{
+  TIndex i;
+  this->GetInputImage( )->TransformPhysicalPointToIndex( p, i );
+  return( this->EvaluateAtIndex( i ) );
+}
+
+// -------------------------------------------------------------------------
+template< class _TGradient >
+typename fpa::Image::Functors::GradientImageFunctionBase< _TGradient >::
+TOutput fpa::Image::Functors::GradientImageFunctionBase< _TGradient >::
+EvaluateAtIndex( const TIndex& i ) const
+{
+  return( this->_Evaluate( i ) );
+}
+
+// -------------------------------------------------------------------------
+template< class _TGradient >
+typename fpa::Image::Functors::GradientImageFunctionBase< _TGradient >::
+TOutput fpa::Image::Functors::GradientImageFunctionBase< _TGradient >::
+EvaluateAtContinuousIndex( const TContIndex& i ) const
+{
+  TPoint p;
+  this->GetInputImage( )->TransformContinuousIndexToPhysicalPoint( i, p );
+  return( this->Evaluate( p ) );
+}
+
+// -------------------------------------------------------------------------
+template< class _TGradient >
+fpa::Image::Functors::GradientImageFunctionBase< _TGradient >::
+GradientImageFunctionBase( )
+  : Superclass( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TGradient >
+fpa::Image::Functors::GradientImageFunctionBase< _TGradient >::
+~GradientImageFunctionBase( )
+{
+}
+
+#endif // __FPA__IMAGE__FUNCTORS__GRADIENTIMAGEFUNCTIONBASE__HXX__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Image/Functors/GulsunTekMedialness.h b/lib/fpa/Image/Functors/GulsunTekMedialness.h
new file mode 100644 (file)
index 0000000..81ea1d5
--- /dev/null
@@ -0,0 +1,80 @@
+#ifndef __FPA__IMAGE__FUNCTORS__GULSUNTEKMEDIALNESS__H__
+#define __FPA__IMAGE__FUNCTORS__GULSUNTEKMEDIALNESS__H__
+
+#include <fpa/Image/Functors/GradientImageFunctionBase.h>
+
+namespace fpa
+{
+  namespace Image
+  {
+    namespace Functors
+    {
+      /**
+       */
+      template< class _TGradient >
+      class GulsunTekMedialness
+        : public GradientImageFunctionBase< _TGradient >
+      {
+      public:
+        typedef GulsunTekMedialness                     Self;
+        typedef GradientImageFunctionBase< _TGradient > Superclass;
+        typedef itk::SmartPointer< Self >               Pointer;
+        typedef itk::SmartPointer< const Self >         ConstPointer;
+
+        itkStaticConstMacro( Dimension, unsigned int, Superclass::Dimension );
+
+        typedef typename Superclass::TOutput TOutput;
+        typedef typename Superclass::TScalar TScalar;
+        typedef typename Superclass::TIndex  TIndex;
+        typedef typename Superclass::TVector TVector;
+        typedef typename Superclass::TPoint  TPoint;
+        typedef typename TIndex::OffsetType  TOffset;
+
+        typedef std::vector< double >  TProfile;
+        typedef std::vector< TOffset > TOffsets;
+
+      public:
+        itkNewMacro( Self );
+        itkTypeMacro( GulsunTekMedialness, GradientImageFunctionBase );
+
+        itkGetConstMacro( MinRadius, double );
+        itkGetConstMacro( MaxRadius, double );
+        itkGetConstMacro( ProfileSampling, unsigned int );
+        itkGetConstMacro( RadialSampling, unsigned int );
+
+        itkSetMacro( MinRadius, double );
+        itkSetMacro( MaxRadius, double );
+        itkSetMacro( ProfileSampling, unsigned int );
+        itkSetMacro( RadialSampling, unsigned int );
+
+      protected:
+        GulsunTekMedialness( );
+        virtual ~GulsunTekMedialness( );
+
+        virtual TOutput _Evaluate( const TIndex& i ) const;
+
+      private:
+        // Purposely not implemented.
+        GulsunTekMedialness( const Self& );
+        void operator=( const Self& );
+
+      protected:
+        double       m_MinRadius;
+        double       m_MaxRadius;
+        unsigned int m_ProfileSampling;
+        unsigned int m_RadialSampling;
+      };
+
+    } // ecapseman
+
+  } // ecapseman
+
+} // ecapseman
+
+#ifndef ITK_MANUAL_INSTANTIATION
+#include <fpa/Image/Functors/GulsunTekMedialness.hxx>
+#endif // ITK_MANUAL_INSTANTIATION
+
+#endif // __FPA__IMAGE__FUNCTORS__GULSUNTEKMEDIALNESS__H__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Image/Functors/GulsunTekMedialness.hxx b/lib/fpa/Image/Functors/GulsunTekMedialness.hxx
new file mode 100644 (file)
index 0000000..9cf48c9
--- /dev/null
@@ -0,0 +1,114 @@
+#ifndef __FPA__IMAGE__FUNCTORS__GULSUNTEKMEDIALNESS__HXX__
+#define __FPA__IMAGE__FUNCTORS__GULSUNTEKMEDIALNESS__HXX__
+
+#include <cmath>
+#include <vnl/vnl_math.h>
+#include <itkLineConstIterator.h>
+
+// -------------------------------------------------------------------------
+template< class _TGradient >
+fpa::Image::Functors::GulsunTekMedialness< _TGradient >::
+GulsunTekMedialness( )
+  : Superclass( ),
+    m_MinRadius( double( 0 ) ),
+    m_MaxRadius( double( 1 ) ),
+    m_ProfileSampling( 4 ),
+    m_RadialSampling( 10 )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TGradient >
+fpa::Image::Functors::GulsunTekMedialness< _TGradient >::
+~GulsunTekMedialness( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TGradient >
+typename fpa::Image::Functors::GulsunTekMedialness< _TGradient >::
+TOutput fpa::Image::Functors::GulsunTekMedialness< _TGradient >::
+_Evaluate( const TIndex& i ) const
+{
+  itk::Object::GlobalWarningDisplayOff( );
+
+  // Various values
+  const _TGradient* in = this->GetInputImage( );
+  double pi2n =
+    double( 2 ) * double( vnl_math::pi ) /
+    double( this->m_ProfileSampling );
+  double rOff = this->m_MaxRadius / double( this->m_RadialSampling - 1 );
+  double optR = double( 0 );
+  TPoint pnt;
+  in->TransformIndexToPhysicalPoint( i, pnt );
+
+  // Main loop
+  for( unsigned int cx = 0; cx < Self::Dimension - 1; cx++ )
+  {
+    for( unsigned int cy = cx + 1; cy < Self::Dimension; cy++ )
+    {
+      TProfile maxProfile( this->m_RadialSampling, double( 0 ) );
+      for( unsigned int p = 0; p < this->m_ProfileSampling; p++ )
+      {
+        double a = pi2n * double( p );
+
+        // Direction of this profile
+        TVector dir;
+        dir.Fill( TScalar( 0 ) );
+        dir[ cx ] = TScalar( std::cos( a ) );
+        dir[ cy ] = TScalar( std::sin( a ) );
+
+        double maxrise = double( 0 );
+        double maxfall = double( -1 );
+        TProfile profile;
+        for( unsigned int r = 0; r < this->m_RadialSampling; r++ )
+        {
+          double radius = double( r ) * rOff;
+          TIndex idx;
+          typename TPoint::VectorType aux;
+          aux.SetVnlVector( dir.GetVnlVector( ) );
+          if(
+            in->TransformPhysicalPointToIndex( pnt + ( aux * radius ), idx )
+            )
+          {
+            TVector g = in->GetPixel( idx );
+            double b = double( g.GetNorm( ) );
+            if( double( g * dir ) < double( 0 ) )
+              b *= double( -1 );
+            maxrise = ( b > maxrise )? b: maxrise;
+            if( radius >= this->m_MinRadius )
+              maxfall = ( b < maxfall )? b: maxfall;
+            profile.push_back( -b - maxrise );
+          }
+          else
+            profile.push_back( double( 0 ) );
+
+        } // rof
+
+        for( unsigned int r = 0; r < this->m_RadialSampling; r++ )
+        {
+          double E = profile[ r ] / -maxfall;
+          E = ( E < double( 0 ) )? double( 0 ): E;
+          E = ( E > double( 1 ) )? double( 1 ): E;
+          maxProfile[ r ] += E;
+
+        } // rof
+
+      } // rof
+
+      for( unsigned int r = 0; r < this->m_RadialSampling; r++ )
+      {
+        double E = maxProfile[ r ] / double( this->m_RadialSampling );
+        optR = ( E > optR )? E: optR;
+
+      } // rof
+
+    } // rof
+
+  } // rof
+  return( TScalar( optR ) );
+}
+
+#endif // __FPA__IMAGE__FUNCTORS__GULSUNTEKMEDIALNESS__HXX__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Image/Functors/MFluxMedialness.h b/lib/fpa/Image/Functors/MFluxMedialness.h
new file mode 100644 (file)
index 0000000..76efda3
--- /dev/null
@@ -0,0 +1,78 @@
+#ifndef __FPA__IMAGE__FUNCTORS__MFLUXMEDIALNESS__H__
+#define __FPA__IMAGE__FUNCTORS__MFLUXMEDIALNESS__H__
+
+#include <fpa/Image/Functors/GradientImageFunctionBase.h>
+
+namespace fpa
+{
+  namespace Image
+  {
+    namespace Functors
+    {
+      /**
+       */
+      template< class _TGradient >
+      class MFluxMedialness
+        : public GradientImageFunctionBase< _TGradient >
+      {
+      public:
+        typedef MFluxMedialness                         Self;
+        typedef GradientImageFunctionBase< _TGradient > Superclass;
+        typedef itk::SmartPointer< Self >               Pointer;
+        typedef itk::SmartPointer< const Self >         ConstPointer;
+
+        itkStaticConstMacro( Dimension, unsigned int, Superclass::Dimension );
+
+        typedef typename Superclass::TOutput TOutput;
+        typedef typename Superclass::TScalar TScalar;
+        typedef typename Superclass::TIndex  TIndex;
+        typedef typename Superclass::TVector TVector;
+        typedef typename Superclass::TPoint  TPoint;
+
+        typedef std::vector< double > TRCandidates;
+
+      public:
+        itkNewMacro( Self );
+        itkTypeMacro( MFluxMedialness, GradientImageFunctionBase );
+
+        itkGetConstMacro( RadiusStep, double );
+        itkGetConstMacro( MinRadius, double );
+        itkGetConstMacro( MaxRadius, double );
+        itkGetConstMacro( RadialSampling, unsigned int );
+
+        itkSetMacro( RadiusStep, double );
+        itkSetMacro( MinRadius, double );
+        itkSetMacro( MaxRadius, double );
+        itkSetMacro( RadialSampling, unsigned int );
+
+      protected:
+        MFluxMedialness( );
+        virtual ~MFluxMedialness( );
+
+        virtual TOutput _Evaluate( const TIndex& i ) const;
+
+      private:
+        // Purposely not implemented.
+        MFluxMedialness( const Self& );
+        void operator=( const Self& );
+
+      protected:
+        double       m_MinRadius;
+        double       m_MaxRadius;
+        unsigned int m_RadialSampling;
+        double       m_RadiusStep;
+      };
+
+    } // ecapseman
+
+  } // ecapseman
+
+} // ecapseman
+
+#ifndef ITK_MANUAL_INSTANTIATION
+#include <fpa/Image/Functors/MFluxMedialness.hxx>
+#endif // ITK_MANUAL_INSTANTIATION
+
+#endif // __FPA__IMAGE__FUNCTORS__MFLUXMEDIALNESS__H__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Image/Functors/MFluxMedialness.hxx b/lib/fpa/Image/Functors/MFluxMedialness.hxx
new file mode 100644 (file)
index 0000000..5e902e5
--- /dev/null
@@ -0,0 +1,159 @@
+#ifndef __FPA__IMAGE__FUNCTORS__MFLUXMEDIALNESS__HXX__
+#define __FPA__IMAGE__FUNCTORS__MFLUXMEDIALNESS__HXX__
+
+#include <cmath>
+#include <vnl/vnl_math.h>
+#include <itkLineConstIterator.h>
+
+// -------------------------------------------------------------------------
+template< class _TGradient >
+fpa::Image::Functors::MFluxMedialness< _TGradient >::
+MFluxMedialness( )
+  : Superclass( ),
+    m_MinRadius( double( 0 ) ),
+    m_MaxRadius( double( 1 ) ),
+    m_RadialSampling( 4 ),
+    m_RadiusStep( double( 1 ) )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TGradient >
+fpa::Image::Functors::MFluxMedialness< _TGradient >::
+~MFluxMedialness( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TGradient >
+typename fpa::Image::Functors::MFluxMedialness< _TGradient >::
+TOutput fpa::Image::Functors::MFluxMedialness< _TGradient >::
+_Evaluate( const TIndex& i ) const
+{
+  itk::Object::GlobalWarningDisplayOff( );
+
+  double pi2n = double( 2 ) * double( vnl_math::pi );
+  pi2n /= int( this->m_RadialSampling );
+  const _TGradient* img = this->GetInputImage( );
+  //const itk::Image::SpacingType& input_spacing = img->GetSpacing( );
+
+  double Flux1 = 0;
+  double Flux2 = 0;
+  double MFlux = 0;
+
+  TRCandidates FluxFinal;
+  TRCandidates radiusGenerated;
+  double dR = double( 0 );
+  double optR = double( 0 );
+  TPoint center;
+  img->TransformIndexToPhysicalPoint( i, center );
+  double radius;
+
+  for( unsigned int cx = 0; cx < Self::Dimension - 1; cx++ )
+  {
+    for( unsigned int cy = cx + 1; cy < Self::Dimension; cy++ )
+    {
+      dR = double( 0 );
+      FluxFinal.clear();
+      radiusGenerated.clear();
+      radius = this->m_MinRadius;
+      while( radius <= this->m_MaxRadius )
+      {
+        MFlux = 0;
+        for( unsigned int I_radial = 0; I_radial < this->m_RadialSampling / 2; I_radial++ )
+        {
+          Flux1 = 0;
+          Flux2 = 0;
+
+          // Direction of first profile
+          typename TPoint::VectorType dir1;
+          dir1.Fill( double( 0 ) );
+          dir1[ cx ] = std::cos( pi2n * double( I_radial ) );
+          dir1[ cy ] = std::sin( pi2n * double( I_radial ) );
+          //dir1 *= (radius);
+
+          TIndex rIdx;
+
+          if ( img->TransformPhysicalPointToIndex( center + (dir1*radius), rIdx ) )
+          {
+            TVector grad_rIdx = img->GetPixel( rIdx );
+            TVector u_i1;
+            u_i1.SetVnlVector( ( center - ( center + dir1 ) ).GetVnlVector( ) );
+            u_i1.Normalize( );
+            // dot product
+            Flux1 = grad_rIdx * u_i1;
+          }
+          else
+          {
+            //if (Self::Dimension==3)
+            //{
+            //std::cout<<"Point Edge x:"<<center[0]+dir1[0] ;
+            //std::cout<<" y:"<<center[1]+dir1[1]<<" z:"<<center[2]+dir1[2]<<std::endl;
+            //}
+            //else if (Self::Dimension==2)
+            //{
+            //std::cout<<"Point Edge x:"<<center[0]+dir1[0] ;
+            //std::cout<<" y:"<<center[1]+dir1[1]<<std::endl;
+            //}
+          }
+
+          // Direction of second profile
+          // pi2n*Iradial + 180°
+          typename TPoint::VectorType dir2;
+          dir2.Fill( double( 0 ) );
+          dir2[ cx ] =  std::cos( (pi2n) * double( I_radial ) + double( vnl_math::pi ));
+          dir2[ cy ] =  std::sin( (pi2n) * double( I_radial ) + double( vnl_math::pi ));
+
+          TIndex rIdx2;
+
+          if ( img->TransformPhysicalPointToIndex( center + (dir2*radius), rIdx2 ) )
+          {
+            TVector grad_rIdx2 = img->GetPixel( rIdx2 );
+            TVector u_i2;
+            u_i2.SetVnlVector( ( center - ( center + dir2 ) ).GetVnlVector( ) );
+            u_i2.Normalize( );
+
+            Flux2 = grad_rIdx2 * u_i2;
+          }
+          else
+          {
+            //if (Self::Dimension==3)
+            //{
+            //std::cout<<"Point Edge x:"<<center[0]+dir2[0] ;
+            //std::cout<<" y:"<<center[1]+dir2[1]<<" z:"<<center[2]+dir2[2]<<std::endl;
+            //}
+            //else if (Self::Dimension==2)
+            //{
+            //std::cout<<"Point Edge x:"<<center[0]+dir2[0] ;
+            //std::cout<<" y:"<<center[1]+dir2[1]<<std::endl;
+            //}
+          }
+
+          MFlux += std::min( Flux1, Flux2 );
+        } // rof
+
+        //std::cout<<Self::Dimension<<" radius:"<<radius<<std::endl;
+        //std::cout<<"Center:"<<center[0]<<" "<<center[1]<<std::endl;
+        //std::cout<<"edge:"<<center[0]+radius*std::cos( pi2n * double( 0 ) )<<std::endl;
+
+        MFlux *= 2;
+        MFlux /= this->m_RadialSampling;
+        FluxFinal.push_back(MFlux);
+        radiusGenerated.push_back(radius);
+
+        radius += this->m_RadiusStep;
+
+      }     //elihw
+
+      dR= *( std::max_element( FluxFinal.begin(), FluxFinal.end() ) );
+      optR= (dR>optR)? dR:optR;
+
+    } // rof
+
+  } // rof
+  return( TScalar(optR) );
+}
+
+#endif // __FPA__IMAGE__FUNCTORS__MFLUXMEDIALNESS__HXX__
+
+// eof - $RCSfile$
index 9196099e76290b1390756c601954976279d3a8f9..0186f7f0fc70cba8ec7ad39aeecab0b8820791d6 100644 (file)
@@ -21,6 +21,7 @@ SET(
   ${CMAKE_CURRENT_SOURCE_DIR}/ExtractEndPointsAndBifurcationsFromMinimumSpanningTree.h
   ${CMAKE_CURRENT_SOURCE_DIR}/RegionGrowThresholdFunction.h
   ${CMAKE_CURRENT_SOURCE_DIR}/InvertCostFunction.h
+  ${CMAKE_CURRENT_SOURCE_DIR}/GradientBaseImageFunctionSource.h
   )
 
 SET(
@@ -44,6 +45,7 @@ SET(
   ${CMAKE_CURRENT_SOURCE_DIR}/ExtractEndPointsAndBifurcationsFromMinimumSpanningTree.cxx
   ${CMAKE_CURRENT_SOURCE_DIR}/RegionGrowThresholdFunction.cxx
   ${CMAKE_CURRENT_SOURCE_DIR}/InvertCostFunction.cxx
+  ${CMAKE_CURRENT_SOURCE_DIR}/GradientBaseImageFunctionSource.cxx
   )
 
 SET(
diff --git a/plugins/fpa/GradientBaseImageFunctionSource.cxx b/plugins/fpa/GradientBaseImageFunctionSource.cxx
new file mode 100644 (file)
index 0000000..129968e
--- /dev/null
@@ -0,0 +1,110 @@
+#include "GradientBaseImageFunctionSource.h"
+#include <cpPlugins/Image.h>
+#include <fpa/Image/Functors/GulsunTekMedialness.h>
+#include <fpa/Image/Functors/FluxMedialness.h>
+#include <fpa/Image/Functors/MFluxMedialness.h>
+#include <cpExtensions/Algorithms/ImageFunctionFilter.h>
+#include <fpa/Image/Functors/GulsunTekMedialness.hxx>
+#include <fpa/Image/Functors/FluxMedialness.hxx>
+#include <fpa/Image/Functors/MFluxMedialness.hxx>
+#include <fpa/Image/Functors/GradientImageFunctionBase.hxx>
+#include <cpExtensions/Algorithms/ImageFunctionFilter.hxx>
+#include <itkImageConstIteratorWithIndex.hxx>
+#include <itkImageIteratorWithIndex.hxx>
+#include <itkImageToImageFilter.hxx>
+#include <itkImageFunction.hxx>
+#include <itkImageRegionConstIteratorWithIndex.hxx>
+
+// -------------------------------------------------------------------------
+fpaPlugins::GradientBaseImageFunctionSource::
+GradientBaseImageFunctionSource( )
+  : Superclass( )
+{
+  this->_AddInput( "Input" );
+  this->_AddOutput< cpPlugins::Image >( "Output" );
+
+  std::vector< std::string > choices;
+  choices.push_back( "Gulsun&Tek" );
+  choices.push_back( "Flux" );
+  choices.push_back( "MFlux" );
+  this->m_Parameters.ConfigureAsChoices( "FunctionType", choices );
+  this->m_Parameters.SetSelectedChoice( "FunctionType", "Gulsun&Tek" );
+
+  this->m_Parameters.ConfigureAsReal( "MinRadius" );
+  this->m_Parameters.ConfigureAsReal( "MaxRadius" );
+  this->m_Parameters.ConfigureAsUint( "ProfileSampling" );
+  this->m_Parameters.ConfigureAsUint( "RadialSampling" );
+  this->m_Parameters.SetReal( "MinRadius", 0 );
+  this->m_Parameters.SetReal( "MaxRadius", 1 );
+  this->m_Parameters.SetUint( "ProfileSampling", 4 );
+  this->m_Parameters.SetUint( "RadialSampling", 10 );
+}
+
+// -------------------------------------------------------------------------
+fpaPlugins::GradientBaseImageFunctionSource::
+~GradientBaseImageFunctionSource( )
+{
+}
+
+// -------------------------------------------------------------------------
+std::string fpaPlugins::GradientBaseImageFunctionSource::
+_GenerateData( )
+{
+  auto image = this->GetInputData( "Input" )->GetITK< itk::DataObject >( );
+  std::string   cpPlugin_Image_Demangle_VectorPixel_AllFloats( r, _GD0, image, itk::CovariantVector, 2 );
+  if( r != "" ) cpPlugin_Image_Demangle_VectorPixel_AllFloats( r, _GD0, image, itk::CovariantVector, 3 );
+  return( r );
+}
+
+// -------------------------------------------------------------------------
+template< class _TImage >
+std::string fpaPlugins::GradientBaseImageFunctionSource::
+_GD0( _TImage* image )
+{
+  typedef fpa::Image::Functors::GulsunTekMedialness< _TImage > _TGT;
+  typedef fpa::Image::Functors::FluxMedialness< _TImage > _TFl;
+  typedef fpa::Image::Functors::MFluxMedialness< _TImage > _TMFl;
+
+  if( image == NULL )
+    return( "GradientBaseImageFunctionSource: Invalid input image." );
+
+  auto ft = this->m_Parameters.GetSelectedChoice( "FunctionType" );
+  if     ( ft == "Gulsun&Tek" ) return( this->_GD1< _TImage, _TGT >( image ) );
+  else if( ft == "Flux" )       return( this->_GD1< _TImage, _TFl >( image ) );
+  else if( ft == "MFlux" )      return( this->_GD1< _TImage, _TMFl >( image ) );
+  else
+    return( "GradientBaseImageFunctionSource: Invalid function type." );
+}
+
+// -------------------------------------------------------------------------
+template< class _TImage, class _TFunction >
+std::string fpaPlugins::GradientBaseImageFunctionSource::
+_GD1( _TImage* image )
+{
+  typedef itk::Image< typename _TFunction::TOutput, _TImage::ImageDimension > _TOutImage;
+  typedef cpExtensions::Algorithms::ImageFunctionFilter< _TImage, _TOutImage, _TFunction > _TFilter;
+
+  _TFilter* filter = this->_CreateITK< _TFilter >( );
+  filter->SetInput( image );
+  _TFunction* function = filter->GetFunction( );
+  if( function == NULL )
+  {
+    filter->SetFunction( _TFunction::New( ) );
+    function = filter->GetFunction( );
+
+  } // fi
+  function->SetMinRadius( this->m_Parameters.GetReal( "MinRadius" ) );
+  function->SetMaxRadius( this->m_Parameters.GetReal( "MaxRadius" ) );
+  /*
+    function->SetProfileSampling( this->m_Parameters.GetUint( "ProfileSampling" ) );
+    function->SetRadialSampling( this->m_Parameters.GetUint( "RadialSampling" ) );
+  */
+
+  filter->Update( );
+
+  // Connect output and finish
+  this->GetOutputData( "Output" )->SetITK( filter->GetOutput( ) );
+  return( "" );
+}
+
+// eof - $RCSfile$
diff --git a/plugins/fpa/GradientBaseImageFunctionSource.h b/plugins/fpa/GradientBaseImageFunctionSource.h
new file mode 100644 (file)
index 0000000..97862a2
--- /dev/null
@@ -0,0 +1,51 @@
+#ifndef __FPAPLUGINS__GRADIENTBASEIMAGEFUNCTIONSOURCE__H__
+#define __FPAPLUGINS__GRADIENTBASEIMAGEFUNCTIONSOURCE__H__
+
+#include <fpa/fpaPlugins_Export.h>
+#include <cpPlugins/ProcessObject.h>
+
+namespace fpaPlugins
+{
+  /**
+   */
+  class fpaPlugins_EXPORT GradientBaseImageFunctionSource
+    : public cpPlugins::ProcessObject
+  {
+  public:
+    typedef GradientBaseImageFunctionSource Self;
+    typedef cpPlugins::ProcessObject        Superclass;
+    typedef itk::SmartPointer< Self >       Pointer;
+    typedef itk::SmartPointer< const Self > ConstPointer;
+
+  public:
+    itkNewMacro( Self );
+    itkTypeMacro(
+      GradientBaseImageFunctionSource, cpPlugins::ProcessObject
+      );
+    cpPlugins_Id_Macro(
+      GradientBaseImageFunctionSource, fpaImageAlgorithmFunctors
+      );
+
+  protected:
+    GradientBaseImageFunctionSource( );
+    virtual ~GradientBaseImageFunctionSource( );
+
+    virtual std::string _GenerateData( ) ITK_OVERRIDE;
+
+    template< class _TImage >
+      inline std::string _GD0( _TImage* image );
+
+    template< class _TImage, class _TFunction >
+      inline std::string _GD1( _TImage* image );
+
+  private:
+    // Purposely not implemented.
+    GradientBaseImageFunctionSource( const Self& other );
+    Self& operator=( const Self& other );
+  };
+
+} // ecapseman
+
+#endif // __FPAPLUGINS__GRADIENTBASEIMAGEFUNCTIONSOURCE__H__
+
+// eof - $RCSfile$