## == Functions, packages and options ==
## =====================================
-#INCLUDE(cmake/BaseConfig.cmake)
+INCLUDE(cmake/BaseConfig.cmake)
+INCLUDE(cmake/KitwareTools.cmake)
#OPTION(USE_cpPlugins "Build cpPlugins-based code" OFF)
#IF(USE_cpPlugins)
# FIND_PACKAGE(cpPlugins)
# MARK_AS_ADVANCED(CLEAR cpPlugins_DIR)
#ENDIF(USE_cpPlugins)
-#INCLUDE(cmake/KitwareTools.cmake)
#INCLUDE(cmake/QtTools.cmake)
#INCLUDE(cmake/Functions.cmake)
--- /dev/null
+## =======================================================================
+## == Force c++11 language version ==
+## == NOTE: It seems that by default on Visual Studio Compiler supports ==
+## == c++11, so it only need to be tested on other OS. ==
+## =======================================================================
+
+IF(NOT MSVC)
+ INCLUDE(CheckCXXCompilerFlag)
+ CHECK_CXX_COMPILER_FLAG("-std=c++11" COMPILER_SUPPORTS_CXX11)
+ IF(COMPILER_SUPPORTS_CXX11)
+ SET(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=c++11")
+ ELSE(COMPILER_SUPPORTS_CXX11)
+ CHECK_CXX_COMPILER_FLAG("-std=c++0x" COMPILER_SUPPORTS_CXX0X)
+ IF(COMPILER_SUPPORTS_CXX0X)
+ SET(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=c++0x")
+ ELSE(COMPILER_SUPPORTS_CXX0X)
+ MESSAGE(
+ FATAL_ERROR
+ "The compiler ${CMAKE_CXX_COMPILER} has no C++11 support."
+ )
+ ENDIF(COMPILER_SUPPORTS_CXX0X)
+ ENDIF(COMPILER_SUPPORTS_CXX11)
+ENDIF(NOT MSVC)
+
+## ===================================================
+## == Prepare header generator to build shared libs ==
+## ===================================================
+
+INCLUDE(GenerateExportHeader)
+
+## ==================================================
+## == Do not allow to build inside the source tree ==
+## ==================================================
+
+IF(PROJECT_BINARY_DIR STREQUAL ${PROJECT_SOURCE_DIR})
+ MESSAGE(FATAL_ERROR "Building in the source tree is not allowed.")
+ENDIF(PROJECT_BINARY_DIR STREQUAL ${PROJECT_SOURCE_DIR})
+
+## =================================================
+## == Where to put targets (executables and libs) ==
+## =================================================
+
+SET(EXECUTABLE_OUTPUT_PATH ${PROJECT_BINARY_DIR})
+SET(LIBRARY_OUTPUT_PATH ${PROJECT_BINARY_DIR})
+MARK_AS_ADVANCED(
+ CMAKE_BACKWARDS_COMPATIBILITY
+ EXECUTABLE_OUTPUT_PATH
+ LIBRARY_OUTPUT_PATH
+ )
+
+## eof - $RCSfile$
--- /dev/null
+# ========================
+# == Find Kitware tools ==
+# ========================
+
+FIND_PACKAGE(ITK REQUIRED)
+INCLUDE(${ITK_USE_FILE})
+
+## eof - $RCSfile$
--- /dev/null
+#include <climits>
+#include <cstdlib>
+#include <fstream>
+#include <iostream>
+#include <string>
+
+#include <itkImage.h>
+#include <itkImageFileReader.h>
+#include <itkImageFileWriter.h>
+
+#include <fpa/Image/MoriFilter.h>
+
+// -------------------------------------------------------------------------
+static const unsigned int Dim = 3;
+typedef short TPixel;
+typedef itk::Image< TPixel, Dim > TImage;
+
+// -------------------------------------------------------------------------
+std::string GetFullPath( const std::string& filename )
+{
+ /* On windows:
+ #include <windows.h>
+ TCHAR full_path[MAX_PATH];
+ GetFullPathName(_T("foo.dat"), MAX_PATH, full_path, NULL);
+ */
+
+ char* buffer = realpath( filename.c_str( ), NULL );
+ std::string path = buffer;
+ std::free( buffer );
+ return( path );
+}
+
+// -------------------------------------------------------------------------
+std::string GetFullPathToDirectory( const std::string& filename )
+{
+ std::string path = GetFullPath( filename );
+ std::size_t found = path.find_last_of( "/\\" );
+ return( path.substr( 0, found + 1 ) );
+}
+
+// -------------------------------------------------------------------------
+int main( int argc, char* argv[] )
+{
+ // Command configuration
+ if( argc < 3 )
+ {
+ std::cerr
+ << "Usage: " << argv[ 0 ] << " input_image output_image"
+ << std::endl;
+ return( 1 );
+
+ } // fi
+ std::string input_image_filename = GetFullPath( argv[ 1 ] );
+ std::string output_image_filename = argv[ 2 ];
+
+ // Try to guess initial seed
+ std::string seed_filename = GetFullPathToDirectory( argv[ 1 ] ) + "seed.txt";
+ std::ifstream seed_str( seed_filename.c_str( ) );
+ if( !seed_str )
+ {
+ std::cerr
+ << "No \"seed.txt\" file in the same input image directory."
+ << std::endl;
+ return( 1 );
+
+ } // fi
+ TImage::IndexType seed;
+ seed_str >> seed[ 0 ] >> seed[ 1 ] >> seed[ 2 ];
+
+ // Read image
+ typedef itk::ImageFileReader< TImage > TReader;
+ TReader::Pointer reader = TReader::New( );
+ reader->SetFileName( input_image_filename );
+
+ // Mori's algorithms
+ typedef fpa::Image::MoriFilter< TImage, TImage > TFilter;
+ TFilter::Pointer filter = TFilter::New( );
+ filter->SetInput( reader->GetOutput( ) );
+ filter->SetSeed( seed );
+ filter->SetThresholdRange( -1024, 0 );
+ filter->SetInsideValue( 1 );
+ filter->SetOutsideValue( 0 );
+
+ // Write results
+ typedef itk::ImageFileWriter< TImage > TWriter;
+ TWriter::Pointer writer = TWriter::New( );
+ writer->SetInput( filter->GetOutput( ) );
+ writer->SetFileName( output_image_filename );
+
+ // Execute pipeline
+ try
+ {
+ writer->Update( );
+ }
+ catch( std::exception& err )
+ {
+ std::cerr << "Error caught: " << err.what( ) << std::endl;
+ return( 1 );
+
+ } // yrt
+ return( 0 );
+}
+
+// eof - $RCSfile$
--- /dev/null
+OPTION(BUILD_EXAMPLES "Build cpPlugins-free examples." OFF)
+IF(BUILD_EXAMPLES)
+ SET(
+ _examples
+ BronchiiInitialSegmentationWithMori
+ )
+ FOREACH(_e ${_examples})
+ ADD_EXECUTABLE(fpa_example_${_e} ${_e}.cxx)
+ TARGET_LINK_LIBRARIES(fpa_example_${_e} ${ITK_LIBRARIES})
+ ENDFOREACH(_e)
+ENDIF(BUILD_EXAMPLES)
+
+## eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Base__Algorithm__h__
+#define __fpa__Base__Algorithm__h__
+
+#include <vector>
+#include <itkFunctionBase.h>
+#include <fpa/Config.h>
+#include <fpa/Base/Events.h>
+#include <fpa/Base/Functors/VertexCostFunctionBase.h>
+
+namespace fpa
+{
+ namespace Base
+ {
+ /**
+ */
+ template< class _TFilter, class _TVertex, class _TOutput >
+ class Algorithm
+ : public _TFilter
+ {
+ public:
+ typedef Algorithm Self;
+ typedef _TFilter Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ typedef _TVertex TVertex;
+ typedef _TOutput TOutput;
+ typedef unsigned long TFrontId;
+
+ // Different functions
+ typedef std::vector< TVertex > TNeighborhood;
+ typedef itk::FunctionBase< TVertex, TNeighborhood > TNeighborhoodFunction;
+ typedef itk::FunctionBase< TOutput, TOutput > TConversionFunction;
+ typedef fpa::Base::Functors::VertexCostFunctionBase< TVertex, TOutput > TVertexFunction;
+
+ // Minigraph to represent collisions
+ typedef std::pair< _TVertex, bool > TCollision;
+ typedef std::vector< TCollision > TCollisionsRow;
+ typedef std::vector< TCollisionsRow > TCollisions;
+
+ protected:
+ struct _TQueueNode
+ {
+ TVertex Vertex;
+ TVertex Parent;
+ TOutput Result;
+ TFrontId FrontId;
+ _TQueueNode( );
+ _TQueueNode( const TVertex& v );
+ _TQueueNode( const TVertex& v, const _TQueueNode& n );
+ };
+
+ public:
+ itkTypeMacro( Self, _TFilter );
+
+ fpa_Base_NewEvent( TStartEvent );
+ fpa_Base_NewEvent( TEndEvent );
+ fpa_Base_NewEvent( TStartLoopEvent );
+ fpa_Base_NewEvent( TEndLoopEvent );
+ fpa_Base_NewEventWithVertex( TPushEvent, TVertex );
+ fpa_Base_NewEventWithVertex( TPopEvent, TVertex );
+ fpa_Base_NewEventWithVertex( TMarkEvent, TVertex );
+
+ itkGetConstMacro( InitResult, _TOutput );
+ itkSetMacro( InitResult, _TOutput );
+
+ itkBooleanMacro( StopAtOneFront );
+ itkGetConstMacro( StopAtOneFront, bool );
+ itkSetMacro( StopAtOneFront, bool );
+
+ itkGetObjectMacro( NeighborhoodFunction, TNeighborhoodFunction );
+ itkGetObjectMacro( ConversionFunction, TConversionFunction );
+ itkGetObjectMacro( VertexFunction, TVertexFunction );
+
+ itkGetConstObjectMacro( NeighborhoodFunction, TNeighborhoodFunction );
+ itkGetConstObjectMacro( ConversionFunction, TConversionFunction );
+ itkGetConstObjectMacro( VertexFunction, TVertexFunction );
+
+ itkSetObjectMacro( NeighborhoodFunction, TNeighborhoodFunction );
+ itkSetObjectMacro( ConversionFunction, TConversionFunction );
+ itkSetObjectMacro( VertexFunction, TVertexFunction );
+
+ public:
+ void ClearSeeds( );
+ void AddSeed( const TVertex& seed, const TOutput& value );
+
+ protected:
+ Algorithm( );
+ virtual ~Algorithm( );
+
+ virtual void GenerateData( ) override;
+
+ // Particular methods
+ virtual bool _ContinueGenerateData( );
+ virtual void _Loop( );
+
+ virtual void _BeforeGenerateData( );
+ virtual void _AfterGenerateData( );
+ virtual void _BeforeLoop( );
+ virtual void _AfterLoop( );
+
+ virtual bool _ValidLoop( ) const;
+ virtual void _UpdateCollisions( const TVertex& a, const TVertex& b );
+
+ virtual _TOutput _GetInputValue( const _TQueueNode& v, const _TQueueNode& p );
+
+ virtual void _InitMarks( ) = 0;
+ virtual void _InitResults( const TOutput& init_value ) = 0;
+ virtual bool _IsMarked( const _TVertex& v ) const = 0;
+ virtual void _Mark( const _TQueueNode& n ) = 0;
+ virtual TFrontId _GetMark( const _TVertex& v ) const = 0;
+ virtual void _UpdateResult( const _TQueueNode& n ) = 0;
+ virtual TOutput _GetResult( const _TVertex& v ) const = 0;
+ virtual unsigned int _GetNumberOfDimensions( ) const = 0;
+
+ virtual bool _UpdateValue( _TQueueNode& v, const _TQueueNode& p ) = 0;
+ virtual unsigned long _QueueSize( ) const = 0;
+ virtual void _QueueClear( ) = 0;
+ virtual void _QueuePush( const _TQueueNode& node ) = 0;
+ virtual _TQueueNode _QueuePop( ) = 0;
+
+ private:
+ // Purposely not implemented.
+ Algorithm( const Self& other );
+ Self& operator=( const Self& other );
+
+ protected:
+ _TOutput m_InitResult;
+ bool m_StopAtOneFront;
+
+ typename TNeighborhoodFunction::Pointer m_NeighborhoodFunction;
+ typename TConversionFunction::Pointer m_ConversionFunction;
+ typename TVertexFunction::Pointer m_VertexFunction;
+
+ std::vector< _TQueueNode > m_Seeds;
+ TCollisions m_Collisions;
+ unsigned int m_NumberOfFronts;
+ };
+
+ } // ecaseman
+
+} // ecaseman
+
+#ifndef ITK_MANUAL_INSTANTIATION
+# include <fpa/Base/Algorithm.hxx>
+#endif // ITK_MANUAL_INSTANTIATION
+
+#endif // __fpa__Base__Algorithm__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Base__Algorithm__hxx__
+#define __fpa__Base__Algorithm__hxx__
+
+#include <queue>
+
+// -------------------------------------------------------------------------
+template < class _TFilter, class _TVertex, class _TOutput >
+fpa::Base::Algorithm< _TFilter, _TVertex, _TOutput >::_TQueueNode::
+_TQueueNode( )
+{
+ this->Vertex.Fill( 0 );
+ this->Parent.Fill( 0 );
+ this->Result = _TOutput( 0 );
+ this->FrontId = 0;
+}
+
+// -------------------------------------------------------------------------
+template < class _TFilter, class _TVertex, class _TOutput >
+fpa::Base::Algorithm< _TFilter, _TVertex, _TOutput >::_TQueueNode::
+_TQueueNode( const _TVertex& v )
+{
+ this->Vertex = v;
+ this->Parent = v;
+ this->Result = _TOutput( 0 );
+ this->FrontId = 0;
+}
+
+// -------------------------------------------------------------------------
+template < class _TFilter, class _TVertex, class _TOutput >
+fpa::Base::Algorithm< _TFilter, _TVertex, _TOutput >::_TQueueNode::
+_TQueueNode( const _TVertex& v, const _TQueueNode& n )
+{
+ this->Vertex = v;
+ this->Parent = n.Vertex;
+ this->Result = n.Result;
+ this->FrontId = n.FrontId;
+}
+
+// -------------------------------------------------------------------------
+template < class _TFilter, class _TVertex, class _TOutput >
+void fpa::Base::Algorithm< _TFilter, _TVertex, _TOutput >::
+ClearSeeds( )
+{
+ this->m_Seeds.clear( );
+ this->Modified( );
+}
+
+// -------------------------------------------------------------------------
+template < class _TFilter, class _TVertex, class _TOutput >
+void fpa::Base::Algorithm< _TFilter, _TVertex, _TOutput >::
+AddSeed( const _TVertex& seed, const TOutput& value )
+{
+ _TQueueNode node( seed );
+ node.FrontId = this->m_Seeds.size( ) + 1;
+ node.Result = value;
+ this->m_Seeds.push_back( node );
+ this->Modified( );
+}
+
+// -------------------------------------------------------------------------
+template < class _TFilter, class _TVertex, class _TOutput >
+fpa::Base::Algorithm< _TFilter, _TVertex, _TOutput >::
+Algorithm( )
+ : Superclass( ),
+ m_InitResult( _TOutput( 0 ) ),
+ m_StopAtOneFront( false )
+{
+}
+
+// -------------------------------------------------------------------------
+template < class _TFilter, class _TVertex, class _TOutput >
+fpa::Base::Algorithm< _TFilter, _TVertex, _TOutput >::
+~Algorithm( )
+{
+}
+
+// -------------------------------------------------------------------------
+template < class _TFilter, class _TVertex, class _TOutput >
+void fpa::Base::Algorithm< _TFilter, _TVertex, _TOutput >::
+GenerateData( )
+{
+ this->InvokeEvent( TStartEvent( ) );
+ this->_BeforeGenerateData( );
+ this->m_NumberOfFronts = this->m_Seeds.size( );
+ if( this->m_NumberOfFronts > 0 )
+ {
+ // Prepare collisions
+ TCollision coll( TVertex( ), false );
+ TCollisionsRow row( this->m_NumberOfFronts, coll );
+ this->m_Collisions.clear( );
+ this->m_Collisions.resize( this->m_NumberOfFronts, row );
+
+ // Put seeds on queue
+ this->_QueueClear( );
+ for(
+ auto nIt = this->m_Seeds.begin( );
+ nIt != this->m_Seeds.end( );
+ ++nIt
+ )
+ this->_QueuePush( *nIt );
+
+ // Init marks and results
+ this->_InitMarks( );
+ this->_InitResults( this->m_InitResult );
+
+ // Main loop
+ do
+ this->_Loop( );
+ while( this->_ContinueGenerateData( ) );
+
+ } // fi
+ this->_AfterGenerateData( );
+ this->InvokeEvent( TEndEvent( ) );
+}
+
+// -------------------------------------------------------------------------
+template < class _TFilter, class _TVertex, class _TOutput >
+bool fpa::Base::Algorithm< _TFilter, _TVertex, _TOutput >::
+_ContinueGenerateData( )
+{
+ return( false );
+}
+
+// -------------------------------------------------------------------------
+template < class _TFilter, class _TVertex, class _TOutput >
+void fpa::Base::Algorithm< _TFilter, _TVertex, _TOutput >::
+_Loop( )
+{
+ this->InvokeEvent( TStartLoopEvent( ) );
+ this->_BeforeLoop( );
+ while( this->_ValidLoop( ) )
+ {
+ _TQueueNode node = this->_QueuePop( );
+ this->InvokeEvent( TPopEvent( node.Vertex, node.FrontId ) );
+
+ // Process actual vertex
+ if( this->_IsMarked( node.Vertex ) )
+ continue;
+ this->_Mark( node );
+ this->_UpdateResult( node );
+ this->InvokeEvent( TMarkEvent( node.Vertex, node.FrontId ) );
+
+ // Add neighbors to queue
+ TNeighborhood neighs = this->m_NeighborhoodFunction->Evaluate( node.Vertex );
+ for( auto it = neighs.begin( ); it != neighs.end( ); ++it )
+ {
+ if( !( this->_IsMarked( *it ) ) )
+ {
+ _TQueueNode neigh( *it, node );
+ if( this->_UpdateValue( neigh, node ) )
+ {
+ this->_QueuePush( neigh );
+ this->InvokeEvent( TPushEvent( node.Vertex, node.FrontId ) );
+
+ } // fi
+ }
+ else
+ this->_UpdateCollisions( node.Vertex, *it );
+
+ } // rof
+
+ } // elihw
+ this->_AfterLoop( );
+ this->InvokeEvent( TEndLoopEvent( ) );
+}
+
+// -------------------------------------------------------------------------
+template < class _TFilter, class _TVertex, class _TOutput >
+void fpa::Base::Algorithm< _TFilter, _TVertex, _TOutput >::
+_BeforeGenerateData( )
+{
+}
+
+// -------------------------------------------------------------------------
+template < class _TFilter, class _TVertex, class _TOutput >
+void fpa::Base::Algorithm< _TFilter, _TVertex, _TOutput >::
+_AfterGenerateData( )
+{
+}
+
+// -------------------------------------------------------------------------
+template < class _TFilter, class _TVertex, class _TOutput >
+void fpa::Base::Algorithm< _TFilter, _TVertex, _TOutput >::
+_BeforeLoop( )
+{
+}
+
+// -------------------------------------------------------------------------
+template < class _TFilter, class _TVertex, class _TOutput >
+void fpa::Base::Algorithm< _TFilter, _TVertex, _TOutput >::
+_AfterLoop( )
+{
+}
+
+// -------------------------------------------------------------------------
+template < class _TFilter, class _TVertex, class _TOutput >
+bool fpa::Base::Algorithm< _TFilter, _TVertex, _TOutput >::
+_ValidLoop( ) const
+{
+ bool r = ( this->_QueueSize( ) > 0 );
+ if( this->m_StopAtOneFront )
+ r &= ( this->m_NumberOfFronts > 0 );
+ return( r );
+}
+
+// -------------------------------------------------------------------------
+template < class _TFilter, class _TVertex, class _TOutput >
+void fpa::Base::Algorithm< _TFilter, _TVertex, _TOutput >::
+_UpdateCollisions( const TVertex& a, const TVertex& b )
+{
+ auto ma = this->_GetMark( a );
+ auto mb = this->_GetMark( b );
+ if( ma == mb || ma == 0 || mb == 0 )
+ return;
+
+ // Mark collision, if it is new
+ ma--; mb--;
+ bool ret = false;
+ bool exists = this->m_Collisions[ ma ][ mb ].second;
+ exists &= this->m_Collisions[ mb ][ ma ].second;
+ if( !exists )
+ {
+ this->m_Collisions[ ma ][ mb ].first = a;
+ this->m_Collisions[ ma ][ mb ].second = true;
+ this->m_Collisions[ mb ][ ma ].first = b;
+ this->m_Collisions[ mb ][ ma ].second = true;
+
+ // Update number of fronts
+ unsigned long N = this->m_Seeds.size( );
+ unsigned long count = 0;
+ std::vector< bool > m( N, false );
+ std::queue< unsigned long > q;
+ q.push( 0 );
+ while( !q.empty( ) )
+ {
+ unsigned long f = q.front( );
+ q.pop( );
+
+ if( m[ f ] )
+ continue;
+ m[ f ] = true;
+ count++;
+
+ for( unsigned int n = 0; n < N; ++n )
+ if( this->m_Collisions[ f ][ n ].second && !m[ n ] )
+ q.push( n );
+
+ } // elihw
+ this->m_NumberOfFronts = N - count;
+
+ } // fi
+}
+
+// -------------------------------------------------------------------------
+template < class _TFilter, class _TVertex, class _TOutput >
+_TOutput fpa::Base::Algorithm< _TFilter, _TVertex, _TOutput >::
+_GetInputValue( const _TQueueNode& v, const _TQueueNode& p )
+{
+ _TOutput res = this->m_InitResult;
+ if( this->m_VertexFunction.IsNotNull( ) )
+ res = this->m_VertexFunction->Evaluate( v.Vertex, p.Vertex );
+ if( this->m_ConversionFunction.IsNotNull( ) )
+ res = this->m_ConversionFunction->Evaluate( res );
+ return( res );
+}
+
+#endif // __fpa__Base__Algorithm__hxx__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Base__Dijkstra__h__
+#define __fpa__Base__Dijkstra__h__
+
+#include <fpa/Base/PriorityQueueAlgorithm.h>
+
+namespace fpa
+{
+ namespace Base
+ {
+ /**
+ */
+ template< class _TSuperclass, class _TMST >
+ class Dijkstra
+ : public fpa::Base::PriorityQueueAlgorithm< _TSuperclass >
+ {
+ public:
+ typedef Dijkstra Self;
+ typedef fpa::Base::PriorityQueueAlgorithm< _TSuperclass > Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ typedef _TMST TMST;
+ typedef typename Superclass::TOutput TOutput;
+ typedef typename Superclass::TVertex TVertex;
+
+ protected:
+ typedef typename Superclass::_TQueueNode _TQueueNode;
+
+ public:
+ itkTypeMacro( Dijkstra, Algorithm );
+
+ public:
+ _TMST* GetMinimumSpanningTree( );
+ const _TMST* GetMinimumSpanningTree( ) const;
+
+ protected:
+ Dijkstra( );
+ virtual ~Dijkstra( );
+
+ virtual void _AfterGenerateData( ) override;
+
+ virtual void _UpdateResult( const _TQueueNode& n ) override;
+ virtual bool _UpdateValue(
+ _TQueueNode& v, const _TQueueNode& p
+ ) override;
+
+ private:
+ // Purposely not defined
+ Dijkstra( const Self& other );
+ Self& operator=( const Self& other );
+
+ protected:
+ unsigned long m_MSTIndex;
+ };
+
+ } // ecapseman
+
+} // ecapseman
+
+#ifndef ITK_MANUAL_INSTANTIATION
+# include <fpa/Base/Dijkstra.hxx>
+#endif // ITK_MANUAL_INSTANTIATION
+
+#endif // __fpa__Base__Dijkstra__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Base__Dijkstra__hxx__
+#define __fpa__Base__Dijkstra__hxx__
+
+#include <algorithm>
+#include <limits>
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass, class _TMST >
+_TMST* fpa::Base::Dijkstra< _TSuperclass, _TMST >::
+GetMinimumSpanningTree( )
+{
+ return(
+ dynamic_cast< _TMST* >(
+ this->itk::ProcessObject::GetOutput( this->m_MSTIndex )
+ )
+ );
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass, class _TMST >
+const _TMST* fpa::Base::Dijkstra< _TSuperclass, _TMST >::
+GetMinimumSpanningTree( ) const
+{
+ return(
+ dynamic_cast< const _TMST* >(
+ this->itk::ProcessObject::GetOutput( this->m_MSTIndex )
+ )
+ );
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass, class _TMST >
+fpa::Base::Dijkstra< _TSuperclass, _TMST >::
+Dijkstra( )
+ : Superclass( )
+{
+ this->m_InitResult = TOutput( 0 );
+ this->m_MSTIndex = this->GetNumberOfRequiredOutputs( );
+ this->SetNumberOfRequiredOutputs( this->m_MSTIndex + 1 );
+ this->itk::ProcessObject::SetNthOutput( this->m_MSTIndex, _TMST::New( ) );
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass, class _TMST >
+fpa::Base::Dijkstra< _TSuperclass, _TMST >::
+~Dijkstra( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass, class _TMST >
+void fpa::Base::Dijkstra< _TSuperclass, _TMST >::
+_AfterGenerateData( )
+{
+ this->Superclass::_AfterGenerateData( );
+
+ auto mst = this->GetMinimumSpanningTree( );
+ mst->ClearSeeds( );
+ for( auto s = this->m_Seeds.begin( ); s != this->m_Seeds.end( ); ++s )
+ mst->AddSeed( s->Vertex );
+ mst->SetCollisions( this->m_Collisions );
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass, class _TMST >
+void fpa::Base::Dijkstra< _TSuperclass, _TMST >::
+_UpdateResult( const _TQueueNode& n )
+{
+ this->Superclass::_UpdateResult( n );
+ this->GetMinimumSpanningTree( )->SetParent( n.Vertex, n.Parent );
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass, class _TMST >
+bool fpa::Base::Dijkstra< _TSuperclass, _TMST >::
+_UpdateValue( _TQueueNode& v, const _TQueueNode& p )
+{
+ v.Result = this->_GetInputValue( p.Vertex, v.Vertex );
+ if( v.Result >= TOutput( 0 ) )
+ {
+ v.Result += p.Result;
+ return( true );
+ }
+ else
+ {
+ v.Result = this->m_InitResult;
+ return( false );
+
+ } // fi
+}
+
+#endif // __fpa__Base__Dijkstra__hxx__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Base__Events__h__
+#define __fpa__Base__Events__h__
+
+// -------------------------------------------------------------------------
+#define fpa_Base_NewEvent( name ) \
+ class name \
+ : public BaseEvent \
+ { \
+ public: \
+ name( ) : BaseEvent( ) { } \
+ virtual ~name( ) { } \
+ const char* GetEventName( ) const \
+ { return( "fpa::Base::##name" ); } \
+ bool CheckEvent( const itk::EventObject* e ) const \
+ { return( dynamic_cast< const name* >( e ) != NULL ); } \
+ itk::EventObject* MakeObject( ) const \
+ { return( new name( ) ); } \
+ };
+
+// -------------------------------------------------------------------------
+#define fpa_Base_NewEventWithVertex( name, type ) \
+ class name \
+ : public BaseEventWithVertex< type > \
+ { \
+ public: \
+ name( ) : BaseEventWithVertex< type >( ) { } \
+ name( const type& v, long fid ) \
+ : BaseEventWithVertex< type >( v, fid ) { } \
+ virtual ~name( ) { } \
+ const char* GetEventName( ) const \
+ { return( "fpa::Base::##name" ); } \
+ bool CheckEvent( const itk::EventObject* e ) const \
+ { return( dynamic_cast< const name* >( e ) != NULL ); } \
+ itk::EventObject* MakeObject( ) const \
+ { return( new name( ) ); } \
+ };
+
+namespace fpa
+{
+ namespace Base
+ {
+ /**
+ */
+ class BaseEvent
+ : public itk::AnyEvent
+ {
+ public:
+ BaseEvent( )
+ : itk::AnyEvent( )
+ { }
+ virtual ~BaseEvent( )
+ { }
+ const char* GetEventName( ) const
+ { return( "fpa::Base::BaseEvent" ); }
+ bool CheckEvent( const itk::EventObject* e ) const
+ { return( dynamic_cast< const BaseEvent* >( e ) != NULL ); }
+ itk::EventObject* MakeObject( ) const
+ { return( new BaseEvent( ) ); }
+ };
+
+ /**
+ */
+ template< class _TVertex >
+ class BaseEventWithVertex
+ : public BaseEvent
+ {
+ public:
+ BaseEventWithVertex( )
+ : BaseEvent( )
+ { }
+ BaseEventWithVertex( const _TVertex& v, long fid )
+ : BaseEvent( ),
+ Vertex( v ),
+ FrontId( fid )
+ { }
+ virtual ~BaseEventWithVertex( )
+ { }
+ const char* GetEventName( ) const
+ { return( "fpa::Base::BaseEventWithVertex" ); }
+ bool CheckEvent( const itk::EventObject* e ) const
+ {
+ return(
+ dynamic_cast< const BaseEventWithVertex< _TVertex >* >( e ) != NULL
+ );
+ }
+ itk::EventObject* MakeObject( ) const
+ { return( new BaseEventWithVertex< _TVertex >( ) ); }
+
+ public:
+ _TVertex Vertex;
+ long FrontId;
+ };
+
+ } // ecapseman
+
+} // ecapseman
+
+#endif // __fpa__Base__Events__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Base__FastMarching__h__
+#define __fpa__Base__FastMarching__h__
+
+#include <fpa/Base/PriorityQueueAlgorithm.h>
+
+namespace fpa
+{
+ namespace Base
+ {
+ /**
+ */
+ template< class _TSuperclass >
+ class FastMarching
+ : public fpa::Base::PriorityQueueAlgorithm< _TSuperclass >
+ {
+ public:
+ typedef FastMarching Self;
+ typedef fpa::Base::PriorityQueueAlgorithm< _TSuperclass > Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ typedef typename Superclass::TFrontId TFrontId;
+ typedef typename Superclass::TOutput TOutput;
+ typedef typename Superclass::TVertex TVertex;
+
+ typedef std::pair< TVertex, double > TFastMarchingNeighbor;
+ typedef std::vector< TFastMarchingNeighbor > TFastMarchingNeighborhood;
+
+ protected:
+ typedef typename Superclass::_TQueueNode _TQueueNode;
+
+ public:
+ itkTypeMacro( FastMarching, Algorithm );
+
+ protected:
+ FastMarching( );
+ virtual ~FastMarching( );
+
+ virtual TFastMarchingNeighborhood _FastMarchingNeighbors( const TVertex& v ) const = 0;
+
+ virtual bool _UpdateValue( _TQueueNode& v, const _TQueueNode& p ) override;
+
+ private:
+ // Purposely not defined
+ FastMarching( const Self& other );
+ Self& operator=( const Self& other );
+ };
+
+ } // ecapseman
+
+} // ecapseman
+
+#ifndef ITK_MANUAL_INSTANTIATION
+# include <fpa/Base/FastMarching.hxx>
+#endif // ITK_MANUAL_INSTANTIATION
+
+#endif // __fpa__Base__FastMarching__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Base__FastMarching__hxx__
+#define __fpa__Base__FastMarching__hxx__
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass >
+fpa::Base::FastMarching< _TSuperclass >::
+FastMarching( )
+ : Superclass( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass >
+fpa::Base::FastMarching< _TSuperclass >::
+~FastMarching( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass >
+bool fpa::Base::FastMarching< _TSuperclass >::
+_UpdateValue( _TQueueNode& v, const _TQueueNode& p )
+{
+ // Compute quadratic coefficients
+ double a = double( 0 );
+ double b = double( 0 );
+ double c = double( 0 );
+ TFastMarchingNeighborhood neighs = this->_FastMarchingNeighbors( v.Vertex );
+ for( unsigned int i = 0; i < neighs.size( ); i += 2 )
+ {
+ double tn = double( this->_GetResult( neighs[ i ].first ) );
+ double tp = double( this->_GetResult( neighs[ i + 1 ].first ) );
+ double dd = double( 0 );
+ double td = this->m_InitResult;
+ if( tn < tp )
+ {
+ td = tn;
+ dd = neighs[ i ].second;
+ }
+ else
+ {
+ td = tp;
+ dd = neighs[ i + 1 ].second;
+
+ } // fi
+
+ if( td < double( this->m_InitResult ) )
+ {
+ dd = double( 1 ) / ( dd * dd );
+ a += dd;
+ b += td * dd;
+ c += td * td * dd;
+
+ } // fi
+
+ } // rof
+ double F = double( this->_GetInputValue( v, p ) );
+ c -= double( 1 ) / ( F * F );
+
+ // Solve quadratic equation
+ double d = ( b * b ) - ( a * c );
+ if( d >= double( 0 ) )
+ {
+ d = std::sqrt( d );
+ b *= double( -1 );
+ double s1 = std::fabs( ( b + d ) / a );
+ double s2 = std::fabs( ( b - d ) / a );
+ v.Result = TOutput( ( s2 < s1 )? s1: s2 );
+ if( v.Result < this->_GetResult( v.Vertex ) )
+ this->_UpdateResult( v );
+ }
+ else
+ {
+ std::cout << std::endl << "-- FM --> " << v.Vertex << " " << d << std::endl;
+ std::exit( 1 );
+ }
+ return( true );
+}
+
+#endif // __fpa__Base__FastMarching__hxx__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Base__Functors__GaussianModel__h__
+#define __fpa__Base__Functors__GaussianModel__h__
+
+#include <fpa/Config.h>
+#include <itkFunctionBase.h>
+#include <cpExtensions/Algorithms/IterativeGaussianModelEstimator.h>
+
+namespace fpa
+{
+ namespace Base
+ {
+ namespace Functors
+ {
+ /**
+ */
+ template< class _TInput, class _TOutput >
+ class GaussianModel
+ : public itk::FunctionBase< _TInput, _TOutput >
+ {
+ public:
+ typedef GaussianModel Self;
+ typedef itk::FunctionBase< _TInput, _TOutput > Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ typedef _TInput TInput;
+ typedef _TOutput TOutput;
+
+ // Model estimator
+ typedef
+ cpExtensions::Algorithms::IterativeGaussianModelEstimator< TOutput, 1 >
+ TModel;
+
+ public:
+ itkNewMacro( Self );
+ itkTypeMacro( GaussianModel, itk::FunctionBase );
+
+ /* TODO
+ itkGetConstMacro( SupportSize, unsigned int );
+ itkGetConstMacro( MinimumCost, TOutput );
+ itkGetObjectMacro( Model, TModel );
+ itkGetConstObjectMacro( Model, TModel );
+ itkSetMacro( SupportSize, unsigned int );
+ itkSetMacro( MinimumCost, TOutput );
+ */
+
+ public:
+ virtual TOutput Evaluate( const TInput& x ) const override;
+
+ protected:
+ GaussianModel( );
+ virtual ~GaussianModel( );
+
+ private:
+ // Purposely not implemented
+ GaussianModel( const Self& other );
+ Self& operator=( const Self& other );
+
+ protected:
+ /* TODO
+ unsigned int m_SupportSize;
+ TOutput m_MinimumCost;
+ typename TModel::Pointer m_Model;
+ */
+ mutable double m_S1;
+ mutable double m_S2;
+ mutable unsigned long m_N;
+ };
+
+ } // ecapseman
+
+ } // ecapseman
+
+} // ecapseman
+
+#ifndef ITK_MANUAL_INSTANTIATION
+# include <fpa/Base/Functors/GaussianModel.hxx>
+#endif // ITK_MANUAL_INSTANTIATION
+
+#endif // __fpa__Base__Functors__GaussianModel__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Base__Functors__GaussianModel__hxx__
+#define __fpa__Base__Functors__GaussianModel__hxx__
+
+// -------------------------------------------------------------------------
+template< class _TInput, class _TOutput >
+typename fpa::Base::Functors::GaussianModel< _TInput, _TOutput >::
+TOutput fpa::Base::Functors::GaussianModel< _TInput, _TOutput >::
+Evaluate( const TInput& x ) const
+{
+ double v = double( x );
+ double d = double( 0 );
+ if( this->m_N > 1 )
+ {
+ double N = double( this->m_N );
+ double m = this->m_S1 / N;
+ double s =
+ ( this->m_S2 - ( ( this->m_S1 * this->m_S1 ) / N ) ) /
+ ( N - double( 1 ) );
+ if( s > double( 0 ) )
+ {
+ d = ( v - m );
+ d *= d;
+ d /= s;
+ d = std::sqrt( d );
+ if( d <= double( 1.5 ) ) // 2sigma
+ {
+ this->m_S1 += v;
+ this->m_S2 += v * v;
+ this->m_N += 1;
+
+ } // fi
+ }
+ else
+ {
+ this->m_S1 += v;
+ this->m_S2 += v * v;
+ this->m_N += 1;
+
+ } // fi
+ }
+ else
+ {
+ this->m_S1 += v;
+ this->m_S2 += v * v;
+ this->m_N += 1;
+
+ } // fi
+ return( d );
+
+ /* TODO
+ TOutput r = TOutput( 0 );
+ if( this->m_Model->GetNumberOfSamples( ) > this->m_SupportSize )
+ {
+ r = this->m_Model->SquaredMahalanobis( TOutput( x ) );
+ if( r <= TOutput( 1 ) )
+ this->m_Model->AddSample( TOutput( x ) );
+ }
+ else
+ {
+ this->m_Model->AddSample( TOutput( x ) );
+ if( this->m_Model->GetNumberOfSamples( ) > 2 )
+ r = this->m_Model->SquaredMahalanobis( TOutput( x ) );
+
+ } // fi
+ if( r < this->m_MinimumCost )
+ return( this->m_MinimumCost );
+ else
+ return( r );
+ */
+}
+
+// -------------------------------------------------------------------------
+template< class _TInput, class _TOutput >
+fpa::Base::Functors::GaussianModel< _TInput, _TOutput >::
+GaussianModel( )
+ : Superclass( )
+ /* TODO
+ ,
+ m_SupportSize( 30 ),
+ m_MinimumCost( TOutput( 1e-10 ) )
+ */
+{
+ /* TODO
+ this->m_Model = TModel::New( );
+ */
+ this->m_S1 = double( 0 );
+ this->m_S2 = double( 0 );
+ this->m_N = 0;
+}
+
+// -------------------------------------------------------------------------
+template< class _TInput, class _TOutput >
+fpa::Base::Functors::GaussianModel< _TInput, _TOutput >::
+~GaussianModel( )
+{
+}
+
+#endif // __fpa__Base__Functors__GaussianModel__hxx__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Base__Functors__Inverse__h__
+#define __fpa__Base__Functors__Inverse__h__
+
+#include <fpa/Config.h>
+#include <itkFunctionBase.h>
+
+namespace fpa
+{
+ namespace Base
+ {
+ namespace Functors
+ {
+ /**
+ */
+ template< class _TInput, class _TOutput >
+ class Inverse
+ : public itk::FunctionBase< _TInput, _TOutput >
+ {
+ public:
+ typedef Inverse Self;
+ typedef itk::FunctionBase< _TInput, _TOutput > Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ typedef _TInput TInput;
+ typedef _TOutput TOutput;
+
+ public:
+ itkNewMacro( Self );
+ itkTypeMacro( Inverse, itk::FunctionBase );
+
+ itkGetConstMacro( NegativeValue, _TOutput );
+ itkSetMacro( NegativeValue, _TOutput );
+
+ public:
+ virtual TOutput Evaluate( const TInput& x ) const override;
+
+ protected:
+ Inverse( );
+ virtual ~Inverse( );
+
+ private:
+ // Purposely not implemented
+ Inverse( const Self& other );
+ Self& operator=( const Self& other );
+
+ protected:
+ _TOutput m_NegativeValue;
+ };
+
+ } // ecapseman
+
+ } // ecapseman
+
+} // ecapseman
+
+#ifndef ITK_MANUAL_INSTANTIATION
+# include <fpa/Base/Functors/Inverse.hxx>
+#endif // ITK_MANUAL_INSTANTIATION
+
+#endif // __fpa__Base__Functors__Inverse__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Base__Functors__Inverse__hxx__
+#define __fpa__Base__Functors__Inverse__hxx__
+
+// -------------------------------------------------------------------------
+template< class _TInput, class _TOutput >
+typename fpa::Base::Functors::Inverse< _TInput, _TOutput >::
+TOutput fpa::Base::Functors::Inverse< _TInput, _TOutput >::
+Evaluate( const TInput& x ) const
+{
+ TInput sign = TInput( ( x < TInput( 0 ) )? -1: 1 );
+ TOutput y = TOutput( 1 ) / ( TOutput( 1 ) + TOutput( x * sign ) );
+ if( sign < TInput( 0 ) )
+ return( y * this->m_NegativeValue );
+ else
+ return( y );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInput, class _TOutput >
+fpa::Base::Functors::Inverse< _TInput, _TOutput >::
+Inverse( )
+ : Superclass( ),
+ m_NegativeValue( _TOutput( -1 ) )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TInput, class _TOutput >
+fpa::Base::Functors::Inverse< _TInput, _TOutput >::
+~Inverse( )
+{
+}
+
+#endif // __fpa__Base__Functors__Inverse__hxx__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Base__Functors__RegionGrow__Base__h__
+#define __fpa__Base__Functors__RegionGrow__Base__h__
+
+#include <itkObject.h>
+#include <itkObjectFactory.h>
+#include <fpa/Base/Functors/VertexCostFunctionBase.h>
+
+namespace fpa
+{
+ namespace Base
+ {
+ namespace Functors
+ {
+ namespace RegionGrow
+ {
+ /**
+ */
+ template< class _TVertex, class _TOutput >
+ class Base
+ : public fpa::Base::Functors::VertexCostFunctionBase< _TVertex, _TOutput >
+ {
+ public:
+ typedef fpa::Base::Functors::VertexCostFunctionBase< _TVertex, _TOutput > Superclass;
+ typedef Base Self;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ typedef _TVertex TVertex;
+ typedef _TOutput TOutput;
+
+ public:
+ itkTypeMacro( Base, itk::Object );
+
+ itkGetConstMacro( InsideValue, _TOutput );
+ itkGetConstMacro( OutsideValue, _TOutput );
+
+ itkSetMacro( InsideValue, _TOutput );
+ itkSetMacro( OutsideValue, _TOutput );
+
+ public:
+ virtual TOutput Evaluate( const TVertex& a, const TVertex& b ) const = 0;
+
+ protected:
+ Base( )
+ : Superclass( ),
+ m_InsideValue( TOutput( 1 ) ),
+ m_OutsideValue( TOutput( 0 ) )
+ { }
+ virtual ~Base( )
+ { }
+
+ private:
+ // Purposely not defined
+ Base( const Self& other );
+ Self& operator=( const Self& other );
+
+ protected:
+ _TOutput m_InsideValue;
+ _TOutput m_OutsideValue;
+ };
+
+ } // ecapseman
+
+ } // ecapseman
+
+ } // ecapseman
+
+} // ecapseman
+
+#endif // __fpa__Base__Functors__RegionGrow__Base__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Base__Functors__RegionGrow__Tautology__h__
+#define __fpa__Base__Functors__RegionGrow__Tautology__h__
+
+#include <fpa/Base/Functors/RegionGrow/Base.h>
+
+namespace fpa
+{
+ namespace Base
+ {
+ namespace Functors
+ {
+ namespace RegionGrow
+ {
+ /**
+ */
+ template< class _TVertex, class _TOutput >
+ class Tautology
+ : public Base< _TVertex, _TOutput >
+ {
+ public:
+ typedef Tautology Self;
+ typedef Base< _TVertex, _TOutput > Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ typedef typename Superclass::TOutput TOutput;
+ typedef typename Superclass::TVertex TVertex;
+
+ public:
+ itkNewMacro( Self );
+ itkTypeMacro( Tautology, Base );
+
+ public:
+ virtual TOutput Evaluate( const TVertex& a, const TVertex& b ) const override
+ {
+ return( this->m_InsideValue );
+ }
+
+ protected:
+ Tautology( )
+ : Superclass( )
+ { }
+ virtual ~Tautology( )
+ { }
+
+ private:
+ // Purposely not defined
+ Tautology( const Self& other );
+ Self& operator=( const Self& other );
+ };
+
+ } // ecapseman
+
+ } // ecapseman
+
+ } // ecapseman
+
+} // ecapseman
+
+#endif // __fpa__Base__Functors__RegionGrow__Tautology__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Base__Functors__VertexCostFunctionBase__h__
+#define __fpa__Base__Functors__VertexCostFunctionBase__h__
+
+#include <itkObject.h>
+#include <itkObjectFactory.h>
+
+namespace fpa
+{
+ namespace Base
+ {
+ namespace Functors
+ {
+ /**
+ */
+ template< class _TVertex, class _TOutput >
+ class VertexCostFunctionBase
+ : public itk::Object
+ {
+ public:
+ typedef VertexCostFunctionBase Self;
+ typedef itk::Object Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ typedef _TVertex TVertex;
+ typedef _TOutput TOutput;
+
+ public:
+ itkTypeMacro( VertexCostFunctionBase, itk::Object );
+
+ public:
+ virtual TOutput Evaluate( const TVertex& a, const TVertex& b ) const = 0;
+
+ protected:
+ VertexCostFunctionBase( )
+ : Superclass( )
+ { }
+ virtual ~VertexCostFunctionBase( )
+ { }
+
+ private:
+ // Purposely not defined
+ VertexCostFunctionBase( const Self& other );
+ Self& operator=( const Self& other );
+ };
+
+ } // ecapseman
+
+ } // ecapseman
+
+} // ecapseman
+
+#endif // __fpa__Base__Functors__VertexCostFunctionBase__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Base__MinimumSpanningTree__h__
+#define __fpa__Base__MinimumSpanningTree__h__
+
+#include <fpa/Config.h>
+#include <vector>
+#include <utility>
+#include <itkObject.h>
+
+namespace fpa
+{
+ namespace Base
+ {
+ /**
+ */
+ template< class _TVertex, class _TPath, class _TSuperclass >
+ class MinimumSpanningTree
+ : public _TSuperclass
+ {
+ public:
+ typedef MinimumSpanningTree Self;
+ typedef _TSuperclass Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ typedef _TVertex TVertex;
+ typedef _TPath TPath;
+ typedef std::pair< _TVertex, bool > TCollision;
+ typedef std::vector< TCollision > TCollisionsRow;
+ typedef std::vector< TCollisionsRow > TCollisions;
+
+ protected:
+ typedef std::vector< unsigned long > _TRow;
+ typedef std::vector< _TRow > _TMatrix;
+
+ public:
+ itkTypeMacro( MinimumSpanningTree, _TSuperclass );
+
+ public:
+ const TCollisions& GetCollisions( ) const;
+ void SetCollisions( const TCollisions& collisions );
+
+ void ClearSeeds( );
+ void AddSeed( const _TVertex& seed );
+
+ virtual _TVertex GetParent( const _TVertex& v ) const = 0;
+ virtual void SetParent( const _TVertex& v, const _TVertex& p ) = 0;
+
+ virtual void GetPath(
+ typename _TPath::Pointer& path, const _TVertex& a
+ ) const;
+ virtual void GetPath(
+ typename _TPath::Pointer& path, const _TVertex& a, const _TVertex& b
+ ) const;
+
+ protected:
+ MinimumSpanningTree( );
+ virtual ~MinimumSpanningTree( );
+
+ private:
+ // Purposely not defined
+ MinimumSpanningTree( const Self& other );
+ Self& operator=( const Self& other );
+
+ protected:
+ TCollisions m_Collisions;
+ _TMatrix m_FrontPaths;
+ std::vector< _TVertex > m_Seeds;
+ };
+
+ } // ecapseman
+
+} // ecapseman
+
+#ifndef ITK_MANUAL_INSTANTIATION
+# include <fpa/Base/MinimumSpanningTree.hxx>
+#endif // ITK_MANUAL_INSTANTIATION
+
+#endif // __fpa__Base__MinimumSpanningTree__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Base__MinimumSpanningTree__hxx__
+#define __fpa__Base__MinimumSpanningTree__hxx__
+
+#include <limits>
+
+// -------------------------------------------------------------------------
+template< class _TVertex, class _TPath, class _TSuperclass >
+const typename
+fpa::Base::MinimumSpanningTree< _TVertex, _TPath, _TSuperclass >::
+TCollisions&
+fpa::Base::MinimumSpanningTree< _TVertex, _TPath, _TSuperclass >::
+GetCollisions( ) const
+{
+ return( this->m_Collisions );
+}
+
+// -------------------------------------------------------------------------
+template< class _TVertex, class _TPath, class _TSuperclass >
+void fpa::Base::MinimumSpanningTree< _TVertex, _TPath, _TSuperclass >::
+SetCollisions( const TCollisions& collisions )
+{
+ static const unsigned long _inf =
+ std::numeric_limits< unsigned long >::max( );
+ if( this->m_Collisions == collisions )
+ return;
+
+ this->m_Collisions = collisions;
+
+ // Prepare a front graph
+ unsigned long N = this->m_Collisions.size( );
+ _TMatrix dist( N, _TRow( N, _inf ) );
+ this->m_FrontPaths = dist;
+ for( unsigned long i = 0; i < N; ++i )
+ {
+ for( unsigned long j = 0; j < N; ++j )
+ {
+ if( this->m_Collisions[ i ][ j ].second )
+ {
+ dist[ i ][ j ] = 1;
+ dist[ j ][ i ] = 1;
+ this->m_FrontPaths[ i ][ j ] = j;
+ this->m_FrontPaths[ j ][ i ] = i;
+
+ } // fi
+
+ } // rof
+ dist[ i ][ i ] = 0;
+ this->m_FrontPaths[ i ][ i ] = i;
+
+ } // rof
+
+ // Use Floyd-Warshall to compute all possible paths between fronts
+ for( unsigned long k = 0; k < N; ++k )
+ {
+ for( unsigned long i = 0; i < N; ++i )
+ {
+ for( unsigned long j = 0; j < N; ++j )
+ {
+ // WARNING: you don't want a numeric overflow!!!
+ unsigned long dik = dist[ i ][ k ];
+ unsigned long dkj = dist[ k ][ j ];
+ unsigned long sum = _inf;
+ if( dik < _inf && dkj < _inf )
+ sum = dik + dkj;
+
+ // Ok, continue Floyd-Warshall
+ if( sum < dist[ i ][ j ] )
+ {
+ dist[ i ][ j ] = sum;
+ this->m_FrontPaths[ i ][ j ] = this->m_FrontPaths[ i ][ k ];
+
+ } // fi
+
+ } // rof
+
+ } // rof
+
+ } // rof
+ this->Modified( );
+}
+
+// -------------------------------------------------------------------------
+template< class _TVertex, class _TPath, class _TSuperclass >
+void fpa::Base::MinimumSpanningTree< _TVertex, _TPath, _TSuperclass >::
+ClearSeeds( )
+{
+ this->m_Seeds.clear( );
+ this->Modified( );
+}
+
+// -------------------------------------------------------------------------
+template< class _TVertex, class _TPath, class _TSuperclass >
+void fpa::Base::MinimumSpanningTree< _TVertex, _TPath, _TSuperclass >::
+AddSeed( const _TVertex& seed )
+{
+ this->m_Seeds.push_back( seed );
+ this->Modified( );
+}
+
+// -------------------------------------------------------------------------
+template< class _TVertex, class _TPath, class _TSuperclass >
+void fpa::Base::MinimumSpanningTree< _TVertex, _TPath, _TSuperclass >::
+GetPath( typename _TPath::Pointer& path, const _TVertex& a ) const
+{
+ std::vector< _TVertex > vertices;
+ _TVertex it = a;
+ _TVertex p = this->GetParent( it );
+ while( it != p )
+ {
+ vertices.push_back( it );
+ it = p;
+ p = this->GetParent( it );
+
+ } // elihw
+ vertices.push_back( it );
+
+ if( path.IsNull( ) )
+ path = _TPath::New( );
+ for( auto v = vertices.begin( ); v != vertices.end( ); ++v )
+ path->AddVertex( *v );
+}
+
+// -------------------------------------------------------------------------
+template< class _TVertex, class _TPath, class _TSuperclass >
+void fpa::Base::MinimumSpanningTree< _TVertex, _TPath, _TSuperclass >::
+GetPath(
+ typename _TPath::Pointer& path, const _TVertex& a, const _TVertex& b
+ ) const
+{
+ static const unsigned long _inf =
+ std::numeric_limits< unsigned long >::max( );
+
+ if( path.IsNull( ) )
+ path = _TPath::New( );
+ typename _TPath::Pointer pa, pb;
+ this->GetPath( pa, a );
+ this->GetPath( pb, b );
+ if( pa->GetSize( ) > 0 && pb->GetSize( ) > 0 )
+ {
+ // Find front identifiers
+ unsigned long ia = _inf, ib = _inf;
+ unsigned long N = this->m_Seeds.size( );
+ for( unsigned long i = 0; i < N; ++i )
+ {
+ if( this->m_Seeds[ i ] == pa->GetVertex( pa->GetSize( ) - 1 ) )
+ ia = i;
+ if( this->m_Seeds[ i ] == pb->GetVertex( pb->GetSize( ) - 1 ) )
+ ib = i;
+
+ } // rof
+
+ // Check if there is a front-jump between given seeds
+ if( ia != ib )
+ {
+ // Compute front path
+ std::vector< long > fpath;
+ fpath.push_back( ia );
+ while( ia != ib )
+ {
+ ia = this->m_FrontPaths[ ia ][ ib ];
+ fpath.push_back( ia );
+
+ } // elihw
+
+ // Continue only if both fronts are connected
+ unsigned int N = fpath.size( );
+ if( N > 0 )
+ {
+ // First path: from start vertex to first collision
+ this->GetPath(
+ path, a,
+ this->m_Collisions[ fpath[ 0 ] ][ fpath[ 1 ] ].first
+ );
+
+ // Intermediary paths
+ for( unsigned int i = 1; i < N - 1; ++i )
+ {
+ typename _TPath::Pointer ipath;
+ this->GetPath(
+ ipath,
+ this->m_Collisions[ fpath[ i ] ][ fpath[ i - 1 ] ].first,
+ this->m_Collisions[ fpath[ i ] ][ fpath[ i + 1 ] ].first
+ );
+ for( long id = 0; id < ipath->GetSize( ); ++id )
+ path->AddVertex( ipath->GetVertex( id ) );
+
+ } // rof
+
+ // Final path: from last collision to end point
+ typename _TPath::Pointer lpath;
+ this->GetPath(
+ lpath,
+ this->m_Collisions[ fpath[ N - 1 ] ][ fpath[ N - 2 ] ].first, b
+ );
+ for( long id = 0; id < lpath->GetSize( ); ++id )
+ path->AddVertex( lpath->GetVertex( id ) );
+
+ } // fi
+ }
+ else
+ {
+ // Ignore common part: find common ancestor
+ long aIt = pa->GetSize( ) - 1;
+ long bIt = pb->GetSize( ) - 1;
+ bool cont = true;
+ while( aIt >= 0 && bIt >= 0 && cont )
+ {
+ cont = ( pa->GetVertex( aIt ) == pb->GetVertex( bIt ) );
+ aIt--;
+ bIt--;
+
+ } // elihw
+ aIt++;
+ bIt++;
+
+ // Glue both parts
+ for( long cIt = 0; cIt <= aIt; ++cIt )
+ path->AddVertex( pa->GetVertex( cIt ) );
+ for( ; bIt >= 0; --bIt )
+ path->AddVertex( pb->GetVertex( bIt ) );
+
+ } // fi
+
+ } // fi
+}
+
+// -------------------------------------------------------------------------
+template< class _TVertex, class _TPath, class _TSuperclass >
+fpa::Base::MinimumSpanningTree< _TVertex, _TPath, _TSuperclass >::
+MinimumSpanningTree( )
+ : Superclass( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TVertex, class _TPath, class _TSuperclass >
+fpa::Base::MinimumSpanningTree< _TVertex, _TPath, _TSuperclass >::
+~MinimumSpanningTree( )
+{
+}
+
+#endif // __fpa__Base__MinimumSpanningTree__hxx__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Base__PriorityQueueAlgorithm__h__
+#define __fpa__Base__PriorityQueueAlgorithm__h__
+
+#include <queue>
+
+namespace fpa
+{
+ namespace Base
+ {
+ /**
+ */
+ template< class _TSuperclass >
+ class PriorityQueueAlgorithm
+ : public _TSuperclass
+ {
+ public:
+ typedef PriorityQueueAlgorithm Self;
+ typedef _TSuperclass Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ protected:
+ typedef typename Superclass::_TQueueNode _TQueueNode;
+ struct _TQueueNodeCompare
+ {
+ bool operator( )( const _TQueueNode& a, const _TQueueNode& b )
+ {
+ return( b.Result < a.Result );
+ }
+ };
+ typedef std::vector< _TQueueNode > _TQueue;
+
+ public:
+ itkTypeMacro( PriorityQueueAlgorithm, Algorithm );
+
+ protected:
+ PriorityQueueAlgorithm( );
+ virtual ~PriorityQueueAlgorithm( );
+
+ virtual unsigned long _QueueSize( ) const override;
+ virtual void _QueueClear( ) override;
+ virtual void _QueuePush( const _TQueueNode& node ) override;
+ virtual _TQueueNode _QueuePop( ) override;
+
+ private:
+ // Purposely not defined
+ PriorityQueueAlgorithm( const Self& other );
+ Self& operator=( const Self& other );
+
+ protected:
+ _TQueue m_Queue;
+ };
+
+ } // ecapseman
+
+} // ecapseman
+
+#ifndef ITK_MANUAL_INSTANTIATION
+# include <fpa/Base/PriorityQueueAlgorithm.hxx>
+#endif // ITK_MANUAL_INSTANTIATION
+
+#endif // __fpa__Base__PriorityQueueAlgorithm__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Base__PriorityQueueAlgorithm__hxx__
+#define __fpa__Base__PriorityQueueAlgorithm__hxx__
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass >
+fpa::Base::PriorityQueueAlgorithm< _TSuperclass >::
+PriorityQueueAlgorithm( )
+ : Superclass( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass >
+fpa::Base::PriorityQueueAlgorithm< _TSuperclass >::
+~PriorityQueueAlgorithm( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass >
+unsigned long fpa::Base::PriorityQueueAlgorithm< _TSuperclass >::
+_QueueSize( ) const
+{
+ return( this->m_Queue.size( ) );
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass >
+void fpa::Base::PriorityQueueAlgorithm< _TSuperclass >::
+_QueueClear( )
+{
+ this->m_Queue.clear( );
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass >
+void fpa::Base::PriorityQueueAlgorithm< _TSuperclass >::
+_QueuePush( const _TQueueNode& node )
+{
+ static _TQueueNodeCompare cmp;
+ this->m_Queue.push_back( node );
+ std::push_heap( this->m_Queue.begin( ), this->m_Queue.end( ), cmp );
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass >
+typename fpa::Base::PriorityQueueAlgorithm< _TSuperclass >::
+_TQueueNode fpa::Base::PriorityQueueAlgorithm< _TSuperclass >::
+_QueuePop( )
+{
+ static _TQueueNodeCompare cmp;
+ std::pop_heap( this->m_Queue.begin( ), this->m_Queue.end( ), cmp );
+ _TQueueNode f = this->m_Queue.back( );
+ this->m_Queue.pop_back( );
+ return( f );
+}
+
+#endif // __fpa__Base__PriorityQueueAlgorithm__hxx__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Base__QueueAlgorithm__h__
+#define __fpa__Base__QueueAlgorithm__h__
+
+#include <queue>
+#include <itkMacro.h>
+#include <itkSmartPointer.h>
+
+namespace fpa
+{
+ namespace Base
+ {
+ /**
+ */
+ template< class _TSuperclass >
+ class QueueAlgorithm
+ : public _TSuperclass
+ {
+ public:
+ typedef QueueAlgorithm Self;
+ typedef _TSuperclass Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ protected:
+ typedef typename Superclass::_TQueueNode _TQueueNode;
+ typedef std::queue< _TQueueNode > _TQueue;
+
+ public:
+ itkTypeMacro( QueueAlgorithm, _TSuperclass );
+
+ protected:
+ QueueAlgorithm( );
+ virtual ~QueueAlgorithm( );
+
+ virtual unsigned long _QueueSize( ) const override;
+ virtual void _QueueClear( ) override;
+ virtual void _QueuePush( const _TQueueNode& node ) override;
+ virtual _TQueueNode _QueuePop( ) override;
+
+ private:
+ // Purposely not defined
+ QueueAlgorithm( const Self& other );
+ Self& operator=( const Self& other );
+
+ protected:
+ _TQueue m_Queue;
+ };
+
+ } // ecapseman
+
+} // ecapseman
+
+#ifndef ITK_MANUAL_INSTANTIATION
+# include <fpa/Base/QueueAlgorithm.hxx>
+#endif // ITK_MANUAL_INSTANTIATION
+
+#endif // __fpa__Base__QueueAlgorithm__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Base__QueueAlgorithm__hxx__
+#define __fpa__Base__QueueAlgorithm__hxx__
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass >
+fpa::Base::QueueAlgorithm< _TSuperclass >::
+QueueAlgorithm( )
+ : Superclass( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass >
+fpa::Base::QueueAlgorithm< _TSuperclass >::
+~QueueAlgorithm( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass >
+unsigned long fpa::Base::QueueAlgorithm< _TSuperclass >::
+_QueueSize( ) const
+{
+ return( this->m_Queue.size( ) );
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass >
+void fpa::Base::QueueAlgorithm< _TSuperclass >::
+_QueueClear( )
+{
+ while( this->m_Queue.size( ) > 0 )
+ this->m_Queue.pop( );
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass >
+void fpa::Base::QueueAlgorithm< _TSuperclass >::
+_QueuePush( const _TQueueNode& node )
+{
+ this->m_Queue.push( node );
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass >
+typename fpa::Base::QueueAlgorithm< _TSuperclass >::
+_TQueueNode fpa::Base::QueueAlgorithm< _TSuperclass >::
+_QueuePop( )
+{
+ _TQueueNode f = this->m_Queue.front( );
+ this->m_Queue.pop( );
+ return( f );
+}
+
+#endif // __fpa__Base__QueueAlgorithm__hxx__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Base__RegionGrow__h__
+#define __fpa__Base__RegionGrow__h__
+
+#include <queue>
+#include <fpa/Base/QueueAlgorithm.h>
+#include <fpa/Base/Functors/RegionGrow/Base.h>
+
+namespace fpa
+{
+ namespace Base
+ {
+ /**
+ */
+ template< class _TSuperclass >
+ class RegionGrow
+ : public fpa::Base::QueueAlgorithm< _TSuperclass >
+ {
+ public:
+ typedef RegionGrow Self;
+ typedef fpa::Base::QueueAlgorithm< _TSuperclass > Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ typedef typename Superclass::TOutput TOutput;
+ typedef typename Superclass::TVertex TVertex;
+ typedef fpa::Base::Functors::RegionGrow::Base< TVertex, TOutput > TGrowFunction;
+
+ protected:
+ typedef typename Superclass::_TQueueNode _TQueueNode;
+
+ public:
+ itkTypeMacro( RegionGrow, Algorithm );
+
+ public:
+ TGrowFunction* GetGrowFunction( );
+ const TGrowFunction* GetGrowFunction( ) const;
+ TOutput GetInsideValue( ) const;
+ TOutput GetOutsideValue( ) const;
+
+ void SetGrowFunction( TGrowFunction* f );
+ void SetInsideValue( const TOutput& v );
+ void SetOutsideValue( const TOutput& v );
+
+ protected:
+ RegionGrow( );
+ virtual ~RegionGrow( );
+
+ virtual bool _UpdateValue(
+ _TQueueNode& v, const _TQueueNode& p
+ ) override;
+ virtual TOutput _GetInputValue( const _TQueueNode& v, const _TQueueNode& p ) override
+ {
+ TOutput res = this->m_InitResult;
+ if( this->m_GrowFunction.IsNotNull( ) )
+ res = this->m_GrowFunction->Evaluate( v.Vertex, p.Vertex );
+ return( res );
+ }
+
+
+ private:
+ // Purposely not defined
+ RegionGrow( const Self& other );
+ Self& operator=( const Self& other );
+
+ protected:
+ typename TGrowFunction::Pointer m_GrowFunction;
+ };
+
+ } // ecapseman
+
+} // ecapseman
+
+#ifndef ITK_MANUAL_INSTANTIATION
+# include <fpa/Base/RegionGrow.hxx>
+#endif // ITK_MANUAL_INSTANTIATION
+
+#endif // __fpa__Base__RegionGrow__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Base__RegionGrow__hxx__
+#define __fpa__Base__RegionGrow__hxx__
+
+#include <fpa/Base/Functors/RegionGrow/Tautology.h>
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass >
+typename fpa::Base::RegionGrow< _TSuperclass >::
+TGrowFunction* fpa::Base::RegionGrow< _TSuperclass >::
+GetGrowFunction( )
+{
+ // TODO: return( dynamic_cast< TGrowFunction* >( this->GetVertexFunction( ) ) );
+ return( this->m_GrowFunction );
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass >
+const typename fpa::Base::RegionGrow< _TSuperclass >::
+TGrowFunction* fpa::Base::RegionGrow< _TSuperclass >::
+GetGrowFunction( ) const
+{
+ /* TODO
+ return(
+ dynamic_cast< const TGrowFunction* >( this->GetVertexFunction( ) )
+ );
+ */
+ return( this->m_GrowFunction );
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass >
+typename fpa::Base::RegionGrow< _TSuperclass >::
+TOutput fpa::Base::RegionGrow< _TSuperclass >::
+GetInsideValue( ) const
+{
+ const TGrowFunction* f = this->GetGrowFunction( );
+ if( f != NULL )
+ return( f->GetInsideValue( ) );
+ else
+ return( this->m_InitResult );
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass >
+typename fpa::Base::RegionGrow< _TSuperclass >::
+TOutput fpa::Base::RegionGrow< _TSuperclass >::
+GetOutsideValue( ) const
+{
+ const TGrowFunction* f = this->GetGrowFunction( );
+ if( f != NULL )
+ return( f->GetOutsideValue( ) );
+ else
+ return( this->m_InitResult );
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass >
+void fpa::Base::RegionGrow< _TSuperclass >::
+SetGrowFunction( TGrowFunction* f )
+{
+ TGrowFunction* old_f = this->GetGrowFunction( );
+ if( old_f != NULL )
+ {
+ f->SetInsideValue( old_f->GetInsideValue( ) );
+ f->SetOutsideValue( old_f->GetOutsideValue( ) );
+
+ } // fi
+ this->m_GrowFunction = f;
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass >
+void fpa::Base::RegionGrow< _TSuperclass >::
+SetInsideValue( const TOutput& v )
+{
+ TGrowFunction* f = this->GetGrowFunction( );
+ if( f != NULL )
+ {
+ f->SetInsideValue( v );
+ this->Modified( );
+
+ } // fi
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass >
+void fpa::Base::RegionGrow< _TSuperclass >::
+SetOutsideValue( const TOutput& v )
+{
+ TGrowFunction* f = this->GetGrowFunction( );
+ if( f != NULL )
+ {
+ f->SetOutsideValue( v );
+ this->Modified( );
+
+ } // fi
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass >
+fpa::Base::RegionGrow< _TSuperclass >::
+RegionGrow( )
+ : Superclass( )
+{
+ typedef fpa::Base::Functors::RegionGrow::Tautology< TVertex, TOutput > _TFunc;
+ this->SetGrowFunction( _TFunc::New( ) );
+ this->m_InitResult = this->GetGrowFunction( )->GetOutsideValue( );
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass >
+fpa::Base::RegionGrow< _TSuperclass >::
+~RegionGrow( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TSuperclass >
+bool fpa::Base::RegionGrow< _TSuperclass >::
+_UpdateValue( _TQueueNode& v, const _TQueueNode& p )
+{
+ v.Result = this->_GetInputValue( v, p );
+ return( v.Result == this->GetInsideValue( ) );
+}
+
+#endif // __fpa__Base__RegionGrow__hxx__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Config__h__
+#define __fpa__Config__h__
+
+/*
+ * =========================================================================
+ * Language related macros
+ * =========================================================================
+ */
+
+#endif // __fpa__Config__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Image__Algorithm__h__
+#define __fpa__Image__Algorithm__h__
+
+#include <itkFunctionBase.h>
+#include <itkImageToImageFilter.h>
+#include <fpa/Base/Algorithm.h>
+#include <fpa/Image/Functors/Base.h>
+
+namespace fpa
+{
+ namespace Image
+ {
+ /**
+ */
+ template< class _TInputImage, class _TOutputImage >
+ class Algorithm
+ : public fpa::Base::Algorithm< itk::ImageToImageFilter< _TInputImage, _TOutputImage >, typename _TInputImage::IndexType, typename _TOutputImage::PixelType >
+ {
+ public:
+ typedef itk::ImageToImageFilter< _TInputImage, _TOutputImage > TFilter;
+ typedef typename _TInputImage::IndexType TVertex;
+ typedef typename _TOutputImage::PixelType TOutput;
+
+ typedef Algorithm Self;
+ typedef fpa::Base::Algorithm< TFilter, TVertex, TOutput > Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ typedef typename Superclass::TFrontId TFrontId;
+ typedef typename Superclass::TNeighborhood TNeighborhood;
+
+ typedef fpa::Image::Functors::Base< itk::ImageBase< _TInputImage::ImageDimension >, typename Superclass::TNeighborhoodFunction > TNeighborhoodFunction;
+ typedef fpa::Image::Functors::Base< _TInputImage, typename Superclass::TVertexFunction > TVertexFunction;
+
+ protected:
+ typedef typename Superclass::_TQueueNode _TQueueNode;
+
+ public:
+ itkTypeMacro( fpa::Image::Algorithm, fpa::Base::Algorithm );
+
+ protected:
+ Algorithm( );
+ virtual ~Algorithm( );
+
+ virtual void _BeforeGenerateData( ) override;
+ virtual void _InitMarks( ) override;
+ virtual void _InitResults( const TOutput& init_value ) override;
+ virtual bool _IsMarked( const TVertex& v ) const override;
+ virtual void _Mark( const _TQueueNode& n ) override;
+ virtual TFrontId _GetMark( const TVertex& v ) const override;
+ virtual void _UpdateResult( const _TQueueNode& n ) override;
+ virtual TOutput _GetResult( const TVertex& v ) const override;
+ virtual unsigned int _GetNumberOfDimensions( ) const override;
+
+ private:
+ // Purposely not defined
+ Algorithm( const Self& other );
+ Self& operator=( const Self& other );
+
+ protected:
+ unsigned int m_MarksIdx;
+ };
+
+ } // ecapseman
+
+} // ecapseman
+
+#ifndef ITK_MANUAL_INSTANTIATION
+# include <fpa/Image/Algorithm.hxx>
+#endif // ITK_MANUAL_INSTANTIATION
+
+#endif // __fpa__Image__Algorithm__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Image__Algorithm__hxx__
+#define __fpa__Image__Algorithm__hxx__
+
+// Send Piotr's code to Anna
+
+#include <itkImage.h>
+#include <fpa/Image/Functors/SimpleNeighborhood.h>
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+fpa::Image::Algorithm< _TInputImage, _TOutputImage >::
+Algorithm( )
+ : Superclass( )
+{
+ typedef fpa::Image::Functors::SimpleNeighborhood< _TInputImage::ImageDimension > _TNeigh;
+ typedef itk::Image< TFrontId, _TInputImage::ImageDimension > _TMarks;
+
+ this->m_MarksIdx = this->GetNumberOfRequiredOutputs( );
+ this->itk::ProcessObject::SetNumberOfRequiredOutputs( this->m_MarksIdx + 1 );
+ typename _TMarks::Pointer marks = _TMarks::New( );
+ this->SetNthOutput( this->m_MarksIdx, marks );
+
+ typename _TNeigh::Pointer neigh = _TNeigh::New( );
+ this->SetNeighborhoodFunction( neigh );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+fpa::Image::Algorithm< _TInputImage, _TOutputImage >::
+~Algorithm( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+void fpa::Image::Algorithm< _TInputImage, _TOutputImage >::
+_BeforeGenerateData( )
+{
+ this->Superclass::_BeforeGenerateData( );
+ this->AllocateOutputs( );
+
+ TNeighborhoodFunction* neighFunc =
+ dynamic_cast< TNeighborhoodFunction* >(
+ this->GetNeighborhoodFunction( )
+ );
+ if( neighFunc == NULL )
+ itkExceptionMacro( << "NeighborhoodFunction not well defined." );
+ neighFunc->SetImage( this->GetInput( ) );
+
+ TVertexFunction* vertexFunc =
+ dynamic_cast< TVertexFunction* >(
+ this->GetVertexFunction( )
+ );
+ if( vertexFunc != NULL )
+ vertexFunc->SetImage( this->GetInput( ) );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+void fpa::Image::Algorithm< _TInputImage, _TOutputImage >::
+_InitMarks( )
+{
+ typedef itk::Image< TFrontId, _TInputImage::ImageDimension > _TMarks;
+ _TMarks* marks =
+ dynamic_cast< _TMarks* >(
+ this->itk::ProcessObject::GetOutput( this->m_MarksIdx )
+ );
+ marks->FillBuffer( 0 );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+void fpa::Image::Algorithm< _TInputImage, _TOutputImage >::
+_InitResults( const TOutput& init_value )
+{
+ this->GetOutput( )->FillBuffer( init_value );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+bool fpa::Image::Algorithm< _TInputImage, _TOutputImage >::
+_IsMarked( const TVertex& v ) const
+{
+ typedef itk::Image< TFrontId, _TInputImage::ImageDimension > _TMarks;
+ const _TMarks* marks =
+ dynamic_cast< const _TMarks* >(
+ this->itk::ProcessObject::GetOutput( this->m_MarksIdx )
+ );
+ return( marks->GetPixel( v ) != 0 );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+void fpa::Image::Algorithm< _TInputImage, _TOutputImage >::
+_Mark( const _TQueueNode& n )
+{
+ typedef itk::Image< TFrontId, _TInputImage::ImageDimension > _TMarks;
+ _TMarks* marks =
+ dynamic_cast< _TMarks* >(
+ this->itk::ProcessObject::GetOutput( this->m_MarksIdx )
+ );
+ marks->SetPixel( n.Vertex, n.FrontId );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+typename fpa::Image::Algorithm< _TInputImage, _TOutputImage >::
+TFrontId fpa::Image::Algorithm< _TInputImage, _TOutputImage >::
+_GetMark( const TVertex& v ) const
+{
+ typedef itk::Image< TFrontId, _TInputImage::ImageDimension > _TMarks;
+ const _TMarks* marks =
+ dynamic_cast< const _TMarks* >(
+ this->itk::ProcessObject::GetOutput( this->m_MarksIdx )
+ );
+ return( marks->GetPixel( v ) );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+void fpa::Image::Algorithm< _TInputImage, _TOutputImage >::
+_UpdateResult( const _TQueueNode& n )
+{
+ this->GetOutput( )->SetPixel( n.Vertex, n.Result );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+typename fpa::Image::Algorithm< _TInputImage, _TOutputImage >::TOutput
+fpa::Image::Algorithm< _TInputImage, _TOutputImage >::
+_GetResult( const TVertex& v ) const
+{
+ if( this->GetOutput( )->GetLargestPossibleRegion( ).IsInside( v ) )
+ return( this->GetOutput( )->GetPixel( v ) );
+ else
+ return( this->m_InitResult );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+unsigned int fpa::Image::Algorithm< _TInputImage, _TOutputImage >::
+_GetNumberOfDimensions( ) const
+{
+ return( _TInputImage::ImageDimension );
+}
+
+#endif // __fpa__Image__Algorithm__hxx__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Image__Dijkstra__h__
+#define __fpa__Image__Dijkstra__h__
+
+#include <fpa/Base/Dijkstra.h>
+#include <fpa/Image/Algorithm.h>
+#include <fpa/Image/MinimumSpanningTree.h>
+#include <fpa/Image/Functors/Base.h>
+
+namespace fpa
+{
+ namespace Image
+ {
+ /**
+ */
+ template< class _TInputImage, class _TOutputImage >
+ class Dijkstra
+ : public fpa::Base::Dijkstra< fpa::Image::Algorithm< _TInputImage, _TOutputImage >, fpa::Image::MinimumSpanningTree< _TInputImage::ImageDimension > >
+ {
+ public:
+ typedef Dijkstra Self;
+ typedef fpa::Image::Algorithm< _TInputImage, _TOutputImage > TAlgorithm;
+ typedef fpa::Image::MinimumSpanningTree< _TInputImage::ImageDimension > TMST;
+ typedef fpa::Base::Dijkstra< TAlgorithm, TMST > Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ typedef typename Superclass::TOutput TOutput;
+ typedef typename Superclass::TVertex TVertex;
+
+ public:
+ itkNewMacro( Self );
+ itkTypeMacro( fpa::Image::Dijkstra, fpa::Base::Dijkstra );
+
+ protected:
+ Dijkstra( );
+ virtual ~Dijkstra( );
+
+ private:
+ // Purposely not defined
+ Dijkstra( const Self& other );
+ Self& operator=( const Self& other );
+ };
+
+ } // ecapseman
+
+} // ecapseman
+
+#ifndef ITK_MANUAL_INSTANTIATION
+# include <fpa/Image/Dijkstra.hxx>
+#endif // ITK_MANUAL_INSTANTIATION
+
+#endif // __fpa__Image__Dijkstra__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Image__Dijkstra__hxx__
+#define __fpa__Image__Dijkstra__hxx__
+
+#include <fpa/Image/Functors/VertexCost.h>
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+fpa::Image::Dijkstra< _TInputImage, _TOutputImage >::
+Dijkstra( )
+ : Superclass( )
+{
+ typedef fpa::Image::Functors::VertexCost< _TInputImage, TOutput > _TCost;
+ typename _TCost::Pointer cost = _TCost::New( );
+ this->SetVertexFunction( cost );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+fpa::Image::Dijkstra< _TInputImage, _TOutputImage >::
+~Dijkstra( )
+{
+}
+
+#endif // __fpa__Image__Dijkstra__hxx__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Image__EndPointsFilter__h__
+#define __fpa__Image__EndPointsFilter__h__
+
+#include <fpa/Image/MinimumSpanningTree.h>
+#include <set>
+
+namespace fpa
+{
+ namespace Image
+ {
+ /**
+ */
+ template< class _TDistanceMap, class _TCostMap >
+ class EndPointsFilter
+ : public itk::Object
+ {
+ public:
+ typedef EndPointsFilter Self;
+ typedef itk::Object Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ typedef _TDistanceMap TDistanceMap;
+ typedef _TCostMap TCostMap;
+ typedef typename TCostMap::IndexType TIndex;
+ typedef MinimumSpanningTree< TCostMap::ImageDimension > TMST;
+
+ typedef itk::Functor::IndexLexicographicCompare< _TCostMap::ImageDimension > TIndexCompare;
+ typedef std::set< TIndex, TIndexCompare > TIndices;
+
+ public:
+ itkNewMacro( Self );
+ itkTypeMacro( fpa::Image::EndPointsFilter, itk::Object );
+
+ itkGetConstObjectMacro( DistanceMap, _TDistanceMap );
+ itkGetConstObjectMacro( CostMap, _TCostMap );
+ itkGetConstObjectMacro( MST, TMST );
+
+ itkSetConstObjectMacro( DistanceMap, _TDistanceMap );
+ itkSetConstObjectMacro( CostMap, _TCostMap );
+ itkSetConstObjectMacro( MST, TMST );
+
+ public:
+ const TIndices& GetBifurcations( ) const;
+ const TIndices& GetEndPoints( ) const;
+
+ void Compute( );
+
+ protected:
+ EndPointsFilter( );
+ virtual ~EndPointsFilter( );
+
+ private:
+ // Purposely not defined
+ EndPointsFilter( const Self& other );
+ Self& operator=( const Self& other );
+
+ protected:
+ typename _TDistanceMap::ConstPointer m_DistanceMap;
+ typename _TCostMap::ConstPointer m_CostMap;
+ typename TMST::ConstPointer m_MST;
+
+ TIndices m_Bifurcations;
+ TIndices m_EndPoints;
+ };
+
+ } // ecapseman
+
+} // ecapseman
+
+#ifndef ITK_MANUAL_INSTANTIATION
+# include <fpa/Image/EndPointsFilter.hxx>
+#endif // ITK_MANUAL_INSTANTIATION
+
+#endif // __fpa__Image__EndPointsFilter__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Image__EndPointsFilter__hxx__
+#define __fpa__Image__EndPointsFilter__hxx__
+
+#include <functional>
+#include <map>
+#include <queue>
+#include <itkImageRegionConstIteratorWithIndex.h>
+
+// -------------------------------------------------------------------------
+template< class _TDistanceMap, class _TCostMap >
+const typename fpa::Image::EndPointsFilter< _TDistanceMap, _TCostMap >::
+TIndices& fpa::Image::EndPointsFilter< _TDistanceMap, _TCostMap >::
+GetBifurcations( ) const
+{
+ return( this->m_Bifurcations );
+}
+
+// -------------------------------------------------------------------------
+template< class _TDistanceMap, class _TCostMap >
+const typename fpa::Image::EndPointsFilter< _TDistanceMap, _TCostMap >::
+TIndices& fpa::Image::EndPointsFilter< _TDistanceMap, _TCostMap >::
+GetEndPoints( ) const
+{
+ return( this->m_EndPoints );
+}
+
+// -------------------------------------------------------------------------
+template< class _TDistanceMap, class _TCostMap >
+void fpa::Image::EndPointsFilter< _TDistanceMap, _TCostMap >::
+Compute( )
+{
+ typedef itk::ImageRegionConstIteratorWithIndex< _TDistanceMap > _TDistMapIt;
+ typedef itk::ImageRegionConstIteratorWithIndex< _TCostMap > _TCostMapIt;
+ typedef std::multimap< double, TIndex, std::greater< double > > _TQueue;
+ typedef typename _TQueue::value_type _TQueueValue;
+ typedef itk::Image< unsigned char, _TCostMap::ImageDimension > _TMarks;
+
+ // Some values
+ typename _TDistanceMap::RegionType region =
+ this->m_DistanceMap->GetRequestedRegion( );
+ typename _TMarks::Pointer marks = _TMarks::New( );
+ marks->SetLargestPossibleRegion( this->m_DistanceMap->GetLargestPossibleRegion( ) );
+ marks->SetRequestedRegion( this->m_DistanceMap->GetRequestedRegion( ) );
+ marks->SetBufferedRegion( this->m_DistanceMap->GetBufferedRegion( ) );
+ marks->SetSpacing( this->m_DistanceMap->GetSpacing( ) );
+ marks->SetOrigin( this->m_DistanceMap->GetOrigin( ) );
+ marks->SetDirection( this->m_DistanceMap->GetDirection( ) );
+ marks->Allocate( );
+ marks->FillBuffer( 0 );
+
+ // Create queue
+ _TQueue queue;
+ _TDistMapIt dIt( this->m_DistanceMap, region );
+ _TCostMapIt cIt( this->m_CostMap, region );
+ dIt.GoToBegin( );
+ cIt.GoToBegin( );
+ for( ; !dIt.IsAtEnd( ) && !cIt.IsAtEnd( ); ++dIt, ++cIt )
+ {
+ double d = double( dIt.Get( ) );
+ if( d > 0 )
+ {
+ double v = double( cIt.Get( ) ) / ( d * d );
+ queue.insert( _TQueueValue( v, dIt.GetIndex( ) ) );
+
+ } // fi
+
+ } // rof
+
+ // BFS from maximum queue
+ while( queue.size( ) > 0 )
+ {
+ // Get node
+ auto nIt = queue.begin( );
+ auto n = *nIt;
+ queue.erase( nIt );
+
+ unsigned char m = marks->GetPixel( n.second );
+ if( ( m & 0x01 ) == 0x01 )
+ continue;
+
+ // Mark it
+ m |= 0x01;
+ marks->SetPixel( n.second, m );
+ this->m_EndPoints.insert( n.second );
+
+ // Get path
+ typename TMST::TPath::Pointer path;
+ this->m_MST->GetPath( path, n.second );
+ for( unsigned long i = 0; i < path->GetSize( ); ++i )
+ {
+ TIndex idx = path->GetVertex( path->GetSize( ) - 1 - i );
+ typename _TCostMap::PointType cnt;
+ this->m_CostMap->TransformIndexToPhysicalPoint( idx, cnt );
+ double d = double( this->m_DistanceMap->GetPixel( idx ) );
+ d *= double( 2 );
+
+ // Mark sphere
+ std::queue< TIndex > q;
+ q.push( idx );
+ while( q.size( ) > 0 )
+ {
+ TIndex v = q.front( );
+ q.pop( );
+ unsigned char m = marks->GetPixel( v );
+ if( ( m & 0x02 ) == 0x02 )
+ continue;
+ m |= 0x03;
+ marks->SetPixel( v, m );
+
+ for( unsigned int x = 0; x < _TCostMap::ImageDimension; ++x )
+ {
+ for( int y = -1; y <= 1; y += 2 )
+ {
+ TIndex w = v;
+ w[ x ] += y;
+ if( region.IsInside( w ) )
+ {
+ typename _TCostMap::PointType p;
+ this->m_CostMap->TransformIndexToPhysicalPoint( w, p );
+ if( cnt.EuclideanDistanceTo( p ) <= d )
+ q.push( w );
+
+ } // fi
+
+ } // rof
+
+ } // rof
+
+ } // elihw
+
+ } // rof
+
+ } // elihw
+}
+
+// -------------------------------------------------------------------------
+template< class _TDistanceMap, class _TCostMap >
+fpa::Image::EndPointsFilter< _TDistanceMap, _TCostMap >::
+EndPointsFilter( )
+ : Superclass( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TDistanceMap, class _TCostMap >
+fpa::Image::EndPointsFilter< _TDistanceMap, _TCostMap >::
+~EndPointsFilter( )
+{
+}
+
+#endif // __fpa__Image__EndPointsFilter__hxx__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Image__FastMarching__h__
+#define __fpa__Image__FastMarching__h__
+
+#include <fpa/Base/FastMarching.h>
+#include <fpa/Image/Algorithm.h>
+
+namespace fpa
+{
+ namespace Image
+ {
+ /**
+ */
+ template< class _TInputImage, class _TOutputImage >
+ class FastMarching
+ : public fpa::Base::FastMarching< fpa::Image::Algorithm< _TInputImage, _TOutputImage > >
+ {
+ public:
+ typedef FastMarching Self;
+ typedef fpa::Image::Algorithm< _TInputImage, _TOutputImage > TAlgorithm;
+ typedef fpa::Base::FastMarching< TAlgorithm > Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ typedef typename Superclass::TOutput TOutput;
+ typedef typename Superclass::TVertex TVertex;
+ typedef typename Superclass::TFastMarchingNeighbor TFastMarchingNeighbor;
+ typedef typename Superclass::TFastMarchingNeighborhood TFastMarchingNeighborhood;
+
+ public:
+ itkNewMacro( Self );
+ itkTypeMacro( fpa::Image::FastMarching, fpa::Base::FastMarching );
+
+ protected:
+ FastMarching( );
+ virtual ~FastMarching( );
+
+ virtual TFastMarchingNeighborhood _FastMarchingNeighbors( const TVertex& v ) const override;
+
+ private:
+ // Purposely not defined
+ FastMarching( const Self& other );
+ Self& operator=( const Self& other );
+ };
+
+ } // ecapseman
+
+} // ecapseman
+
+#ifndef ITK_MANUAL_INSTANTIATION
+# include <fpa/Image/FastMarching.hxx>
+#endif // ITK_MANUAL_INSTANTIATION
+
+#endif // __fpa__Image__FastMarching__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Image__FastMarching__hxx__
+#define __fpa__Image__FastMarching__hxx__
+
+#include <fpa/Image/Functors/VertexCost.h>
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+fpa::Image::FastMarching< _TInputImage, _TOutputImage >::
+FastMarching( )
+ : Superclass( )
+{
+ typedef fpa::Image::Functors::VertexCost< _TInputImage, TOutput > _TCost;
+ typename _TCost::Pointer cost = _TCost::New( );
+ this->SetVertexFunction( cost );
+ this->m_InitResult = std::numeric_limits< TOutput >::max( );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+fpa::Image::FastMarching< _TInputImage, _TOutputImage >::
+~FastMarching( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+typename fpa::Image::FastMarching< _TInputImage, _TOutputImage >::
+TFastMarchingNeighborhood
+fpa::Image::FastMarching< _TInputImage, _TOutputImage >::
+_FastMarchingNeighbors( const TVertex& v ) const
+{
+ TFastMarchingNeighborhood neighs;
+ typename _TInputImage::SpacingType spac = this->GetInput( )->GetSpacing( );
+ for( unsigned int d = 0; d < _TInputImage::ImageDimension; ++d )
+ {
+ for( int i = -1; i <= 1; i += 2 )
+ {
+ TVertex n = v;
+ n[ d ] += i;
+ neighs.push_back( TFastMarchingNeighbor( n, spac[ d ] ) );
+
+ } // rof
+
+ } // rof
+ return( neighs );
+}
+
+#endif // __fpa__Image__FastMarching__hxx__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Image__Functors__Base__h__
+#define __fpa__Image__Functors__Base__h__
+
+#include <fpa/Config.h>
+#include <itkFunctionBase.h>
+#include <itkImageBase.h>
+
+namespace fpa
+{
+ namespace Image
+ {
+ namespace Functors
+ {
+ /**
+ */
+ template< class _TImage, class _TSuperclass >
+ class Base
+ : public _TSuperclass
+ {
+ public:
+ typedef Base Self;
+ typedef _TSuperclass Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ typedef _TImage TImage;
+ typedef itk::ImageBase< TImage::ImageDimension > TImageBase;
+
+ public:
+ itkTypeMacro( Base, itk::FunctionBase );
+
+ itkGetConstObjectMacro( Image, TImage );
+ itkSetConstObjectMacro( Image, TImage );
+
+ protected:
+ Base( ) : Superclass( ) { }
+ virtual ~Base( ) { }
+
+ private:
+ // Purposely not implemented
+ Base( const Self& other );
+ Self& operator=( const Self& other );
+
+ protected:
+ typename TImage::ConstPointer m_Image;
+ };
+
+ } // ecapseman
+
+ } // ecapseman
+
+} // ecapseman
+
+#endif // __fpa__Image__Functors__Base__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Image__Functors__RegionGrow__BinaryThreshold__h__
+#define __fpa__Image__Functors__RegionGrow__BinaryThreshold__h__
+
+#include <fpa/Image/Functors/Base.h>
+#include <fpa/Base/Functors/RegionGrow/Base.h>
+
+namespace fpa
+{
+ namespace Image
+ {
+ namespace Functors
+ {
+ namespace RegionGrow
+ {
+ /**
+ */
+ template< class _TImage, class _TOutput >
+ class BinaryThreshold
+ : public fpa::Image::Functors::Base< _TImage, fpa::Base::Functors::RegionGrow::Base< typename _TImage::IndexType, _TOutput > >
+ {
+ public:
+ typedef _TImage TImage;
+ typedef _TOutput TOutput;
+ typedef typename TImage::IndexType TIndex;
+ typedef typename TImage::PixelType TPixel;
+
+ typedef fpa::Base::Functors::RegionGrow::Base< TIndex, TOutput > TBase;
+ typedef fpa::Image::Functors::Base< TImage, TBase > Superclass;
+ typedef BinaryThreshold Self;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ public:
+ itkNewMacro( Self );
+ itkTypeMacro( BinaryThreshold, Base );
+
+ itkGetConstMacro( Lower, TPixel );
+ itkGetConstMacro( Upper, TPixel );
+ itkSetMacro( Lower, TPixel );
+ itkSetMacro( Upper, TPixel );
+
+ public:
+ virtual TOutput Evaluate(
+ const TIndex& a, const TIndex& b
+ ) const override;
+
+ protected:
+ BinaryThreshold( );
+ virtual ~BinaryThreshold( );
+
+ private:
+ // Purposely not implemented
+ BinaryThreshold( const Self& other );
+ Self& operator=( const Self& other );
+
+ protected:
+ TPixel m_Lower;
+ TPixel m_Upper;
+ };
+
+ } // ecapseman
+
+ } // ecapseman
+
+ } // ecapseman
+
+} // ecapseman
+
+#ifndef ITK_MANUAL_INSTANTIATION
+# include <fpa/Image/Functors/RegionGrow/BinaryThreshold.hxx>
+#endif // ITK_MANUAL_INSTANTIATION
+
+#endif // __fpa__Image__Functors__RegionGrow__BinaryThreshold__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Image__Functors__RegionGrow__BinaryThreshold__hxx__
+#define __fpa__Image__Functors__RegionGrow__BinaryThreshold__hxx__
+
+#include <limits>
+
+// -------------------------------------------------------------------------
+template< class _TImage, class _TOutput >
+_TOutput
+fpa::Image::Functors::RegionGrow::BinaryThreshold< _TImage, _TOutput >::
+Evaluate( const TIndex& a, const TIndex& b ) const
+{
+ const _TImage* im =
+ dynamic_cast< const _TImage* >( this->m_Image.GetPointer( ) );
+ if( im != NULL )
+ {
+ TPixel v = im->GetPixel( b );
+ return(
+ ( this->m_Lower < v && v < this->m_Upper )?
+ this->m_InsideValue:
+ this->m_OutsideValue
+ );
+ }
+ else
+ return( this->m_OutsideValue );
+}
+
+// -------------------------------------------------------------------------
+template< class _TImage, class _TOutput >
+fpa::Image::Functors::RegionGrow::BinaryThreshold< _TImage, _TOutput >::
+BinaryThreshold( )
+ : Superclass( )
+{
+ this->m_Upper = std::numeric_limits< TPixel >::max( );
+ if( std::numeric_limits< TPixel >::is_integer )
+ this->m_Lower = std::numeric_limits< TPixel >::min( );
+ else
+ this->m_Lower = -this->m_Upper;
+}
+
+// -------------------------------------------------------------------------
+template< class _TImage, class _TOutput >
+fpa::Image::Functors::RegionGrow::BinaryThreshold< _TImage, _TOutput >::
+~BinaryThreshold( )
+{
+}
+
+#endif // __fpa__Image__Functors__RegionGrow__BinaryThreshold__hxx__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Image__Functors__SimpleNeighborhood__h__
+#define __fpa__Image__Functors__SimpleNeighborhood__h__
+
+#include <vector>
+#include <fpa/Image/Functors/Base.h>
+#include <itkFunctionBase.h>
+#include <itkImageBase.h>
+
+namespace fpa
+{
+ namespace Image
+ {
+ namespace Functors
+ {
+ /**
+ */
+ template< unsigned int _VDim >
+ class SimpleNeighborhood
+ : public fpa::Image::Functors::Base< itk::ImageBase< _VDim >, itk::FunctionBase< itk::Index< _VDim >, std::vector< itk::Index< _VDim > > > >
+ {
+ public:
+ typedef itk::ImageBase< _VDim > TImage;
+ typedef typename TImage::IndexType TIndex;
+ typedef typename TIndex::OffsetType TOffset;
+ typedef std::vector< TIndex > TOutput;
+ typedef itk::FunctionBase< TIndex, TOutput > TBaseFunctor;
+ typedef fpa::Image::Functors::Base< TImage, TBaseFunctor > Superclass;
+ typedef SimpleNeighborhood Self;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ public:
+ itkNewMacro( Self );
+ itkTypeMacro( SimpleNeighborhood, Base );
+
+ itkGetConstMacro( Order, unsigned int );
+ itkSetMacro( Order, unsigned int );
+
+ public:
+ virtual TOutput Evaluate( const TIndex& center ) const override;
+
+ protected:
+ SimpleNeighborhood( );
+ virtual ~SimpleNeighborhood( );
+
+ void _1stCombination( ) const;
+ void _2ndCombination( ) const;
+
+ private:
+ // Purposely not implemented
+ SimpleNeighborhood( const Self& other );
+ Self& operator=( const Self& other );
+
+ protected:
+ unsigned int m_Order;
+ mutable std::vector< TOffset > m_Offsets;
+ };
+
+ } // ecapseman
+
+ } // ecapseman
+
+} // ecapseman
+
+#ifndef ITK_MANUAL_INSTANTIATION
+# include <fpa/Image/Functors/SimpleNeighborhood.hxx>
+#endif // ITK_MANUAL_INSTANTIATION
+
+#endif // __fpa__Image__Functors__SimpleNeighborhood__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Image__Functors__SimpleNeighborhood__hxx__
+#define __fpa__Image__Functors__SimpleNeighborhood__hxx__
+
+// -------------------------------------------------------------------------
+template< unsigned int _VDim >
+typename fpa::Image::Functors::SimpleNeighborhood< _VDim >::
+TOutput fpa::Image::Functors::SimpleNeighborhood< _VDim >::
+Evaluate( const TIndex& center ) const
+{
+ if( this->m_Offsets.size( ) == 0 )
+ {
+ if( this->m_Order == 1 )
+ this->_1stCombination( );
+ else
+ this->_2ndCombination( );
+
+ } // fi
+
+ TOutput res;
+ typename TImage::RegionType reg = this->m_Image->GetRequestedRegion( );
+ for( int i = 0; i < this->m_Offsets.size( ); ++i )
+ {
+ TIndex idx = center + this->m_Offsets[ i ];
+ if( reg.IsInside( idx ) && idx != center )
+ res.push_back( idx );
+
+ } // rof
+ return( res );
+}
+
+// -------------------------------------------------------------------------
+template< unsigned int _VDim >
+fpa::Image::Functors::SimpleNeighborhood< _VDim >::
+SimpleNeighborhood( )
+ : Superclass( ),
+ m_Order( 1 )
+{
+}
+
+// -------------------------------------------------------------------------
+template< unsigned int _VDim >
+ fpa::Image::Functors::SimpleNeighborhood< _VDim >::
+~SimpleNeighborhood( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< unsigned int _VDim >
+void fpa::Image::Functors::SimpleNeighborhood< _VDim >::
+_1stCombination( ) const
+{
+ for( int d = 0; d < TImage::ImageDimension; ++d )
+ {
+ typename TIndex::OffsetType off;
+ off.Fill( 0 );
+ for( int i = -1; i <= 1; i += 2 )
+ {
+ off[ d ] = i;
+ this->m_Offsets.push_back( off );
+
+ } // rof
+
+ } // rof
+}
+
+// -------------------------------------------------------------------------
+template< unsigned int _VDim >
+void fpa::Image::Functors::SimpleNeighborhood< _VDim >::
+_2ndCombination( ) const
+{
+ std::vector< std::vector< std::vector< int > > > M;
+ M.push_back( std::vector< std::vector< int > >( ) );
+
+ std::vector< int > base;
+ base.push_back( -1 );
+ base.push_back( 0 );
+ base.push_back( 1 );
+ int dim = TImage::ImageDimension;
+
+ M.push_back( std::vector< std::vector< int > >( ) );
+ for( int j = 0; j < base.size( ); ++j )
+ M[ 1 ].push_back( std::vector< int >( 1, base[ j ] ) );
+
+ for( int i = 2; i <= dim; ++i )
+ {
+ M.push_back( std::vector< std::vector< int > >( ) );
+ for( int j = 0; j < base.size( ); ++j )
+ {
+ for( int k = 0; k < M[ i - 1 ].size( ); ++k )
+ {
+ M[ i ].push_back( std::vector< int >( 1, base[ j ] ) );
+ for( int l = 0; l < M[ i - 1 ][ k ].size( ); ++l )
+ M[ i ].back( ).push_back( M[ i - 1 ][ k ][ l ] );
+
+ } // rof
+
+ } // rof
+
+ } // rof
+
+ for( int i = 0; i < M[ dim ].size( ); ++i )
+ {
+ TOffset off;
+ for( int j = 0; j < M[ dim ][ i ].size( ); ++j )
+ off[ j ] = M[ dim ][ i ][ j ];
+ this->m_Offsets.push_back( off );
+
+ } // rof
+}
+
+#endif // __fpa__Image__Functors__SimpleNeighborhood__hxx__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Image__Functors__VertexCost__h__
+#define __fpa__Image__Functors__VertexCost__h__
+
+#include <fpa/Image/Functors/Base.h>
+#include <fpa/Base/Functors/VertexCostFunctionBase.h>
+
+namespace fpa
+{
+ namespace Image
+ {
+ namespace Functors
+ {
+ /**
+ */
+ template< class _TImage, class _TOutput >
+ class VertexCost
+ : public fpa::Image::Functors::Base< _TImage, fpa::Base::Functors::VertexCostFunctionBase< typename _TImage::IndexType, _TOutput > >
+ {
+ public:
+ typedef _TImage TImage;
+ typedef typename TImage::IndexType TIndex;
+ typedef _TOutput TOutput;
+
+ typedef fpa::Base::Functors::VertexCostFunctionBase< TIndex, TOutput > TBaseFunctor;
+ typedef fpa::Image::Functors::Base< TImage, TBaseFunctor > Superclass;
+ typedef VertexCost Self;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ public:
+ itkNewMacro( Self );
+ itkTypeMacro( VertexCost, Base );
+
+ itkBooleanMacro( UseImageSpacing );
+ itkGetConstMacro( UseImageSpacing, bool );
+ itkSetMacro( UseImageSpacing, bool );
+
+ public:
+ virtual TOutput Evaluate(
+ const TIndex& a, const TIndex& b
+ ) const override;
+
+ protected:
+ VertexCost( );
+ virtual ~VertexCost( );
+
+ private:
+ // Purposely not implemented
+ VertexCost( const Self& other );
+ Self& operator=( const Self& other );
+
+ protected:
+ bool m_UseImageSpacing;
+ };
+
+ } // ecapseman
+
+ } // ecapseman
+
+} // ecapseman
+
+#ifndef ITK_MANUAL_INSTANTIATION
+# include <fpa/Image/Functors/VertexCost.hxx>
+#endif // ITK_MANUAL_INSTANTIATION
+
+#endif // __fpa__Image__Functors__VertexCost__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Image__Functors__VertexCost__hxx__
+#define __fpa__Image__Functors__VertexCost__hxx__
+
+// -------------------------------------------------------------------------
+template< class _TImage, class _TOutput >
+typename fpa::Image::Functors::VertexCost< _TImage, _TOutput >::
+TOutput fpa::Image::Functors::VertexCost< _TImage, _TOutput >::
+Evaluate( const TIndex& a, const TIndex& b ) const
+{
+ const _TImage* im =
+ dynamic_cast< const _TImage* >( this->m_Image.GetPointer( ) );
+ TOutput coeff = TOutput( 1 );
+ if( this->m_UseImageSpacing )
+ {
+ typename _TImage::PointType pa, pb;
+ im->TransformIndexToPhysicalPoint( a, pa );
+ im->TransformIndexToPhysicalPoint( b, pb );
+ coeff = TOutput( pa.EuclideanDistanceTo( pb ) );
+
+ } // fi
+ return( TOutput( im->GetPixel( a ) ) * coeff );
+}
+
+// -------------------------------------------------------------------------
+template< class _TImage, class _TOutput >
+fpa::Image::Functors::VertexCost< _TImage, _TOutput >::
+VertexCost( )
+ : Superclass( ),
+ m_UseImageSpacing( false )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TImage, class _TOutput >
+fpa::Image::Functors::VertexCost< _TImage, _TOutput >::
+~VertexCost( )
+{
+}
+
+#endif // __fpa__Image__Functors__VertexCost__hxx__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Image__MinimumSpanningTree__h__
+#define __fpa__Image__MinimumSpanningTree__h__
+
+#include <fpa/Base/MinimumSpanningTree.h>
+#include <cpExtensions/DataStructures/PolyLineParametricPath.h>
+#include <itkImage.h>
+
+namespace fpa
+{
+ namespace Image
+ {
+ /**
+ */
+ template< unsigned int _VDim >
+ class MinimumSpanningTree
+ : public fpa::Base::MinimumSpanningTree< itk::Index< _VDim >, cpExtensions::DataStructures::PolyLineParametricPath< _VDim >, itk::Image< itk::Offset< _VDim >, _VDim > >
+ {
+ public:
+ typedef itk::Index< _VDim > TVertex;
+ typedef cpExtensions::DataStructures::PolyLineParametricPath< _VDim > TPath;
+ typedef itk::Offset< _VDim > TOffset;
+ typedef itk::Image< TOffset, _VDim > TImage;
+ typedef fpa::Base::MinimumSpanningTree< TVertex, TPath, TImage > Superclass;
+ typedef MinimumSpanningTree Self;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ public:
+ itkNewMacro( Self );
+ itkTypeMacro(
+ fpa::Image::MinimumSpanningTree, fpa::Base::MinimumSpanningTree
+ );
+
+ public:
+ virtual TVertex GetParent( const TVertex& v ) const override;
+ virtual void SetParent( const TVertex& v, const TVertex& p ) override;
+
+ virtual void GetPath(
+ typename TPath::Pointer& path, const TVertex& a
+ ) const override;
+ virtual void GetPath(
+ typename TPath::Pointer& path, const TVertex& a, const TVertex& b
+ ) const override;
+
+ protected:
+ MinimumSpanningTree( );
+ virtual ~MinimumSpanningTree( );
+
+ private:
+ // Purposely not defined
+ MinimumSpanningTree( const Self& other );
+ Self& operator=( const Self& other );
+ };
+
+ } // ecapseman
+
+} // ecapseman
+
+#ifndef ITK_MANUAL_INSTANTIATION
+# include <fpa/Image/MinimumSpanningTree.hxx>
+#endif // ITK_MANUAL_INSTANTIATION
+
+#endif // __fpa__Image__MinimumSpanningTree__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Image__MinimumSpanningTree__hxx__
+#define __fpa__Image__MinimumSpanningTree__hxx__
+
+// -------------------------------------------------------------------------
+template< unsigned int _VDim >
+typename fpa::Image::MinimumSpanningTree< _VDim >::
+TVertex fpa::Image::MinimumSpanningTree< _VDim >::
+GetParent( const TVertex& v ) const
+{
+ TVertex p = v + this->GetPixel( v );
+ return( p );
+}
+
+// -------------------------------------------------------------------------
+template< unsigned int _VDim >
+void fpa::Image::MinimumSpanningTree< _VDim >::
+SetParent( const TVertex& v, const TVertex& p )
+{
+ this->SetPixel( v, p - v );
+}
+
+// -------------------------------------------------------------------------
+template< unsigned int _VDim >
+void fpa::Image::MinimumSpanningTree< _VDim >::
+GetPath( typename TPath::Pointer& path, const TVertex& a ) const
+{
+ if( path.GetPointer( ) == NULL )
+ path = TPath::New( );
+ path->SetReferenceImage( this );
+ this->Superclass::GetPath( path, a );
+}
+
+// -------------------------------------------------------------------------
+template< unsigned int _VDim >
+void fpa::Image::MinimumSpanningTree< _VDim >::
+GetPath(
+ typename TPath::Pointer& path, const TVertex& a, const TVertex& b
+ ) const
+{
+ if( path.GetPointer( ) == NULL )
+ path = TPath::New( );
+ path->SetReferenceImage( this );
+ this->Superclass::GetPath( path, a, b );
+}
+
+// -------------------------------------------------------------------------
+template< unsigned int _VDim >
+fpa::Image::MinimumSpanningTree< _VDim >::
+MinimumSpanningTree( )
+ : Superclass( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< unsigned int _VDim >
+fpa::Image::MinimumSpanningTree< _VDim >::
+~MinimumSpanningTree( )
+{
+}
+
+#endif // __fpa__Image__MinimumSpanningTree__hxx__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Image__MoriRegionGrow__h__
+#define __fpa__Image__MoriRegionGrow__h__
+
+#include <itkImageToImageFilter.h>
+#include <itkBinaryThresholdImageFilter.h>
+#include <fpa/Image/MoriRegionGrowHelper.h>
+
+namespace fpa
+{
+ namespace Image
+ {
+ /**
+ */
+ template< class _TInputImage, class _TOutputImage, class _TAuxPixel = unsigned short >
+ class MoriRegionGrow
+ : public itk::ImageToImageFilter< _TInputImage, _TOutputImage >
+ {
+ public:
+ typedef MoriRegionGrow Self;
+ typedef itk::ImageToImageFilter< _TInputImage, _TOutputImage > Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ typedef itk::Image< _TAuxPixel, _TInputImage::ImageDimension > TAuxImage;
+ typedef fpa::Image::MoriRegionGrowHelper< _TInputImage, TAuxImage > THelper;
+ typedef itk::BinaryThresholdImageFilter< TAuxImage, _TOutputImage > TThreshold;
+
+ typedef typename _TInputImage::IndexType TIndex;
+ typedef typename _TInputImage::PixelType TInputPixel;
+ typedef typename _TOutputImage::PixelType TOutputPixel;
+
+ public:
+ itkNewMacro( Self );
+ itkTypeMacro( fpa::Image::MoriRegionGrow, itk::ImageToImageFilter );
+
+ itkGetConstMacro( Seed, TIndex );
+ itkSetMacro( Seed, TIndex );
+
+ public:
+ TAuxImage* GetAuxiliaryImage( );
+ const TAuxImage* GetAuxiliaryImage( ) const;
+
+ TInputPixel GetLower( ) const;
+ TInputPixel GetUpper( ) const;
+ TInputPixel GetStep( ) const;
+ TOutputPixel GetInsideValue( ) const;
+ TOutputPixel GetOutsideValue( ) const;
+
+ void SetLower( const TInputPixel& v );
+ void SetUpper( const TInputPixel& v );
+ void SetStep( const TInputPixel& v );
+ void SetInsideValue( const TOutputPixel& v );
+ void SetOutsideValue( const TOutputPixel& v );
+
+
+ /* TODO
+ itkGetConstMacro( Lower, TPixel );
+ itkGetConstMacro( Upper, TPixel );
+ itkGetConstMacro( Step, TPixel );
+ itkGetConstMacro( Sensitivity, double );
+ itkSetMacro( Lower, TPixel );
+ itkSetMacro( Upper, TPixel );
+ itkSetMacro( Step, TPixel );
+ itkSetMacro( Sensitivity, double );
+ */
+
+ protected:
+ MoriRegionGrow( );
+ virtual ~MoriRegionGrow( );
+
+ virtual void GenerateData( ) override;
+
+ private:
+ // Purposely not defined
+ MoriRegionGrow( const Self& other );
+ Self& operator=( const Self& other );
+
+ protected:
+ typename THelper::Pointer m_Helper;
+ typename TThreshold::Pointer m_Threshold;
+ TIndex m_Seed;
+ };
+
+ } // ecapseman
+
+} // ecapseman
+
+#ifndef ITK_MANUAL_INSTANTIATION
+# include <fpa/Image/MoriRegionGrow.hxx>
+#endif // ITK_MANUAL_INSTANTIATION
+
+#endif // __fpa__Image__MoriRegionGrow__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Image__MoriRegionGrow__hxx__
+#define __fpa__Image__MoriRegionGrow__hxx__
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
+typename
+fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
+TAuxImage*
+fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
+GetAuxiliaryImage( )
+{
+ return( this->m_Helper->GetOutput( ) );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
+const typename
+fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
+TAuxImage*
+fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
+GetAuxiliaryImage( ) const
+{
+ return( this->m_Helper->GetOutput( ) );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
+typename
+fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
+TInputPixel
+fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
+GetLower( ) const
+{
+ return( this->m_Helper->GetLower( ) );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
+typename
+fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
+TInputPixel
+fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
+GetUpper( ) const
+{
+ return( this->m_Helper->GetUpper( ) );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
+typename
+fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
+TInputPixel
+fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
+GetStep( ) const
+{
+ return( this->m_Helper->GetStep( ) );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
+typename
+fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
+TOutputPixel
+fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
+GetInsideValue( ) const
+{
+ return( this->m_Threshold->GetInsideValue( ) );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
+typename
+fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
+TOutputPixel
+fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
+GetOutsideValue( ) const
+{
+ return( this->m_Threshold->GetInsideValue( ) );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
+void fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
+SetLower( const TInputPixel& v )
+{
+ this->m_Helper->SetLower( v );
+ this->Modified( );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
+void fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
+SetUpper( const TInputPixel& v )
+{
+ this->m_Helper->SetUpper( v );
+ this->Modified( );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
+void fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
+SetStep( const TInputPixel& v )
+{
+ this->m_Helper->SetStep( v );
+ this->Modified( );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
+void fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
+SetInsideValue( const TOutputPixel& v )
+{
+ this->m_Threshold->SetInsideValue( v );
+ this->Modified( );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
+void fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
+SetOutsideValue( const TOutputPixel& v )
+{
+ this->m_Threshold->SetOutsideValue( v );
+ this->Modified( );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
+fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
+MoriRegionGrow( )
+ : Superclass( )
+{
+ this->m_Helper = THelper::New( );
+ this->m_Threshold = TThreshold::New( );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
+fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
+~MoriRegionGrow( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
+void fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
+GenerateData( )
+{
+ this->m_Helper->ClearSeeds( );
+ this->m_Helper->AddSeed( this->m_Seed, 0 );
+ this->m_Helper->SetInput( this->GetInput( ) );
+ this->m_Helper->Update( );
+
+ this->m_Threshold->SetInput( this->m_Helper->GetOutput( ) );
+ this->m_Threshold->SetLowerThreshold( std::numeric_limits< _TAuxPixel >::min( ) );
+ this->m_Threshold->SetUpperThreshold( this->m_Helper->GetOptimumThreshold( ) );
+ this->m_Threshold->Update( );
+ this->GetOutput( )->Graft( this->m_Threshold->GetOutput( ) );
+}
+
+#endif // __fpa__Image__MoriRegionGrow__hxx__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Image__MoriRegionGrowHelper__h__
+#define __fpa__Image__MoriRegionGrowHelper__h__
+
+#include <utility>
+#include <vector>
+#include <fpa/Image/RegionGrow.h>
+#include <fpa/Image/Functors/RegionGrow/BinaryThreshold.h>
+
+namespace fpa
+{
+ namespace Image
+ {
+ /**
+ */
+ template< class _TInputImage, class _TOutputImage >
+ class MoriRegionGrowHelper
+ : public fpa::Image::RegionGrow< _TInputImage, _TOutputImage >
+ {
+ public:
+ typedef MoriRegionGrowHelper Self;
+ typedef fpa::Image::RegionGrow< _TInputImage, _TOutputImage > Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ typedef typename Superclass::TOutput TOutput;
+ typedef typename Superclass::TVertex TVertex;
+ typedef typename Superclass::TGrowFunction TGrowFunction;
+ typedef typename _TInputImage::PixelType TPixel;
+ typedef fpa::Image::Functors::RegionGrow::BinaryThreshold< _TInputImage, TOutput > TBinThresholdFunction;
+
+ protected:
+ typedef typename Superclass::_TQueueNode _TQueueNode;
+ typedef typename Superclass::_TQueue _TQueue;
+
+ typedef std::pair< TPixel, unsigned long > TCurveData;
+ typedef std::vector< TCurveData > TCurve;
+
+ public:
+ itkNewMacro( Self );
+ itkTypeMacro( fpa::Image::MoriRegionGrowHelper, fpa::Image::RegionGrow );
+
+ itkGetConstMacro( Lower, TPixel );
+ itkGetConstMacro( Upper, TPixel );
+ itkGetConstMacro( Step, TPixel );
+ itkGetConstMacro( OptimumThreshold, typename _TOutputImage::PixelType );
+
+ itkSetMacro( Lower, TPixel );
+ itkSetMacro( Upper, TPixel );
+ itkSetMacro( Step, TPixel );
+
+
+ protected:
+ MoriRegionGrowHelper( );
+ virtual ~MoriRegionGrowHelper( );
+
+ virtual bool _ContinueGenerateData( ) override;
+ virtual void _BeforeGenerateData( ) override;
+ virtual void _AfterGenerateData( ) override;
+ virtual void _BeforeLoop( ) override;
+ virtual void _AfterLoop( ) override;
+ virtual bool _UpdateValue( _TQueueNode& v, const _TQueueNode& p ) override;
+ virtual void _UpdateResult( const _TQueueNode& n ) override;
+
+ private:
+ // Purposely not defined
+ MoriRegionGrowHelper( const Self& other );
+ Self& operator=( const Self& other );
+
+ protected:
+ TPixel m_Lower;
+ TPixel m_Upper;
+ TPixel m_Step;
+ typename _TOutputImage::PixelType m_OptimumThreshold;
+
+ _TQueue m_NextQueue;
+ unsigned long m_ActualCount;
+ TCurve m_Curve;
+ };
+
+ } // ecapseman
+
+} // ecapseman
+
+#ifndef ITK_MANUAL_INSTANTIATION
+# include <fpa/Image/MoriRegionGrowHelper.hxx>
+#endif // ITK_MANUAL_INSTANTIATION
+
+#endif // __fpa__Image__MoriRegionGrowHelper__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Image__MoriRegionGrowHelper__hxx__
+#define __fpa__Image__MoriRegionGrowHelper__hxx__
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+fpa::Image::MoriRegionGrowHelper< _TInputImage, _TOutputImage >::
+MoriRegionGrowHelper( )
+ : Superclass( ),
+ m_Step( TPixel( 1 ) )
+{
+ this->m_Upper = std::numeric_limits< TPixel >::max( );
+ if( std::numeric_limits< TPixel >::is_integer )
+ this->m_Lower = std::numeric_limits< TPixel >::min( );
+ else
+ this->m_Lower = -this->m_Upper;
+ typename TBinThresholdFunction::Pointer functor =
+ TBinThresholdFunction::New( );
+ this->SetGrowFunction( functor );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+fpa::Image::MoriRegionGrowHelper< _TInputImage, _TOutputImage >::
+~MoriRegionGrowHelper( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+bool fpa::Image::MoriRegionGrowHelper< _TInputImage, _TOutputImage >::
+_ContinueGenerateData( )
+{
+ TBinThresholdFunction* functor =
+ dynamic_cast< TBinThresholdFunction* >( this->GetGrowFunction( ) );
+ TPixel u = functor->GetUpper( );
+
+ // Update flooding data
+ this->m_Curve.push_back( TCurveData( u, this->m_ActualCount ) );
+
+ // Update thresholds
+ if( u < this->m_Upper )
+ {
+ u += this->m_Step;
+ if( u > this->m_Upper )
+ u = this->m_Upper;
+ functor->SetUpper( u );
+ this->m_Queue = this->m_NextQueue;
+ while( this->m_NextQueue.size( ) > 0 )
+ this->m_NextQueue.pop( );
+ return( true );
+ }
+ else
+ return( false );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+void fpa::Image::MoriRegionGrowHelper< _TInputImage, _TOutputImage >::
+_BeforeGenerateData( )
+{
+ this->Superclass::_BeforeGenerateData( );
+
+ while( this->m_NextQueue.size( ) > 0 )
+ this->m_NextQueue.pop( );
+ this->m_OptimumThreshold = ( typename _TOutputImage::PixelType )( 0 );
+ this->m_ActualCount = 0;
+ this->m_Curve.clear( );
+ TBinThresholdFunction* functor =
+ dynamic_cast< TBinThresholdFunction* >( this->GetGrowFunction( ) );
+ functor->SetLower( this->m_Lower );
+ functor->SetUpper( this->m_Lower + this->m_Step );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+void fpa::Image::MoriRegionGrowHelper< _TInputImage, _TOutputImage >::
+_AfterGenerateData( )
+{
+ typedef typename _TOutputImage::PixelType _TOut;
+
+ this->Superclass::_AfterGenerateData( );
+ while( this->m_NextQueue.size( ) > 0 )
+ this->m_NextQueue.pop( );
+
+ // Find optimum threshold by dicotomy
+ unsigned long l = 0;
+ unsigned long r = this->m_Curve.size( ) - 1;
+ while( ( r - l ) > 1 )
+ {
+ unsigned long m = ( r + l ) >> 1;
+ double vm = double( this->m_Curve[ m ].second );
+ double dl = vm - double( this->m_Curve[ l ].second );
+ double dr = double( this->m_Curve[ r ].second ) - vm;
+ if( dl > dr )
+ r = m;
+ else
+ l = m;
+
+ } // elihw
+ this->m_OptimumThreshold = _TOut( ( r + l ) >> 1 );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+void fpa::Image::MoriRegionGrowHelper< _TInputImage, _TOutputImage >::
+_BeforeLoop( )
+{
+ this->Superclass::_BeforeLoop( );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+void fpa::Image::MoriRegionGrowHelper< _TInputImage, _TOutputImage >::
+_AfterLoop( )
+{
+ this->Superclass::_AfterLoop( );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+bool fpa::Image::MoriRegionGrowHelper< _TInputImage, _TOutputImage >::
+_UpdateValue( _TQueueNode& v, const _TQueueNode& p )
+{
+ typedef typename _TOutputImage::PixelType _TOut;
+
+ bool ret = this->Superclass::_UpdateValue( v, p );
+ v.Result = _TOut( this->m_Curve.size( ) + 1 );
+ if( !ret )
+ this->m_NextQueue.push( v );
+ return( ret );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+void fpa::Image::MoriRegionGrowHelper< _TInputImage, _TOutputImage >::
+_UpdateResult( const _TQueueNode& n )
+{
+ this->Superclass::_UpdateResult( n );
+ this->m_ActualCount += 1;
+}
+
+#endif // __fpa__Image__MoriRegionGrowHelper__hxx__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Image__RegionGrow__h__
+#define __fpa__Image__RegionGrow__h__
+
+#include <fpa/Base/RegionGrow.h>
+#include <fpa/Image/Algorithm.h>
+
+namespace fpa
+{
+ namespace Image
+ {
+ /**
+ */
+ template< class _TInputImage, class _TOutputImage >
+ class RegionGrow
+ : public fpa::Base::RegionGrow< fpa::Image::Algorithm< _TInputImage, _TOutputImage > >
+ {
+ public:
+ typedef fpa::Image::Algorithm< _TInputImage, _TOutputImage > TAlgorithm;
+ typedef RegionGrow Self;
+ typedef fpa::Base::RegionGrow< TAlgorithm > Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ typedef typename Superclass::TOutput TOutput;
+ typedef typename Superclass::TVertex TVertex;
+
+ typedef fpa::Image::Functors::Base< _TInputImage, typename Superclass::TGrowFunction > TGrowFunction;
+
+ public:
+ itkNewMacro( Self );
+ itkTypeMacro( fpa::Image::RegionGrow, fpa::Base::RegionGrow );
+
+ protected:
+ RegionGrow( );
+ virtual ~RegionGrow( );
+
+ virtual void _BeforeGenerateData( ) override;
+
+ private:
+ // Purposely not defined
+ RegionGrow( const Self& other );
+ Self& operator=( const Self& other );
+ };
+
+ } // ecapseman
+
+} // ecapseman
+
+#ifndef ITK_MANUAL_INSTANTIATION
+# include <fpa/Image/RegionGrow.hxx>
+#endif // ITK_MANUAL_INSTANTIATION
+
+#endif // __fpa__Image__RegionGrow__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Image__RegionGrow__hxx__
+#define __fpa__Image__RegionGrow__hxx__
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+fpa::Image::RegionGrow< _TInputImage, _TOutputImage >::
+RegionGrow( )
+ : Superclass( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+fpa::Image::RegionGrow< _TInputImage, _TOutputImage >::
+~RegionGrow( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TOutputImage >
+void fpa::Image::RegionGrow< _TInputImage, _TOutputImage >::
+_BeforeGenerateData( )
+{
+ this->Superclass::_BeforeGenerateData( );
+ this->m_InitResult = this->GetOutsideValue( );
+ TGrowFunction* grow =
+ dynamic_cast< TGrowFunction* >( this->GetGrowFunction( ) );
+ if( grow != NULL )
+ grow->SetImage( this->GetInput( ) );
+}
+
+#endif // __fpa__Image__RegionGrow__hxx__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Image__SkeletonFilter__h__
+#define __fpa__Image__SkeletonFilter__h__
+
+#include <map>
+#include <vector>
+#include <fpa/Image/Dijkstra.h>
+#include <cpExtensions/DataStructures/Skeleton.h>
+#include <itkSimpleDataObjectDecorator.h>
+
+namespace fpa
+{
+ namespace Image
+ {
+ /**
+ */
+ template< class _TImage >
+ class SkeletonFilter
+ : public Dijkstra< _TImage, _TImage >
+ {
+ public:
+ typedef SkeletonFilter Self;
+ typedef Dijkstra< _TImage, _TImage > Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ typedef _TImage TImage;
+ itkStaticConstMacro(
+ Dimension, unsigned int, TImage::ImageDimension
+ );
+
+ typedef typename Superclass::TMST TMST;
+ typedef typename TImage::IndexType TIndex;
+ typedef itk::Image< unsigned char, Self::Dimension > TMarks;
+ typedef cpExtensions::DataStructures::Skeleton< Self::Dimension > TSkeleton;
+
+ protected:
+ typedef typename Superclass::_TQueueNode _TQueueNode;
+ typedef std::multimap< double, TIndex, std::greater< double > > _TSkeletonQueue;
+ typedef std::map< TIndex, TIndex, typename TIndex::LexicographicCompare > _TAdjacencies;
+
+ public:
+ itkNewMacro( Self );
+ itkTypeMacro( fpa::Image::SkeletonFilter, itk::Object );
+
+ public:
+ TSkeleton* GetSkeleton( );
+ TMarks* GetMarks( );
+
+ protected:
+ SkeletonFilter( );
+ virtual ~SkeletonFilter( );
+
+ virtual void _BeforeGenerateData( ) override;
+ virtual void _UpdateResult( const _TQueueNode& n ) override;
+ virtual void _AfterGenerateData( ) override;
+
+ void _EndPoints( std::vector< TIndex >& end_points, _TAdjacencies& A );
+ void _Skeleton( const std::vector< TIndex >& end_points, _TAdjacencies& A );
+ void _MarkSphere( const TIndex& idx );
+
+ private:
+ // Purposely not defined
+ SkeletonFilter( const Self& other );
+ Self& operator=( const Self& other );
+
+ protected:
+ _TSkeletonQueue m_SkeletonQueue;
+
+ unsigned long m_SkeletonIdx;
+ unsigned long m_MarksIdx;
+ };
+
+ } // ecapseman
+
+} // ecapseman
+
+#ifndef ITK_MANUAL_INSTANTIATION
+# include <fpa/Image/SkeletonFilter.hxx>
+#endif // ITK_MANUAL_INSTANTIATION
+
+#endif // __fpa__Image__SkeletonFilter__h__
+
+// eof - $RCSfile$
--- /dev/null
+#ifndef __fpa__Image__SkeletonFilter__hxx__
+#define __fpa__Image__SkeletonFilter__hxx__
+
+#include <fpa/Base/Functors/Inverse.h>
+#include <fpa/Image/Functors/SimpleNeighborhood.h>
+#include <itkImageRegionIteratorWithIndex.h>
+
+// -------------------------------------------------------------------------
+#define fpa_Image_SkeletonFilter_OutputMacro( o_n, o_t ) \
+ template< class _TImage > \
+ typename fpa::Image::SkeletonFilter< _TImage >:: \
+ o_t* fpa::Image::SkeletonFilter< _TImage >:: \
+ Get##o_n( ) \
+ { \
+ return( \
+ dynamic_cast< o_t* >( \
+ this->itk::ProcessObject::GetOutput( this->m_##o_n##Idx ) \
+ ) \
+ ); \
+ }
+
+fpa_Image_SkeletonFilter_OutputMacro( Skeleton, TSkeleton );
+fpa_Image_SkeletonFilter_OutputMacro( Marks, TMarks );
+
+// -------------------------------------------------------------------------
+template< class _TImage >
+fpa::Image::SkeletonFilter< _TImage >::
+SkeletonFilter( )
+ : Superclass( )
+{
+ typedef typename _TImage::PixelType _TPixel;
+ typedef fpa::Image::Functors::SimpleNeighborhood< _TImage::ImageDimension > _TNeighFunc;
+ typedef fpa::Base::Functors::Inverse< _TPixel, _TPixel > _TInvFunc;
+
+ unsigned int nOutputs = this->GetNumberOfRequiredOutputs( );
+ this->SetNumberOfRequiredOutputs( nOutputs + 2 );
+ this->m_SkeletonIdx = nOutputs;
+ this->m_MarksIdx = nOutputs + 1;
+
+ typename TSkeleton::Pointer skeleton = TSkeleton::New( );
+ typename TMarks::Pointer marks = TMarks::New( );
+ this->SetNthOutput( this->m_SkeletonIdx, skeleton.GetPointer( ) );
+ this->SetNthOutput( this->m_MarksIdx, marks.GetPointer( ) );
+
+ typename _TNeighFunc::Pointer nfunc = _TNeighFunc::New( );
+ nfunc->SetOrder( 2 );
+ this->SetNeighborhoodFunction( nfunc );
+
+ typename _TInvFunc::Pointer ifunc = _TInvFunc::New( );
+ this->SetConversionFunction( ifunc );
+}
+
+// -------------------------------------------------------------------------
+template< class _TImage >
+fpa::Image::SkeletonFilter< _TImage >::
+~SkeletonFilter( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TImage >
+void fpa::Image::SkeletonFilter< _TImage >::
+_BeforeGenerateData( )
+{
+ this->Superclass::_BeforeGenerateData( );
+ this->m_SkeletonQueue.clear( );
+}
+
+// -------------------------------------------------------------------------
+template< class _TImage >
+void fpa::Image::SkeletonFilter< _TImage >::
+_UpdateResult( const _TQueueNode& n )
+{
+ typedef typename _TSkeletonQueue::value_type _TSkeletonQueueValue;
+
+ this->Superclass::_UpdateResult( n );
+ double d = double( this->GetInput( )->GetPixel( n.Vertex ) );
+ if( d >= double( 0 ) )
+ {
+ // Update skeleton candidates
+ d += double( 1e-5 );
+ double v = double( n.Result ) / ( d * d );
+ this->m_SkeletonQueue.insert( _TSkeletonQueueValue( v, n.Vertex ) );
+
+ } // fi
+}
+
+// -------------------------------------------------------------------------
+template< class _TImage >
+void fpa::Image::SkeletonFilter< _TImage >::
+_AfterGenerateData( )
+{
+ this->Superclass::_AfterGenerateData( );
+
+ _TAdjacencies A;
+ std::vector< TIndex > end_points;
+ this->_EndPoints( end_points, A );
+ this->_Skeleton( end_points, A );
+}
+
+// -------------------------------------------------------------------------
+template< class _TImage >
+void fpa::Image::SkeletonFilter< _TImage >::
+_EndPoints( std::vector< TIndex >& end_points, _TAdjacencies& A )
+{
+ auto marks = this->GetMarks( );
+ auto mst = this->GetMinimumSpanningTree( );
+ auto spac = marks->GetSpacing( );
+
+ // Some values
+ marks->SetLargestPossibleRegion( mst->GetLargestPossibleRegion( ) );
+ marks->SetRequestedRegion( mst->GetRequestedRegion( ) );
+ marks->SetBufferedRegion( mst->GetBufferedRegion( ) );
+ marks->SetSpacing( mst->GetSpacing( ) );
+ marks->SetOrigin( mst->GetOrigin( ) );
+ marks->SetDirection( mst->GetDirection( ) );
+ marks->Allocate( );
+ marks->FillBuffer( 0 );
+
+ // BFS from maximum queue
+ while( this->m_SkeletonQueue.size( ) > 0 )
+ {
+ // Get node
+ auto nIt = this->m_SkeletonQueue.begin( );
+ auto n = *nIt;
+ this->m_SkeletonQueue.erase( nIt );
+
+ // Mark it and update end-points
+ unsigned char m = marks->GetPixel( n.second );
+ if( m != 0 )
+ continue;
+ marks->SetPixel( n.second, 1 );
+ end_points.push_back( n.second );
+
+ // Mark path
+ TIndex it = n.second;
+ TIndex p = mst->GetParent( it );
+ while( it != p )
+ {
+ this->_MarkSphere( it );
+ it = p;
+ p = mst->GetParent( it );
+
+ } // elihw
+ this->_MarkSphere( it );
+ A[ n.second ] = it;
+
+ } // elihw
+}
+
+// -------------------------------------------------------------------------
+template< class _TImage >
+void fpa::Image::SkeletonFilter< _TImage >::
+_Skeleton( const std::vector< TIndex >& end_points, _TAdjacencies& A )
+{
+ typedef itk::Image< unsigned long, Self::Dimension > _TTagsImage;
+ typedef typename TMST::TPath _TPath;
+
+ auto mst = this->GetMinimumSpanningTree( );
+ auto skeleton = this->GetSkeleton( );
+
+ // Tag branches
+ typename _TTagsImage::Pointer tags = _TTagsImage::New( );
+ tags->SetLargestPossibleRegion( mst->GetLargestPossibleRegion( ) );
+ tags->SetRequestedRegion( mst->GetRequestedRegion( ) );
+ tags->SetBufferedRegion( mst->GetBufferedRegion( ) );
+ tags->Allocate( );
+ tags->FillBuffer( 0 );
+ for( auto eIt = end_points.begin( ); eIt != end_points.end( ); ++eIt )
+ {
+ TIndex it = *eIt;
+ TIndex p = mst->GetParent( it );
+ while( it != p )
+ {
+ tags->SetPixel( it, tags->GetPixel( it ) + 1 );
+ it = p;
+ p = mst->GetParent( it );
+
+ } // elihw
+ tags->SetPixel( it, tags->GetPixel( it ) + 1 );
+
+ } // rof
+
+ // Build paths (branches)
+ for( auto eIt = end_points.begin( ); eIt != end_points.end( ); ++eIt )
+ {
+ TIndex it = *eIt;
+ TIndex p = mst->GetParent( it );
+ TIndex sIdx = *eIt;
+ typename _TPath::Pointer path = _TPath::New( );
+ path->SetReferenceImage( mst );
+ while( it != p )
+ {
+ if( tags->GetPixel( sIdx ) != tags->GetPixel( it ) )
+ {
+ // Ok, a new branch has been added
+ path->AddVertex( it );
+ skeleton->AddBranch( path );
+
+ // Update a new starting index
+ path = _TPath::New( );
+ path->SetReferenceImage( mst );
+ sIdx = it;
+ }
+ else
+ path->AddVertex( it );
+ it = p;
+ p = mst->GetParent( it );
+
+ } // elihw
+
+ // Finally, add last branch
+ path->AddVertex( it );
+ skeleton->AddBranch( path );
+
+ } // rof
+}
+
+// -------------------------------------------------------------------------
+template< class _TImage >
+void fpa::Image::SkeletonFilter< _TImage >::
+_MarkSphere( const TIndex& idx )
+{
+ typedef itk::ImageRegionIteratorWithIndex< TMarks > _TMarksIt;
+
+ static const double _eps = std::sqrt( double( Self::Dimension + 1 ) );
+ auto input = this->GetInput( );
+ auto marks = this->GetMarks( );
+ auto spac = input->GetSpacing( );
+ auto region = input->GetRequestedRegion( );
+
+ typename _TImage::PointType cnt;
+ input->TransformIndexToPhysicalPoint( idx, cnt );
+ double r = double( input->GetPixel( idx ) ) * _eps;
+
+ TIndex i0, i1;
+ for( unsigned int d = 0; d < Self::Dimension; ++d )
+ {
+ long off = long( std::ceil( r / double( spac[ d ] ) ) );
+ if( off < 3 )
+ off = 3;
+ i0[ d ] = idx[ d ] - off;
+ i1[ d ] = idx[ d ] + off;
+
+ if( i0[ d ] < region.GetIndex( )[ d ] )
+ i0[ d ] = region.GetIndex( )[ d ];
+
+ if( i1[ d ] >= region.GetIndex( )[ d ] + region.GetSize( )[ d ] )
+ i1[ d ] = region.GetIndex( )[ d ] + region.GetSize( )[ d ] - 1;
+
+ } // rof
+
+ typename _TImage::SizeType size;
+ for( unsigned int d = 0; d < Self::Dimension; ++d )
+ size[ d ] = i1[ d ] - i0[ d ] + 1;
+
+ typename _TImage::RegionType neighRegion;
+ neighRegion.SetIndex( i0 );
+ neighRegion.SetSize( size );
+
+ _TMarksIt mIt( marks, neighRegion );
+ for( mIt.GoToBegin( ); !mIt.IsAtEnd( ); ++mIt )
+ {
+ TIndex w = mIt.GetIndex( );
+ typename _TImage::PointType p;
+ marks->TransformIndexToPhysicalPoint( w, p );
+ mIt.Set( 1 );
+
+ } // rof
+}
+
+#endif // __fpa__Image__SkeletonFilter__hxx__
+
+// eof - $RCSfile$
git://git.creatis.insa-lyon.fr / FrontAlgorithms.git / commitdiff