{
this->Superclass::_AfterGenerateData( );
+ _TAdjacencies A;
std::vector< TIndex > end_points;
- this->_EndPoints( end_points );
- this->_Skeleton( 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 )
+_EndPoints( std::vector< TIndex >& end_points, _TAdjacencies& A )
{
- typedef itk::ImageRegionIteratorWithIndex< TMarks > _TMarksIt;
-
- static const double _eps = std::sqrt( double( Self::Dimension + 1 ) );
- auto input = this->GetInput( );
- auto mst = this->GetMinimumSpanningTree( );
auto marks = this->GetMarks( );
+ auto mst = this->GetMinimumSpanningTree( );
+ auto spac = marks->GetSpacing( );
// Some values
- marks->SetLargestPossibleRegion( input->GetLargestPossibleRegion( ) );
- marks->SetRequestedRegion( input->GetRequestedRegion( ) );
- marks->SetBufferedRegion( input->GetBufferedRegion( ) );
- marks->SetSpacing( input->GetSpacing( ) );
- marks->SetOrigin( input->GetOrigin( ) );
- marks->SetDirection( input->GetDirection( ) );
+ 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 );
- auto spac = marks->GetSpacing( );
// BFS from maximum queue
- auto region = input->GetRequestedRegion( );
while( this->m_SkeletonQueue.size( ) > 0 )
{
// Get node
continue;
marks->SetPixel( n.second, 1 );
end_points.push_back( n.second );
+ std::cout << this->m_SkeletonQueue.size( ) << std::endl;
// Mark path
TIndex it = n.second;
TIndex p = mst->GetParent( it );
while( it != p )
{
- // TODO: tags->SetPixel( it, tags->GetPixel( it ) + 1 );
+ this->_MarkSphere( it );
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
- */
+ this->_MarkSphere( it );
+ A[ n.second ] = it;
} // elihw
}
// -------------------------------------------------------------------------
template< class _TImage >
void fpa::Image::SkeletonFilter< _TImage >::
-_Skeleton( const std::vector< TIndex >& end_points )
+_Skeleton( const std::vector< TIndex >& end_points, _TAdjacencies& A )
{
- typedef itk::Image< unsigned long, Self::Dimension > _TTagsImage;
- typedef typename TMST::TPath _TPath;
+ std::cout << end_points.size( ) << " " << A.size( ) << std::endl;
+
+ /* TODO
+ 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
+ */
+}
- 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 )
+// -------------------------------------------------------------------------
+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 )
{
- TIndex it = *eIt;
- TIndex p = mst->GetParent( it );
- while( it != p )
- {
- tags->SetPixel( it, tags->GetPixel( it ) + 1 );
- it = p;
- p = mst->GetParent( it );
+ long off = long( std::ceil( r / double( spac[ d ] ) ) );
+ if( off < 3 )
+ off = 3;
+ i0[ d ] = idx[ d ] - off;
+ i1[ d ] = idx[ d ] + off;
- } // elihw
- tags->SetPixel( it, tags->GetPixel( it ) + 1 );
+ 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
- // 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 );
+ typename _TImage::SizeType size;
+ for( unsigned int d = 0; d < Self::Dimension; ++d )
+ size[ d ] = i1[ d ] - i0[ d ] + 1;
- } // elihw
+ typename _TImage::RegionType neighRegion;
+ neighRegion.SetIndex( i0 );
+ neighRegion.SetSize( size );
- // Finally, add last branch
- path->AddVertex( it );
- skeleton->AddBranch( path );
+ _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
}
git://git.creatis.insa-lyon.fr / FrontAlgorithms.git / blobdiff