--- /dev/null
+#ifndef __fpa__Image__SkeletonFilter__h__
+#define __fpa__Image__SkeletonFilter__h__
+
+#include <itkProcessObject.h>
+#include <fpa/Image/MinimumSpanningTree.h>
+#include <cpExtensions/DataStructures/Skeleton.h>
+#include <set>
+
+namespace fpa
+{
+ namespace Image
+ {
+ /**
+ */
+ template< class _TDistanceMap, class _TCostMap >
+ class SkeletonFilter
+ : public itk::ProcessObject
+ {
+ public:
+ typedef SkeletonFilter Self;
+ typedef itk::ProcessObject Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
+
+#ifdef ITK_USE_CONCEPT_CHECKING
+ itkConceptMacro(
+ SameDimension,
+ ( itk::Concept::SameDimension< _TDistanceMap::ImageDimension, _TCostMap::ImageDimension > )
+ );
+#endif
+ typedef _TDistanceMap TDistanceMap;
+ typedef _TCostMap TCostMap;
+ typedef typename TCostMap::IndexType TIndex;
+ typedef MinimumSpanningTree< TCostMap::ImageDimension > TMST;
+ typedef
+ cpExtensions::DataStructures::Skeleton< TCostMap::ImageDimension >
+ TSkeleton;
+
+ typedef
+ itk::Functor::IndexLexicographicCompare< _TCostMap::ImageDimension >
+ TIndexCompare;
+ typedef std::set< TIndex, TIndexCompare > TIndicesData;
+ typedef itk::SimpleDataObjectDecorator< TIndicesData > TIndices;
+
+ public:
+ itkNewMacro( Self );
+ itkTypeMacro( fpa::Image::SkeletonFilter, itk::Object );
+
+ public:
+ _TDistanceMap* GetDistanceMap( );
+ const _TDistanceMap* GetDistanceMap( ) const;
+ void SetDistanceMap( _TDistanceMap* dmap );
+
+ _TCostMap* GetCostMap( );
+ const _TCostMap* GetCostMap( ) const;
+ void SetCostMap( _TCostMap* cmap );
+
+ TMST* GetMinimumSpanningTree( );
+ const TMST* GetMinimumSpanningTree( ) const;
+ void SetMinimumSpanningTree( TMST* mst );
+
+ TIndices* GetEndPoints( );
+ const TIndices* GetEndPoints( ) const;
+
+ TIndices* GetBifurcations( );
+ const TIndices* GetBifurcations( ) const;
+
+ TSkeleton* GetSkeleton( );
+ const TSkeleton* GetSkeleton( ) const;
+
+ protected:
+ SkeletonFilter( );
+ virtual ~SkeletonFilter( );
+
+ virtual void GenerateData( ) fpa_OVERRIDE;
+
+ void _EndPoints(
+ const TDistanceMap* dmap,
+ const TCostMap* cmap,
+ const TMST* mst,
+ TIndices* end_points
+ );
+
+ void _Skeleton(
+ const TDistanceMap* dmap,
+ const TCostMap* cmap,
+ const TMST* mst,
+ TIndices* end_points
+ TIndices* bifurcations,
+ TSkeleton* skeleton
+ );
+
+ private:
+ // Purposely not defined
+ SkeletonFilter( const Self& other );
+ Self& operator=( const Self& other );
+ };
+
+ } // 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__
+
+// -------------------------------------------------------------------------
+#define __fpa__Image__SkeletonFilter__InputMacro( input_name, input_type, input_id ) \
+ template< class _TDistanceMap, class _TCostMap > \
+ typename fpa::Image::SkeletonFilter< _TDistanceMap, _TCostMap >:: \
+ input_type* fpa::Image::SkeletonFilter< _TDistanceMap, _TCostMap >:: \
+ Get##input_name( ) \
+ { \
+ return( \
+ dynamic_cast< input_type* >( \
+ this->Superclass::GetInput( input_id ) \
+ ) \
+ ); \
+ } \
+ template< class _TDistanceMap, class _TCostMap > \
+ const \
+ typename fpa::Image::SkeletonFilter< _TDistanceMap, _TCostMap >:: \
+ input_type* \
+ fpa::Image::SkeletonFilter< _TDistanceMap, _TCostMap >:: \
+ Get##input_name( ) const \
+ { \
+ return( \
+ dynamic_cast< const input_type* >( \
+ this->Superclass::GetInput( input_id ) \
+ ) \
+ ); \
+ } \
+ template< class _TDistanceMap, class _TCostMap > \
+ void fpa::Image::SkeletonFilter< _TDistanceMap, _TCostMap >:: \
+ Set##input_name( input_type* input ) \
+ { \
+ this->Superclass::SetInput( input_id, input ); \
+ }
+
+__fpa__Image__SkeletonFilter__InputMacro( DistanceMap, TDistanceMap, 0 );
+__fpa__Image__SkeletonFilter__InputMacro( CostMap, TCostMap, 1 );
+__fpa__Image__SkeletonFilter__InputMacro( MinimumSpanningTree, TMST, 2 );
+
+// -------------------------------------------------------------------------
+#define __fpa__Image__SkeletonFilter__OutputMacro( output_name, output_type, output_id ) \
+ template< class _TDistanceMap, class _TCostMap > \
+ typename fpa::Image::SkeletonFilter< _TDistanceMap, _TCostMap >:: \
+ output_type* fpa::Image::SkeletonFilter< _TDistanceMap, _TCostMap >:: \
+ Get##output_name( ) \
+ { \
+ return( \
+ dynamic_cast< output_type* >( \
+ this->Superclass::GetOutput( output_id ) \
+ ) \
+ ); \
+ } \
+ template< class _TDistanceMap, class _TCostMap > \
+ const typename \
+ fpa::Image::SkeletonFilter< _TDistanceMap, _TCostMap >:: \
+ output_type* fpa::Image::SkeletonFilter< _TDistanceMap, _TCostMap >:: \
+ Get##output_name( ) const \
+ { \
+ return( \
+ dynamic_cast< const output_type* >( \
+ this->Superclass::GetOutput( output_id ) \
+ ) \
+ ); \
+ }
+
+__fpa__Image__SkeletonFilter__OutputMacro( Skeleton, TSkeleton, 0 );
+__fpa__Image__SkeletonFilter__OutputMacro( EndPoints, TIndices, 1 );
+__fpa__Image__SkeletonFilter__OutputMacro( Bifurcations, TIndices, 2 );
+
+// -------------------------------------------------------------------------
+template< class _TDistanceMap, class _TCostMap >
+fpa::Image::SkeletonFilter< _TDistanceMap, _TCostMap >::
+SkeletonFilter( )
+ : Superclass( )
+{
+ this->SetNumbeOfRequiredInputs( 3 );
+ this->SetNumbeOfRequiredOutputs( 3 );
+
+ typename TIndices::Pointer end_points = TIndices::New( );
+ typename TIndices::Pointer bifurcations = TIndices::New( );
+ typename TSkeleton::Pointer skeleton = TSkeleton::New( );
+ this->SetNthOutput( 0, skeleton );
+ this->SetNthOutput( 1, end_points );
+ this->SetNthOutput( 2, bifurcations );
+}
+
+// -------------------------------------------------------------------------
+template< class _TDistanceMap, class _TCostMap >
+fpa::Image::SkeletonFilter< _TDistanceMap, _TCostMap >::
+~SkeletonFilter( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TDistanceMap, class _TCostMap >
+void fpa::Image::SkeletonFilter< _TDistanceMap, _TCostMap >::
+GenerateData( )
+{
+ // 0. I/O objects
+ const TDistanceMap* dmap = this->GetDistanceMap( );
+ const TCostMap* cmap = this->GetCostMap( );
+ const TMST* mst = this->GetMinimumSpanningTree( );
+ TIndices* end_points = this->GetEndPoints( );
+ TIndices* bifurcations = this->GetBifurcations( );
+ TSkeleton* skeleton = this->GetSkeleton( );
+
+ // 1. Check input correspondance
+ // TODO
+
+ // 2. Detect end-points
+ this->_EndPoints( dmap, cmap, mst, end_points );
+
+ // 3. Build skeleton and keep track of bifurcations
+ this->_Skeleton( dmap, cmap, mst, end_points, bifurcations, skeleton );
+}
+
+// -------------------------------------------------------------------------
+template< class _TDistanceMap, class _TCostMap >
+void fpa::Image::SkeletonFilter< _TDistanceMap, _TCostMap >::
+_EndPoints(
+ const TDistanceMap* dmap,
+ const TCostMap* cmap,
+ const TMST* mst,
+ TIndices* end_points
+ )
+{
+ 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 _TMarks::Pointer marks = _TMarks::New( );
+ marks->SetLargestPossibleRegion( dmap->GetLargestPossibleRegion( ) );
+ marks->SetRequestedRegion( dmap->GetRequestedRegion( ) );
+ marks->SetBufferedRegion( dmap->GetBufferedRegion( ) );
+ marks->SetSpacing( dmap->GetSpacing( ) );
+ marks->SetOrigin( dmap->GetOrigin( ) );
+ marks->SetDirection( dmap->GetDirection( ) );
+ marks->Allocate( );
+ marks->FillBuffer( 0 );
+
+ // Create queue
+ _TQueue queue;
+ _TDistMapIt dIt( dmap, dmap->GetRequestedRegion( ) );
+ _TCostMapIt cIt( cmap, cmap->GetRequestedRegion( ) );
+ 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 and update end-points
+ m |= 0x01;
+ marks->SetPixel( n.second, m );
+ end_points->Get( ).insert( n.second );
+
+ // Get path
+ typename TMST::TPath::Pointer path;
+ 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;
+ cmap->TransformIndexToPhysicalPoint( idx, cnt );
+ double d = double( dmap->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;
+ cmap->TransformIndexToPhysicalPoint( w, p );
+ if( cnt.EuclideanDistanceTo( p ) <= d )
+ q.push( w );
+
+ } // fi
+
+ } // rof
+
+ } // rof
+
+ } // elihw
+
+ } // rof
+
+ } // elihw
+}
+
+// -------------------------------------------------------------------------
+template< class _TDistanceMap, class _TCostMap >
+void fpa::Image::SkeletonFilter< _TDistanceMap, _TCostMap >::
+_Skeleton(
+ const TDistanceMap* dmap,
+ const TCostMap* cmap,
+ const TMST* mst,
+ TIndices* end_points
+ TIndices* bifurcations,
+ TSkeleton* skeleton
+ )
+{
+}
+
+#endif // __fpa__Image__SkeletonFilter__hxx__
+
+// eof - $RCSfile$
git://git.creatis.insa-lyon.fr / FrontAlgorithms.git / commitdiff