RegionGrow_Tautology
RegionGrow_BinaryThreshold
RegionGrow_Mori
+ Dijkstra_Maurer
#CreateMoriInputImage
#BronchiiInitialSegmentationWithMori
#BronchiiInitialSegmentationWithBinaryThresholdRegionGrow
--- /dev/null
+#include <itkImage.h>
+#include <itkImageFileReader.h>
+#include <itkImageFileWriter.h>
+#include <itkMinimumMaximumImageCalculator.h>
+#include <itkSignedMaurerDistanceMapImageFilter.h>
+
+#include <fpa/Image/Dijkstra.h>
+#include <fpa/Image/Functors/VertexIdentity.h>
+#include <fpa/Base/Functors/InvertValue.h>
+
+// -------------------------------------------------------------------------
+static const unsigned int VDim = 2;
+typedef short TPixel;
+typedef double TScalar;
+typedef itk::Image< TPixel, VDim > TImage;
+typedef itk::Image< TScalar, VDim > TScalarImage;
+typedef itk::ImageFileReader< TImage > TReader;
+typedef itk::ImageFileWriter< TScalarImage > TWriter;
+typedef fpa::Image::Dijkstra< TScalarImage, TScalarImage > TFilter;
+typedef itk::MinimumMaximumImageCalculator< TImage > TMinMax;
+typedef itk::SignedMaurerDistanceMapImageFilter< TImage, TScalarImage > TDMap;
+
+typedef fpa::Image::Functors::VertexIdentity< TScalarImage, TScalar > TVertexFunc;
+typedef fpa::Base::Functors::InvertValue< TScalar, TScalar > TValueFunc;
+
+typedef TFilter::TMST TMST;
+typedef itk::ImageFileWriter< TMST > TMSTWriter;
+
+// -------------------------------------------------------------------------
+int main( int argc, char* argv[] )
+{
+ // Get arguments
+ if( argc < 5 + VDim )
+ {
+ std::cerr
+ << "Usage: " << argv[ 0 ]
+ << " input_image output_image output_mst stop_at_one_front";
+ for( unsigned int i = 0; i < VDim; ++i )
+ std::cerr << " s_" << i;
+ std::cerr << " ..." << std::endl;
+ return( 1 );
+
+ } // fi
+ std::string input_image_filename = argv[ 1 ];
+ std::string output_image_filename = argv[ 2 ];
+ std::string output_mst_filename = argv[ 3 ];
+ bool stop_at_one_front = ( std::atoi( argv[ 4 ] ) == 1 );
+
+ TReader::Pointer reader = TReader::New( );
+ reader->SetFileName( input_image_filename );
+ try
+ {
+ reader->Update( );
+ }
+ catch( std::exception& err )
+ {
+ std::cerr << "Error caught: " << err.what( ) << std::endl;
+ return( 1 );
+
+ } // yrt
+
+ TMinMax::Pointer minmax = TMinMax::New( );
+ minmax->SetImage( reader->GetOutput( ) );
+ minmax->Compute( );
+
+ TDMap::Pointer dmap = TDMap::New( );
+ dmap->SetInput( reader->GetOutput( ) );
+ dmap->SetBackgroundValue( minmax->GetMinimum( ) );
+ dmap->InsideIsPositiveOn( );
+ dmap->UseImageSpacingOn( );
+ dmap->SquaredDistanceOff( );
+
+ TFilter::Pointer filter = TFilter::New( );
+ filter->SetInput( dmap->GetOutput( ) );
+ filter->SetStopAtOneFront( stop_at_one_front );
+
+ for( int i = 5; i < argc; i += VDim )
+ {
+ if( i + VDim <= argc )
+ {
+ TImage::IndexType seed;
+ for( int j = 0; j < VDim; ++j )
+ seed[ j ] = std::atoi( argv[ i + j ] );
+ filter->AddSeed( seed );
+
+ } // fi
+
+ } // rof
+
+ TVertexFunc::Pointer vertex_func = TVertexFunc::New( );
+ filter->SetFunctor( vertex_func );
+
+ TValueFunc::Pointer value_func = TValueFunc::New( );
+ value_func->SetAlpha( 1 );
+ value_func->SetBeta( 1 );
+ filter->SetFunctor( value_func );
+
+ TWriter::Pointer writer = TWriter::New( );
+ writer->SetInput( filter->GetOutput( ) );
+ writer->SetFileName( output_image_filename );
+
+ TMSTWriter::Pointer mst_writer = TMSTWriter::New( );
+ mst_writer->SetInput( filter->GetMinimumSpanningTree( ) );
+ mst_writer->SetFileName( output_mst_filename );
+
+ try
+ {
+ writer->Update( );
+ mst_writer->Update( );
+ }
+ catch( std::exception& err )
+ {
+ std::cerr << "ERROR: " << err.what( ) << std::endl;
+ return( 1 );
+
+ } // yrt
+ return( 0 );
+}
+
+// eof - $RCSfile$
#define __fpa__Base__Dijkstra__h__
#include <itkFunctionBase.h>
+#include <fpa/Base/MinimumSpanningTree.h>
+#include <fpa/Base/Functors/VertexParentBase.h>
namespace fpa
{
{
/**
*/
- template< class _TFilter, class _TMarksInterface, class _TSeedsInterface >
+ template< class _TFilter, class _TMarksInterface, class _TSeedsInterface, class _TMST >
class Dijkstra
: public _TFilter,
public _TMarksInterface,
typedef _TFilter Superclass;
typedef _TMarksInterface TMarksInterface;
typedef _TSeedsInterface TSeedsInterface;
+ typedef _TMST TMST;
typedef itk::SmartPointer< Self > Pointer;
typedef itk::SmartPointer< const Self > ConstPointer;
typedef typename Superclass::TVertices TVertices;
typedef itk::FunctionBase< TInputValue, TOutputValue > TIntensityFunctor;
- typedef itk::FunctionBase< TVertex, TOutputValue > TVertexFunctor;
+ typedef fpa::Base::Functors::VertexParentBase< TVertex, TOutputValue > TVertexFunctor;
protected:
struct _TNode
TVertex Vertex;
TVertex Parent;
TOutputValue Cost;
- _TNode( const TVertex& v, const TVertex& p )
+ unsigned long FrontId;
+ _TNode( const TVertex& v, const TVertex& p, const unsigned long& fId )
{
this->Vertex = v;
this->Parent = p;
+ this->FrontId = fId;
this->Cost = TOutputValue( 0 );
}
bool operator<( const _TNode& b ) const
{
- return( this->Cost < b.Cost );
+ return( b.Cost < this->Cost );
}
};
itkTypeMacro( Dijkstra, TFilter );
public:
+ TMST* GetMinimumSpanningTree( );
+ const TMST* GetMinimumSpanningTree( ) const;
+
const TIntensityFunctor* GetIntensityFunctor( ) const;
const TVertexFunctor* GetVertexFunctor( ) const;
protected:
typename TIntensityFunctor::Pointer m_IntensityFunctor;
typename TVertexFunctor::Pointer m_VertexFunctor;
+ unsigned long m_MSTIndex;
};
} // ecapseman
#define __fpa__Base__Dijkstra__hxx__
// -------------------------------------------------------------------------
-template< class _TFilter, class _TMarksInterface, class _TSeedsInterface >
+template< class _TFilter, class _TMarksInterface, class _TSeedsInterface, class _TMST >
+typename
+fpa::Base::Dijkstra< _TFilter, _TMarksInterface, _TSeedsInterface, _TMST >::
+TMST* fpa::Base::Dijkstra< _TFilter, _TMarksInterface, _TSeedsInterface, _TMST >::
+GetMinimumSpanningTree( )
+{
+ return(
+ dynamic_cast< TMST* >(
+ this->itk::ProcessObject::GetOutput( this->m_MSTIndex )
+ )
+ );
+}
+
+// -------------------------------------------------------------------------
+template< class _TFilter, class _TMarksInterface, class _TSeedsInterface, class _TMST >
+const typename
+fpa::Base::Dijkstra< _TFilter, _TMarksInterface, _TSeedsInterface, _TMST >::
+TMST* fpa::Base::Dijkstra< _TFilter, _TMarksInterface, _TSeedsInterface, _TMST >::
+GetMinimumSpanningTree( ) const
+{
+ return(
+ dynamic_cast< const TMST* >(
+ this->itk::ProcessObject::GetOutput( this->m_MSTIndex )
+ )
+ );
+}
+
+// -------------------------------------------------------------------------
+template< class _TFilter, class _TMarksInterface, class _TSeedsInterface, class _TMST >
const typename
-fpa::Base::Dijkstra< _TFilter, _TMarksInterface, _TSeedsInterface >::
+fpa::Base::Dijkstra< _TFilter, _TMarksInterface, _TSeedsInterface, _TMST >::
TIntensityFunctor*
-fpa::Base::Dijkstra< _TFilter, _TMarksInterface, _TSeedsInterface >::
+fpa::Base::Dijkstra< _TFilter, _TMarksInterface, _TSeedsInterface, _TMST >::
GetIntensityFunctor( ) const
{
return( this->m_IntensityFunctor );
}
// -------------------------------------------------------------------------
-template< class _TFilter, class _TMarksInterface, class _TSeedsInterface >
+template< class _TFilter, class _TMarksInterface, class _TSeedsInterface, class _TMST >
const typename
-fpa::Base::Dijkstra< _TFilter, _TMarksInterface, _TSeedsInterface >::
+fpa::Base::Dijkstra< _TFilter, _TMarksInterface, _TSeedsInterface, _TMST >::
TVertexFunctor*
-fpa::Base::Dijkstra< _TFilter, _TMarksInterface, _TSeedsInterface >::
+fpa::Base::Dijkstra< _TFilter, _TMarksInterface, _TSeedsInterface, _TMST >::
GetVertexFunctor( ) const
{
return( this->m_VertexFunctor );
}
// -------------------------------------------------------------------------
-template< class _TFilter, class _TMarksInterface, class _TSeedsInterface >
-void fpa::Base::Dijkstra< _TFilter, _TMarksInterface, _TSeedsInterface >::
+template< class _TFilter, class _TMarksInterface, class _TSeedsInterface, class _TMST >
+void fpa::Base::Dijkstra< _TFilter, _TMarksInterface, _TSeedsInterface, _TMST >::
SetFunctor( TIntensityFunctor* functor )
{
if( this->m_IntensityFunctor.GetPointer( ) != functor )
}
// -------------------------------------------------------------------------
-template< class _TFilter, class _TMarksInterface, class _TSeedsInterface >
-void fpa::Base::Dijkstra< _TFilter, _TMarksInterface, _TSeedsInterface >::
+template< class _TFilter, class _TMarksInterface, class _TSeedsInterface, class _TMST >
+void fpa::Base::Dijkstra< _TFilter, _TMarksInterface, _TSeedsInterface, _TMST >::
SetFunctor( TVertexFunctor* functor )
{
if( this->m_VertexFunctor.GetPointer( ) != functor )
}
// -------------------------------------------------------------------------
-template< class _TFilter, class _TMarksInterface, class _TSeedsInterface >
-fpa::Base::Dijkstra< _TFilter, _TMarksInterface, _TSeedsInterface >::
+template< class _TFilter, class _TMarksInterface, class _TSeedsInterface, class _TMST >
+fpa::Base::Dijkstra< _TFilter, _TMarksInterface, _TSeedsInterface, _TMST >::
Dijkstra( )
: Superclass( ),
_TMarksInterface( this ),
_TSeedsInterface( this )
{
+ this->m_MSTIndex = this->GetNumberOfRequiredOutputs( );
+ this->SetNumberOfRequiredOutputs( this->m_MSTIndex + 1 );
+ this->itk::ProcessObject::SetNthOutput( this->m_MSTIndex, TMST::New( ) );
}
// -------------------------------------------------------------------------
-template< class _TFilter, class _TMarksInterface, class _TSeedsInterface >
-fpa::Base::Dijkstra< _TFilter, _TMarksInterface, _TSeedsInterface >::
+template< class _TFilter, class _TMarksInterface, class _TSeedsInterface, class _TMST >
+fpa::Base::Dijkstra< _TFilter, _TMarksInterface, _TSeedsInterface, _TMST >::
~Dijkstra( )
{
}
// -------------------------------------------------------------------------
-template< class _TFilter, class _TMarksInterface, class _TSeedsInterface >
-void fpa::Base::Dijkstra< _TFilter, _TMarksInterface, _TSeedsInterface >::
+template< class _TFilter, class _TMarksInterface, class _TSeedsInterface, class _TMST >
+void fpa::Base::Dijkstra< _TFilter, _TMarksInterface, _TSeedsInterface, _TMST >::
GenerateData( )
{
+ // Init objects
+ this->_ConfigureOutputs( TOutputValue( 0 ) );
+ this->_InitMarks( this->GetNumberOfSeeds( ) );
+ TMST* mst = this->GetMinimumSpanningTree( );
+
+ // Init queue
+ std::vector< _TNode > q;
+ unsigned long frontId = 1;
+ for( TVertex seed: this->GetSeeds( ) )
+ q.push_back( _TNode( seed, seed, frontId++ ) );
+
+ // Main loop
+ while( q.size( ) > 0 )
+ {
+ // Get next candidate
+ std::pop_heap( q.begin( ), q.end( ) );
+ _TNode node = q.back( );
+ q.pop_back( );
+ if( this->_IsMarked( node.Vertex ) )
+ continue;
+ this->_Mark( node.Vertex, node.FrontId );
+
+ // Ok, pixel lays inside region
+ this->_SetOutputValue( node.Vertex, node.Cost );
+ mst->SetParent( node.Vertex, node.Parent );
+
+ // Add neighborhood
+ TVertices neighbors = this->_GetNeighbors( node.Vertex );
+ for( TVertex neigh: neighbors )
+ {
+ if( this->_IsMarked( neigh ) )
+ {
+ // Invoke stop at collisions
+ unsigned long nColl = this->_Collisions( node.Vertex, neigh );
+ if(
+ this->StopAtOneFront( ) &&
+ this->GetNumberOfSeeds( ) > 1 &&
+ nColl == 1
+ )
+ q.clear( );
+ }
+ else
+ {
+ // Compute new cost
+ TOutputValue ncost =
+ this->m_VertexFunctor->Evaluate( neigh, node.Vertex );
+ if( this->m_IntensityFunctor.IsNotNull( ) )
+ ncost = this->m_IntensityFunctor->Evaluate( ncost );
+
+ // This algorithm only supports positive values
+ if( ncost >= TOutputValue( 0 ) )
+ {
+ // Insert new node
+ _TNode nn( neigh, node.Vertex, node.FrontId );
+ nn.Cost = node.Cost + ncost;
+ q.push_back( nn );
+ std::push_heap( q.begin( ), q.end( ) );
+
+ } // fi
+
+ } // fi
+
+ } // rof
+
+ } // elihw
+ this->_FreeMarks( );
}
#endif // __fpa__Base__Dijkstra__hxx__
--- /dev/null
+// =========================================================================
+// @author Leonardo Florez Valencia
+// @email florez-l@javeriana.edu.co
+// =========================================================================
+
+#ifndef __fpa__Base__Functors__InvertValue__h__
+#define __fpa__Base__Functors__InvertValue__h__
+
+#include <cmath>
+#include <itkFunctionBase.h>
+
+namespace fpa
+{
+ namespace Base
+ {
+ namespace Functors
+ {
+ /**
+ */
+ template< class _TInputValue, class _TOutputValue >
+ class InvertValue
+ : public itk::FunctionBase< _TInputValue, _TOutputValue >
+ {
+ public:
+ typedef InvertValue Self;
+ typedef itk::FunctionBase< _TInputValue, _TOutputValue > Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ typedef _TInputValue TInputValue;
+ typedef _TOutputValue TOutputValue;
+
+ public:
+ itkNewMacro( Self );
+ itkTypeMacro( InvertValue, itk::FunctionBase );
+
+ itkGetConstMacro( Alpha, double );
+ itkGetConstMacro( Beta, double );
+ itkSetMacro( Alpha, double );
+ itkSetMacro( Beta, double );
+
+ public:
+ virtual TOutputValue Evaluate( const TInputValue& a ) const override
+ {
+ double d = this->m_Alpha;
+ d += std::pow( double( a ), this->m_Beta );
+ double x = -1;
+ if( std::fabs( d ) > double( 0 ) )
+ x =
+ double( 1 ) /
+ ( this->m_Alpha + std::pow( double( a ), this->m_Beta ) );
+ return( TOutputValue( x ) );
+ }
+
+ protected:
+ InvertValue( )
+ : Superclass( ),
+ m_Alpha( double( 1 ) ),
+ m_Beta( double( 1 ) )
+ { }
+ virtual ~InvertValue( ) { }
+
+ private:
+ InvertValue( const Self& other );
+ Self& operator=( const Self& other );
+
+ protected:
+ double m_Alpha;
+ double m_Beta;
+ };
+
+ } // ecapseman
+
+ } // ecapseman
+
+} // ecapseman
+
+#endif // __fpa__Base__Functors__InvertValue__h__
+
+// eof - $RCSfile$
--- /dev/null
+// =========================================================================
+// @author Leonardo Florez Valencia
+// @email florez-l@javeriana.edu.co
+// =========================================================================
+
+#ifndef __fpa__Base__Functors__VertexParentBase__h__
+#define __fpa__Base__Functors__VertexParentBase__h__
+
+#include <itkObject.h>
+#include <itkObjectFactory.h>
+
+namespace fpa
+{
+ namespace Base
+ {
+ namespace Functors
+ {
+ /**
+ */
+ template< class _TVertex, class _TOutputValue >
+ class VertexParentBase
+ : public itk::Object
+ {
+ public:
+ typedef VertexParentBase Self;
+ typedef itk::Object Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ typedef _TVertex TVertex;
+ typedef _TOutputValue TOutputValue;
+
+ public:
+ itkTypeMacro( VertexParentBase, TFilter );
+
+ public:
+ virtual TOutputValue Evaluate( const TVertex& v, const TVertex& p ) const = 0;
+
+ protected:
+ VertexParentBase( )
+ : Superclass( )
+ { }
+ virtual ~VertexParentBase( ) { }
+
+ private:
+ VertexParentBase( const Self& other );
+ Self& operator=( const Self& other );
+ };
+
+ } // ecapseman
+
+ } // ecapseman
+
+} // ecapseman
+
+#endif // __fpa__Base__Functors__VertexParentBase__h__
+
+// eof - $RCSfile$
--- /dev/null
+// =========================================================================
+// @author Leonardo Florez Valencia
+// @email florez-l@javeriana.edu.co
+// =========================================================================
+
+#ifndef __fpa__Base__MinimumSpanningTree__h__
+#define __fpa__Base__MinimumSpanningTree__h__
+
+#include <vector>
+
+namespace fpa
+{
+ namespace Base
+ {
+ /**
+ */
+ template< class _TVertex, class _Superclass >
+ class MinimumSpanningTree
+ : public _Superclass
+ {
+ public:
+ typedef MinimumSpanningTree Self;
+ typedef _Superclass Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ typedef _TVertex TVertex;
+ typedef std::pair< TVertex, bool > TCollision;
+ typedef std::vector< TCollision > TCollisionsRow;
+ typedef std::vector< TCollisionsRow > TCollisions;
+ typedef std::vector< TVertex > TVertices;
+
+ protected:
+ typedef std::vector< unsigned long > _TRow;
+ typedef std::vector< _TRow > _TMatrix;
+
+ public:
+ itkTypeMacro( fpa::Base::MinimumSpanningTree, _Superclass );
+
+ 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 TVertices GetPath( const TVertex& a ) const;
+ virtual TVertices GetPath( const TVertex& a, const TVertex& b ) const;
+
+ protected:
+ MinimumSpanningTree( );
+ virtual ~MinimumSpanningTree( );
+
+ private:
+ 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
+// =========================================================================
+// @author Leonardo Florez Valencia
+// @email florez-l@javeriana.edu.co
+// =========================================================================
+
+#ifndef __fpa__Base__MinimumSpanningTree__hxx__
+#define __fpa__Base__MinimumSpanningTree__hxx__
+
+// -------------------------------------------------------------------------
+template< class _TVertex, class _Superclass >
+const typename fpa::Base::MinimumSpanningTree< _TVertex, _Superclass >::
+TCollisions& fpa::Base::MinimumSpanningTree< _TVertex, _Superclass >::
+GetCollisions( ) const
+{
+ return( this->m_Collisions );
+}
+
+// -------------------------------------------------------------------------
+template< class _TVertex, class _Superclass >
+void fpa::Base::MinimumSpanningTree< _TVertex, _Superclass >::
+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 _Superclass >
+void fpa::Base::MinimumSpanningTree< _TVertex, _Superclass >::
+ClearSeeds( )
+{
+ this->m_Seeds.clear( );
+ if( this->m_DataObject != NULL )
+ this->m_DataObject->Modified( );
+}
+
+// -------------------------------------------------------------------------
+template< class _TVertex, class _Superclass >
+void fpa::Base::MinimumSpanningTree< _TVertex, _Superclass >::
+AddSeed( const _TVertex& seed )
+{
+ this->m_Seeds.push_back( seed );
+ if( this->m_DataObject != NULL )
+ this->m_DataObject->Modified( );
+}
+
+// -------------------------------------------------------------------------
+template< class _TVertex, class _Superclass >
+typename fpa::Base::MinimumSpanningTree< _TVertex, _Superclass >::
+TVertices fpa::Base::MinimumSpanningTree< _TVertex, _Superclass >::
+GetPath( const _TVertex& a ) const
+{
+ TVertices 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 );
+ return( vertices );
+}
+
+// -------------------------------------------------------------------------
+template< class _TVertex, class _Superclass >
+typename fpa::Base::MinimumSpanningTree< _TVertex, _Superclass >::
+TVertices fpa::Base::MinimumSpanningTree< _TVertex, _Superclass >::
+GetPath( const _TVertex& a, const _TVertex& b ) const
+{
+ static const unsigned long _inf =
+ std::numeric_limits< unsigned long >::max( );
+
+ TVertices vertices;
+ TVertices pa = this->GetPath( a );
+ TVertices pb = this->GetPath( b );
+ if( pa.size( ) > 0 && pb.size( ) > 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[ pa.size( ) - 1 ] )
+ ia = i;
+ if( this->m_Seeds[ i ] == pb[ pb.size( ) - 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
+ vertices = this->GetPath(
+ a, this->m_Collisions[ fpath[ 0 ] ][ fpath[ 1 ] ].first
+ );
+
+ // Intermediary paths
+ for( unsigned int i = 1; i < N - 1; ++i )
+ {
+ TVertices ipath =
+ this->GetPath(
+ this->m_Collisions[ fpath[ i ] ][ fpath[ i - 1 ] ].first,
+ this->m_Collisions[ fpath[ i ] ][ fpath[ i + 1 ] ].first
+ );
+ for( long id = 0; id < ipath.size( ); ++id )
+ vertices.push_back( ipath[ id ] );
+
+ } // rof
+
+ // Final path: from last collision to end point
+ TVertices lpath =
+ this->GetPath(
+ this->m_Collisions[ fpath[ N - 1 ] ][ fpath[ N - 2 ] ].first, b
+ );
+ for( long id = 0; id < lpath.size( ); ++id )
+ vertices.push_back( lpath[ id ] );
+
+ } // fi
+ }
+ else
+ {
+ // Ignore common part: find common ancestor
+ long aIt = pa.size( ) - 1;
+ long bIt = pb.size( ) - 1;
+ bool cont = true;
+ while( aIt >= 0 && bIt >= 0 && cont )
+ {
+ cont = ( pa[ aIt ] == pb[ bIt ] );
+ aIt--;
+ bIt--;
+
+ } // elihw
+ aIt++;
+ bIt++;
+
+ // Glue both parts
+ for( long cIt = 0; cIt <= aIt; ++cIt )
+ vertices.push_back( pa[ cIt ] );
+ for( ; bIt >= 0; --bIt )
+ vertices.push_back( pb[ bIt ] );
+
+ } // fi
+
+ } // fi
+ return( vertices );
+}
+
+// -------------------------------------------------------------------------
+template< class _TVertex, class _Superclass >
+fpa::Base::MinimumSpanningTree< _TVertex, _Superclass >::
+MinimumSpanningTree( )
+ : Superclass( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TVertex, class _Superclass >
+fpa::Base::MinimumSpanningTree< _TVertex, _Superclass >::
+~MinimumSpanningTree( )
+{
+}
+
+#endif // __fpa__Base__MinimumSpanningTree__hxx__
+
+// eof - $RCSfile$
--- /dev/null
+// =========================================================================
+// @author Leonardo Florez Valencia
+// @email florez-l@javeriana.edu.co
+// =========================================================================
+
+#ifndef __fpa__Image__Dijkstra__h__
+#define __fpa__Image__Dijkstra__h__
+
+#include <fpa/Base/Dijkstra.h>
+#include <fpa/Base/SeedsInterface.h>
+#include <fpa/Image/MarksInterface.h>
+#include <fpa/Image/Filter.h>
+#include <fpa/Image/MinimumSpanningTree.h>
+#include <fpa/Image/Functors/VertexParentBase.h>
+
+namespace fpa
+{
+ namespace Image
+ {
+ /**
+ */
+ template< class _TInputImage, class _TOutputImage >
+ class Dijkstra
+ : public fpa::Base::Dijkstra< fpa::Image::Filter< _TInputImage, _TOutputImage >, fpa::Image::MarksInterface< _TInputImage::ImageDimension >, fpa::Base::SeedsInterface< typename _TInputImage::IndexType, typename _TInputImage::IndexType::LexicographicCompare >, fpa::Image::MinimumSpanningTree< _TInputImage::ImageDimension > >
+ {
+ public:
+ // Interfaces
+ typedef fpa::Image::Filter< _TInputImage, _TOutputImage > TFilter;
+ typedef fpa::Image::MarksInterface< _TInputImage::ImageDimension > TMarksInterface;
+ typedef fpa::Base::SeedsInterface< typename _TInputImage::IndexType, typename _TInputImage::IndexType::LexicographicCompare > TSeedsInterface;
+ typedef fpa::Image::MinimumSpanningTree< _TInputImage::ImageDimension > TMST;
+
+ // Smart pointers
+ typedef Dijkstra Self;
+ typedef fpa::Base::Dijkstra< TFilter, TMarksInterface, TSeedsInterface, TMST > Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+ typedef typename TFilter::TInputImage TInputImage;
+ typedef typename TFilter::TOutputValue TOutputValue;
+ typedef typename TFilter::TVertex TVertex;
+
+ public:
+ itkNewMacro( Self );
+ itkTypeMacro( fpa::Image::Dijkstra, fpa::Base::Dijkstra );
+
+ protected:
+ Dijkstra( ) : Superclass( ) { }
+ virtual ~Dijkstra( ) { }
+
+ virtual void _ConfigureOutputs( const TOutputValue& init_value ) override
+ {
+ this->Superclass::_ConfigureOutputs( init_value );
+
+ typename TVertex::OffsetType o;
+ o.Fill( 0 );
+ const TInputImage* input = this->GetInput( );
+ TMST* mst = this->GetMinimumSpanningTree( );
+ mst->CopyInformation( input );
+ mst->SetBufferedRegion( input->GetRequestedRegion( ) );
+ mst->Allocate( );
+ mst->FillBuffer( o );
+ }
+
+ virtual void GenerateData( ) override
+ {
+ // Configure functors with input image
+ typedef typename TFilter::TInputImage _TInputIage;
+ typedef typename TFilter::TOutputValue _TOutputValue;
+ typedef fpa::Image::Functors::VertexParentBase< _TInputImage, _TOutputValue > _TVFunc;
+ _TVFunc* vfunc =
+ dynamic_cast< _TVFunc* >( this->m_VertexFunctor.GetPointer( ) );
+ if( vfunc != NULL )
+ vfunc->SetImage( this->GetInput( ) );
+
+ // Ok, continue
+ this->Superclass::GenerateData( );
+ }
+
+ private:
+ Dijkstra( const Self& other );
+ Self& operator=( const Self& other );
+ };
+
+ } // ecapseman
+
+} // ecapseman
+
+#endif // __fpa__Image__Dijkstra__h__
+
+// eof - $RCSfile$
--- /dev/null
+// =========================================================================
+// @author Leonardo Florez Valencia
+// @email florez-l@javeriana.edu.co
+// =========================================================================
+
+#ifndef __fpa__Image__Functors__VertexIdentity__h__
+#define __fpa__Image__Functors__VertexIdentity__h__
+
+#include <fpa/Image/Functors/VertexParentBase.h>
+
+namespace fpa
+{
+ namespace Image
+ {
+ namespace Functors
+ {
+ /**
+ */
+ template< class _TInputImage, class _TOutputValue >
+ class VertexIdentity
+ : public fpa::Image::Functors::VertexParentBase< _TInputImage, _TOutputValue >
+ {
+ public:
+ typedef _TInputImage TInputImage;
+ typedef _TOutputValue TOutputValue;
+ typedef VertexIdentity Self;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+ typedef fpa::Image::Functors::VertexParentBase< TInputImage, TOutputValue > Superclass;
+
+ typedef typename Superclass::TVertex TVertex;
+
+ public:
+ itkNewMacro( Self );
+ itkTypeMacro(
+ fpa::Image::Functors::VertexIdentity,
+ fpa::Image::Functors::VertexParentBase
+ );
+
+ public:
+ virtual TOutputValue Evaluate(
+ const TVertex& a, const TVertex& p
+ ) const override
+ {
+ return( TOutputValue( this->m_Image->GetPixel( a ) ) );
+ }
+
+ protected:
+ VertexIdentity( )
+ : Superclass( )
+ { }
+ virtual ~VertexIdentity( ) { }
+
+ private:
+ VertexIdentity( const Self& other );
+ Self& operator=( const Self& other );
+ };
+
+ } // ecapseman
+
+ } // ecapseman
+
+} // ecapseman
+
+#endif // __fpa__Image__Functors__VertexIdentity__h__
+
+// eof - $RCSfile$
--- /dev/null
+// =========================================================================
+// @author Leonardo Florez Valencia
+// @email florez-l@javeriana.edu.co
+// =========================================================================
+
+#ifndef __fpa__Image__Functors__VertexParentBase__h__
+#define __fpa__Image__Functors__VertexParentBase__h__
+
+#include <fpa/Base/Functors/VertexParentBase.h>
+
+namespace fpa
+{
+ namespace Image
+ {
+ namespace Functors
+ {
+ /**
+ */
+ template< class _TInputImage, class _TOutputValue >
+ class VertexParentBase
+ : public fpa::Base::Functors::VertexParentBase< typename _TInputImage::IndexType, _TOutputValue >
+ {
+ public:
+ typedef _TInputImage TInputImage;
+ typedef _TOutputValue TOutputValue;
+ typedef typename TInputImage::IndexType TVertex;
+ typedef VertexParentBase Self;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+ typedef fpa::Base::Functors::VertexParentBase< TVertex, TOutputValue > Superclass;
+
+ public:
+ itkTypeMacro(
+ fpa::Image::Functors::VertexParentBase,
+ fpa::Base::Functors::VertexParentBase
+ );
+
+ itkGetConstObjectMacro( Image, TInputImage );
+ itkSetConstObjectMacro( Image, TInputImage );
+
+ protected:
+ VertexParentBase( )
+ : Superclass( )
+ { }
+ virtual ~VertexParentBase( ) { }
+
+ private:
+ VertexParentBase( const Self& other );
+ Self& operator=( const Self& other );
+
+ protected:
+ typename TInputImage::ConstPointer m_Image;
+ };
+
+ } // ecapseman
+
+ } // ecapseman
+
+} // ecapseman
+
+#endif // __fpa__Image__Functors__VertexParentBase__h__
+
+// eof - $RCSfile$
--- /dev/null
+// =========================================================================
+// @author Leonardo Florez Valencia
+// @email florez-l@javeriana.edu.co
+// =========================================================================
+
+#ifndef __fpa__Image__MinimumSpanningTree__h__
+#define __fpa__Image__MinimumSpanningTree__h__
+
+#include <fpa/Base/MinimumSpanningTree.h>
+#include <itkImage.h>
+
+namespace fpa
+{
+ namespace Image
+ {
+ /**
+ */
+ template< unsigned int _VDim >
+ class MinimumSpanningTree
+ : public fpa::Base::MinimumSpanningTree< itk::Index< _VDim >, itk::Image< itk::Offset< _VDim >, _VDim > >
+ {
+ public:
+ typedef itk::Index< _VDim > TVertex;
+ typedef itk::Image< itk::Offset< _VDim >, _VDim > TBaseImage;
+
+ typedef MinimumSpanningTree Self;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+ typedef fpa::Base::MinimumSpanningTree< TVertex, TBaseImage > Superclass;
+
+ typedef typename Superclass::TCollision TCollision;
+ typedef typename Superclass::TCollisionsRow TCollisionsRow;
+ typedef typename Superclass::TCollisions TCollisions;
+ typedef typename Superclass::TVertices TVertices;
+
+ public:
+ itkNewMacro( Self );
+ itkTypeMacro(
+ fpa::Image::MinimumSpanningTree,
+ fpa::Base::MinimumSpanningTree
+ );
+
+ public:
+ virtual TVertex GetParent( const TVertex& v ) const override
+ {
+ return( v + this->GetPixel( v ) );
+ }
+ virtual void SetParent( const TVertex& v, const TVertex& p ) override
+ {
+ this->SetPixel( v, p - v );
+ }
+
+ protected:
+ MinimumSpanningTree( )
+ : Superclass( )
+ { }
+ virtual ~MinimumSpanningTree( )
+ { }
+
+ private:
+ MinimumSpanningTree( const Self& other );
+ Self& operator=( const Self& other );
+ };
+
+ } // ecapseman
+
+} // ecapseman
+
+#endif // __fpa__Image__MinimumSpanningTree__h__
+
+// eof - $RCSfile$