...

[FrontAlgorithms.git] / lib / fpa / Image / SkeletonFilter.hxx

diff --git a/lib/fpa/Image/SkeletonFilter.hxx b/lib/fpa/Image/SkeletonFilter.hxx
index 50a035453cfd8660aec6cb4f6a2d19c485d0043e..224829b93f5a37119a5ceea5f3421b26814b9f43 100644 (file)
--- a/lib/fpa/Image/SkeletonFilter.hxx
+++ b/lib/fpa/Image/SkeletonFilter.hxx
@@ -174,24 +174,76 @@ GenerateData( )
     mst, dijkstra->GetSkeletonQueue( )
     );
 
-  // Create branches
-  TSkeleton* sk = this->GetOutput( );
+  // Compute symbolic branches
+  typedef std::map< TIndex, TIndex, typename TIndex::LexicographicCompare > _TTags;
+  _TTags tags, branches;
   typename std::vector< TIndex >::const_iterator eIt = end_points.begin( );
-  std::map< TIndex, unsigned long, typename TIndex::LexicographicCompare > tags;
-  unsigned long t = 1;
   for( ; eIt != end_points.end( ); ++eIt )
   {
     // Tag path
     TIndex it = *eIt;
     TIndex p = mst->GetParent( it );
-    while( it != p )
+    typename _TTags::iterator bIt = tags.end( );
+    while( it != p && bIt == tags.end( ) )
     {
-      tags[ it ] = t;
+      typename _TTags::iterator tIt = tags.find( it );
+      if( tIt != tags.end( ) )
+      {
+        // Ok, a bifurcation point has been found
+        // branch1: tIt->second <-> it (ok)
+        branches[ tIt->second ] = it;
+
+        // branch2: *eit <-> it (ok)
+        branches[ *eIt ] = it;
+
+        // branch3: it <-> until next bifurcation
+        bIt = tIt;
+      }
+      else
+        tags[ it ] = *eIt;
       it = p;
       p = mst->GetParent( it );
 
     } // elihw
-    // TODO: More on it
+    if( bIt != tags.end( ) )
+    {
+      TIndex pTag = bIt->second;
+      TIndex nTag = bIt->first;
+      it = bIt->first;
+      p = it;
+      while( tags[ it ] == pTag )
+      {
+        tags[ it ] = nTag;
+        p = it;
+        it = mst->GetParent( it );
+
+      } // elihw
+      branches[ bIt->first ] = p;
+    }
+    else
+    {
+      tags[ it ] = *eIt;
+      branches[ *eIt ] = it;
+
+    } // fi
+
+  } // rof
+
+  // Fill full branches
+  typedef typename _TMST::TVertices _TVertices;
+  typedef typename TSkeleton::TPath _TPath;
+
+  TSkeleton* sk = this->GetOutput( );
+  typename _TTags::const_iterator bIt = branches.begin( );
+  for( ; bIt != branches.end( ); ++bIt )
+  {
+    _TVertices v = mst->GetPath( bIt->first, bIt->second );
+    typename _TPath::Pointer path = _TPath::New( );
+    path->SetReferenceImage( this->GetInput( ) );
+    typename _TVertices::const_reverse_iterator vIt = v.rbegin( );
+    for( ; vIt != v.rend( ); ++vIt )
+      path->AddVertex( *vIt );
+    sk->AddBranch( path );
 
   } // rof
 }
@@ -296,84 +348,6 @@ _EndPoints(
   } // rof
 }
 
-
-
-/* TODO
-   this->m_DistanceMap = TDistanceMap::New( );
-   this->m_DistanceMap->InsideIsPositiveOn( );
-   this->m_DistanceMap->SquaredDistanceOff( );
-   this->m_DistanceMap->UseImageSpacingOn( );
-
-   template< class _TImage, class _TScalar >
-   void fpa::Image::SkeletonFilter< _TImage, _TScalar >::
-   _Skeleton( const std::vector< TVertex >& end_points, _TAdjacencies& A )
-   {
-   typedef typename TSkeleton::TPath _TPath;
-   typedef itk::Image< unsigned long, TImage::ImageDimension > _TTagsImage;
-
-   TMST* mst = this->GetMinimumSpanningTree( );
-   TSkeleton* 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 );
-   typename std::vector< TVertex >::const_iterator eIt = end_points.begin( );
-   for( ; eIt != end_points.end( ); ++eIt )
-   {
-   TVertex it = *eIt;
-   TVertex 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)
-   eIt = end_points.begin( );
-   for( ; eIt != end_points.end( ); ++eIt )
-   {
-   TVertex it = *eIt;
-   TVertex p = mst->GetParent( it );
-   TVertex sIdx = it;
-   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
-   }
-*/
-
 #endif // __fpa__Image__SkeletonFilter__hxx__
 
 // eof - $RCSfile$