#ifndef __fpa__Image__SkeletonFilter__h__
#define __fpa__Image__SkeletonFilter__h__
+#include <vector>
#include <fpa/Image/Dijkstra.h>
#include <cpExtensions/DataStructures/Skeleton.h>
#include <itkSimpleDataObjectDecorator.h>
-#include <functional>
-#include <map>
-#include <set>
-
-/* TODO
- #include <itkProcessObject.h>
- #include <fpa/Image/MinimumSpanningTree.h>
-*/
namespace fpa
{
typedef typename Superclass::TMST TMST;
typedef typename TImage::IndexType TIndex;
-
- typedef
- itk::Functor::IndexLexicographicCompare< Self::Dimension > TIndexCmp;
- typedef std::set< TIndex, TIndexCmp > TIndicesData;
- typedef itk::SimpleDataObjectDecorator< TIndicesData > TIndices;
- typedef itk::Image< unsigned char, Self::Dimension > TMarks;
- typedef
- cpExtensions::DataStructures::Skeleton< Self::Dimension > TSkeleton;
+ 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::multimap< double, TIndex, std::greater< double > > _TSkeletonQueue;
public:
itkNewMacro( Self );
itkTypeMacro( fpa::Image::SkeletonFilter, itk::Object );
public:
- TIndices* GetEndPoints( );
- TIndices* GetBifurcations( );
TSkeleton* GetSkeleton( );
TMarks* GetMarks( );
virtual void _UpdateResult( const _TQueueNode& n ) fpa_OVERRIDE;
virtual void _AfterGenerateData( ) fpa_OVERRIDE;
- void _EndPoints( );
- void _Skeleton( );
-
- /* TODO
- virtual void GenerateData( ) fpa_OVERRIDE;
- const TDistanceMap* dmap,
- const TCostMap* cmap,
- const TMST* mst,
- TIndicesData& end_points
- );
- const TDistanceMap* dmap,
- const TCostMap* cmap,
- const TMST* mst,
- const TIndicesData& end_points,
- TIndicesData& bifurcations,
- TSkeleton* skeleton
- );
- */
+ void _EndPoints( std::vector< TIndex >& end_points );
+ void _Skeleton( const std::vector< TIndex >& end_points );
private:
// Purposely not defined
protected:
_TSkeletonQueue m_SkeletonQueue;
- unsigned long m_EndPointsIdx;
- unsigned long m_BifurcationsIdx;
unsigned long m_SkeletonIdx;
unsigned long m_MarksIdx;
};
#ifndef __fpa__Image__SkeletonFilter__hxx__
#define __fpa__Image__SkeletonFilter__hxx__
-/* TODO
- #include <queue>
- #include <itkImageRegionConstIteratorWithIndex.h>
-*/
#include <fpa/Base/Functors/Inverse.h>
#include <fpa/Image/Functors/SimpleNeighborhood.h>
#include <itkImageRegionIteratorWithIndex.h>
}
fpa_Image_SkeletonFilter_OutputMacro( Skeleton, TSkeleton );
-fpa_Image_SkeletonFilter_OutputMacro( EndPoints, TIndices );
-fpa_Image_SkeletonFilter_OutputMacro( Bifurcations, TIndices );
fpa_Image_SkeletonFilter_OutputMacro( Marks, TMarks );
// -------------------------------------------------------------------------
typedef fpa::Base::Functors::Inverse< _TPixel, _TPixel > _TInvFunc;
unsigned int nOutputs = this->GetNumberOfRequiredOutputs( );
- this->SetNumberOfRequiredOutputs( nOutputs + 4 );
- this->m_EndPointsIdx = nOutputs;
- this->m_BifurcationsIdx = nOutputs + 1;
- this->m_SkeletonIdx = nOutputs + 2;
- this->m_MarksIdx = nOutputs + 3;
-
- typename TIndices::Pointer end_points = TIndices::New( );
- typename TIndices::Pointer bifurcations = TIndices::New( );
+ 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_EndPointsIdx, end_points.GetPointer( ) );
- this->SetNthOutput( this->m_BifurcationsIdx, bifurcations.GetPointer( ) );
this->SetNthOutput( this->m_SkeletonIdx, skeleton.GetPointer( ) );
this->SetNthOutput( this->m_MarksIdx, marks.GetPointer( ) );
typedef typename _TSkeletonQueue::value_type _TSkeletonQueueValue;
this->Superclass::_UpdateResult( n );
-
- // Update skeleton candidates
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 ) );
_AfterGenerateData( )
{
this->Superclass::_AfterGenerateData( );
- this->_EndPoints( );
- this->_Skeleton( );
+
+ std::vector< TIndex > end_points;
+ this->_EndPoints( end_points );
+ this->_Skeleton( end_points );
}
// -------------------------------------------------------------------------
template< class _TImage >
void fpa::Image::SkeletonFilter< _TImage >::
-_EndPoints( )
+_EndPoints( std::vector< TIndex >& end_points )
{
typedef itk::ImageRegionIteratorWithIndex< TMarks > _TMarksIt;
- /* TODO
- typedef itk::ImageRegionConstIteratorWithIndex< _TImage > _TDistMapIt;
- typedef itk::ImageRegionConstIteratorWithIndex< _TOutput > _TOutputIt;
- typedef std::multimap< double, TIndex, std::greater< double > > _TQueue;
- typedef typename _TQueue::value_type _TQueueValue;
- */
-
static const double _eps = std::sqrt( double( Self::Dimension + 1 ) );
auto input = this->GetInput( );
auto mst = this->GetMinimumSpanningTree( );
auto marks = this->GetMarks( );
- auto& end_points = this->GetEndPoints( )->Get( );
// Some values
marks->SetLargestPossibleRegion( input->GetLargestPossibleRegion( ) );
marks->SetDirection( input->GetDirection( ) );
marks->Allocate( );
marks->FillBuffer( 0 );
+ auto spac = marks->GetSpacing( );
// BFS from maximum queue
auto region = input->GetRequestedRegion( );
if( m != 0 )
continue;
marks->SetPixel( n.second, 1 );
- end_points.insert( n.second );
+ end_points.push_back( n.second );
// Mark path
- auto spac = marks->GetSpacing( );
- typename TMST::TPath::Pointer path;
- mst->GetPath( path, n.second );
- for( unsigned long i = 0; i < path->GetSize( ); ++i )
+ TIndex it = n.second;
+ TIndex p = mst->GetParent( it );
+ while( it != p )
{
- TIndex idx = path->GetVertex( i );
- 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
-
- } // rof
+ // TODO: tags->SetPixel( it, tags->GetPixel( it ) + 1 );
+ it = p;
+ p = mst->GetParent( it );
+
+ } // elihw
+ // TODO: tags->SetPixel( it, tags->GetPixel( it ) + 1 );
+
+ /* TODO: use mst directly rather than computing paths
+ typename TMST::TPath::Pointer path;
+ mst->GetPath( path, n.second );
+ for( unsigned long i = 0; i < path->GetSize( ); ++i )
+ {
+ TIndex idx = path->GetVertex( i );
+ 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
+
+ } // rof
+ */
} // elihw
}
// -------------------------------------------------------------------------
template< class _TImage >
void fpa::Image::SkeletonFilter< _TImage >::
-_Skeleton( )
+_Skeleton( const std::vector< TIndex >& end_points )
{
+ typedef itk::Image< unsigned long, Self::Dimension > _TTagsImage;
+ typedef typename TMST::TPath _TPath;
+
auto mst = this->GetMinimumSpanningTree( );
auto skeleton = this->GetSkeleton( );
- auto& end_points = this->GetEndPoints( )->Get( );
- typedef typename TMST::TPath _TPath;
+ // 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 )
{
- typename _TPath::Pointer path;
- mst->GetPath( path, *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
}
-// -------------------------------------------------------------------------
-/* TODO
- template< class _TImage >
- void fpa::Image::SkeletonFilter< _TImage >::
- _Skeleton(
- const TDistanceMap* dmap,
- const TCostMap* cmap,
- const TMST* mst,
- const TIndicesData& end_points,
- TIndicesData& bifurcations,
- TSkeleton* skeleton
- )
- {
- }
-*/
-
#endif // __fpa__Image__SkeletonFilter__hxx__
// eof - $RCSfile$
git://git.creatis.insa-lyon.fr / FrontAlgorithms.git / commitdiff