Submission to GRETSI

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

diff --git a/lib/fpa/Image/DijkstraWithSphereBacktracking.hxx b/lib/fpa/Image/DijkstraWithSphereBacktracking.hxx
index 7d50735014422505f4e995753a0e0b54f4fb93c6..69836a95e4f49a00b336ad5b6f58d3beb1a06bf1 100644 (file)
--- a/lib/fpa/Image/DijkstraWithSphereBacktracking.hxx
+++ b/lib/fpa/Image/DijkstraWithSphereBacktracking.hxx
@@ -61,7 +61,7 @@ _BeforeMainLoop( )
   for( unsigned int d = 1; d < I::ImageDimension; ++d )
     max_spac =
       ( max_spac < double( spac[ d ] ) )? double( spac[ d ] ): max_spac;
-  max_spac *= double( 3 );
+  max_spac *= double( 30 );
 
   // Correct seeds
   for( unsigned int s = 0; s < this->m_Seeds.size( ); ++s )
@@ -99,7 +99,8 @@ _AfterMainLoop( )
 {
   // Finish base algorithm
   this->Superclass::_AfterMainLoop( );
-  this->m_FinalTree.clear( );
+  this->m_FullTree.clear( );
+  this->m_ReducedTree.clear( );
   this->m_EndPoints.clear( );
   this->m_BifurcationPoints.clear( );
   if( this->m_Candidates.size( ) == 0 )
@@ -161,7 +162,7 @@ _AfterMainLoop( )
     {
       if( start )
       {
-        if( this->m_FinalTree.find( max_vertex ) == this->m_FinalTree.end( ) )
+        if( this->m_FullTree.find( max_vertex ) == this->m_FullTree.end( ) )
         {
           // Mark a sphere around current point as visited
           double dist = std::sqrt( double( input->GetPixel( max_vertex ) ) );
@@ -173,7 +174,7 @@ _AfterMainLoop( )
 
           // Next vertex in current path
           this->InvokeEvent( TBacktrackingEvent( max_vertex, this->m_NumberOfBranches ) );
-          this->m_FinalTree[ max_vertex ] =
+          this->m_FullTree[ max_vertex ] =
             TTreeNode( this->_Parent( max_vertex ), this->m_NumberOfBranches );
           // std::cout << "New: " << this->m_NumberOfBranches << std::endl;
         }
@@ -183,7 +184,7 @@ _AfterMainLoop( )
           this->m_BifurcationPoints.push_back( max_vertex );
           this->m_NumberOfBranches++;
 
-          this->m_FinalTree[ max_vertex ] =
+          this->m_FullTree[ max_vertex ] =
             TTreeNode( this->_Parent( max_vertex ), this->m_NumberOfBranches );
           // std::cout << "Bifurcation: " << this->m_NumberOfBranches << std::endl;
 
@@ -211,7 +212,7 @@ _AfterMainLoop( )
 
           // Next vertex in current path
           this->InvokeEvent( TBacktrackingEvent( max_vertex, this->m_NumberOfBranches ) );
-          this->m_FinalTree[ max_vertex ] =
+          this->m_FullTree[ max_vertex ] =
             TTreeNode( this->_Parent( max_vertex ), this->m_NumberOfBranches );
           change = true;
           // std::cout << "Change: " << this->m_NumberOfBranches << std::endl;
@@ -221,7 +222,7 @@ _AfterMainLoop( )
           // A bifurcation point has been reached!
           // TODO: this->m_BifurcationPoints.push_back( max_vertex );
           this->m_NumberOfBranches++;
-          this->m_FinalTree[ max_vertex ] =
+          this->m_FullTree[ max_vertex ] =
             TTreeNode( this->_Parent( max_vertex ), this->m_NumberOfBranches );
           // std::cout << "Change bifurcation: " << this->m_NumberOfBranches << std::endl;
 
@@ -233,6 +234,7 @@ _AfterMainLoop( )
     } while( max_vertex != this->_Parent( max_vertex ) );
     if( start || change )
       this->m_NumberOfBranches++;
+    this->m_FullTree[ max_vertex ] = TTreeNode( max_vertex, this->m_NumberOfBranches );
 
     this->InvokeEvent( TEndBacktrackingEvent( ) );
 
@@ -240,8 +242,8 @@ _AfterMainLoop( )
 
   std::map< TMark, unsigned long > histo;
   for(
-    typename TTree::iterator treeIt = this->m_FinalTree.begin( );
-    treeIt != this->m_FinalTree.end( );
+    typename TTree::iterator treeIt = this->m_FullTree.begin( );
+    treeIt != this->m_FullTree.end( );
     ++treeIt
     )
     histo[ treeIt->second.second ]++;
@@ -257,8 +259,8 @@ _AfterMainLoop( )
   this->m_NumberOfBranches = changes.size( );
 
   for(
-    typename TTree::iterator treeIt = this->m_FinalTree.begin( );
-    treeIt != this->m_FinalTree.end( );
+    typename TTree::iterator treeIt = this->m_FullTree.begin( );
+    treeIt != this->m_FullTree.end( );
     ++treeIt
     )
   {
@@ -267,6 +269,52 @@ _AfterMainLoop( )
 
   } // fi
 
+  // Construct reduced tree
+  for(
+    typename TVertices::const_iterator eIt = this->m_EndPoints.begin( );
+    eIt != this->m_EndPoints.end( );
+    ++eIt
+    )
+  {
+    typename TTree::const_iterator tIt = this->m_FullTree.find( *eIt );
+    if( tIt != this->m_FullTree.end( ) )
+    {
+      TMark id = tIt->second.second;
+      do
+      {
+        tIt = this->m_FullTree.find( tIt->second.first );
+        if( tIt == this->m_FullTree.end( ) )
+          break;
+
+      } while( tIt->second.second == id );
+      this->m_ReducedTree[ *eIt ] = TTreeNode( tIt->first, id );
+
+    } // fi
+
+  } // rof
+
+  for(
+    typename TVertices::const_iterator bIt = this->m_BifurcationPoints.begin( );
+    bIt != this->m_BifurcationPoints.end( );
+    ++bIt
+    )
+  {
+    typename TTree::const_iterator tIt = this->m_FullTree.find( *bIt );
+    if( tIt != this->m_FullTree.end( ) )
+    {
+      TMark id = tIt->second.second;
+      do
+      {
+        tIt = this->m_FullTree.find( tIt->second.first );
+        if( tIt == this->m_FullTree.end( ) )
+          break;
+
+      } while( tIt->second.second == id );
+      this->m_ReducedTree[ *bIt ] = TTreeNode( tIt->first, id );
+
+    } // fi
+
+  } // rof
 }
 
 // -------------------------------------------------------------------------