]> Creatis software - FrontAlgorithms.git/blobdiff - lib/fpa/Base/MinimumSpanningTree.hxx
Windows compilation (1/2)
[FrontAlgorithms.git] / lib / fpa / Base / MinimumSpanningTree.hxx
index 9d0939ddf74b2e469c4416c103f1f6445ee96214..8746005f8a776b11f8021a65a55d86341b65eea1 100644 (file)
@@ -1,28 +1,36 @@
-#ifndef __FPA__BASE__MINIMUMSPANNINGTREE__HXX__
-#define __FPA__BASE__MINIMUMSPANNINGTREE__HXX__
+#ifndef __fpa__Base__MinimumSpanningTree__hxx__
+#define __fpa__Base__MinimumSpanningTree__hxx__
 
 #include <limits>
 
 // -------------------------------------------------------------------------
-template< class V, class B >
-const unsigned long fpa::Base::MinimumSpanningTree< V, B >::INF_VALUE =
-  std::numeric_limits< unsigned long >::max( ) >> 1;
+template< class _TVertex, class _TSuperclass >
+const typename fpa::Base::MinimumSpanningTree< _TVertex, _TSuperclass >::
+TCollisions& fpa::Base::MinimumSpanningTree< _TVertex, _TSuperclass >::
+GetCollisions( ) const
+{
+  return( this->m_Collisions );
+}
 
 // -------------------------------------------------------------------------
-template< class V, class B >
-void fpa::Base::MinimumSpanningTree< V, B >::
+template< class _TVertex, class _TSuperclass >
+void fpa::Base::MinimumSpanningTree< _TVertex, _TSuperclass >::
 SetCollisions( const TCollisions& collisions )
 {
+  static const unsigned long _inf =
+    std::numeric_limits< unsigned long >::max( );
+  if( this->m_Collisions == collisions )
+    return;
+
   this->m_Collisions = collisions;
 
-  // Prepare fronts graph using Floyd-Warshall
-  unsigned long nSeeds = this->m_Collisions.size( );
-  _TMatrix dist( nSeeds, _TRow( nSeeds, Self::INF_VALUE ) );
+  // Prepare a front graph
+  unsigned long N = this->m_Collisions.size( );
+  _TMatrix dist( N, _TRow( N, _inf ) );
   this->m_FrontPaths = dist;
-
-  for( unsigned long i = 0; i < nSeeds; ++i )
+  for( unsigned long i = 0; i < N; ++i )
   {
-    for( unsigned long j = 0; j < nSeeds; ++j )
+    for( unsigned long j = 0; j < N; ++j )
     {
       if( this->m_Collisions[ i ][ j ].second )
       {
@@ -38,15 +46,25 @@ SetCollisions( const TCollisions& collisions )
     this->m_FrontPaths[ i ][ i ] = i;
 
   } // rof
-  for( unsigned long k = 0; k < nSeeds; ++k )
+
+  // Use Floyd-Warshall to compute all possible paths between fronts
+  for( unsigned long k = 0; k < N; ++k )
   {
-    for( unsigned long i = 0; i < nSeeds; ++i )
+    for( unsigned long i = 0; i < N; ++i )
     {
-      for( unsigned long j = 0; j < nSeeds; ++j )
+      for( unsigned long j = 0; j < N; ++j )
       {
-        if( ( dist[ i ][ k ] + dist[ k ][ j ] ) < dist[ i ][ j ] )
+        // WARNING: you don't want a numeric overflow!!!
+        unsigned long dik = dist[ i ][ k ];
+        unsigned long dkj = dist[ k ][ j ];
+        unsigned long sum = _inf;
+        if( dik < _inf && dkj < _inf )
+          sum = dik + dkj;
+
+        // Ok, continue Floyd-Warshall
+        if( sum < dist[ i ][ j ] )
         {
-          dist[ i ][ j ] = dist[ i ][ k ] + dist[ k ][ j ];
+          dist[ i ][ j ] = sum;
           this->m_FrontPaths[ i ][ j ] = this->m_FrontPaths[ i ][ k ];
 
         } // fi
@@ -60,190 +78,155 @@ SetCollisions( const TCollisions& collisions )
 }
 
 // -------------------------------------------------------------------------
-template< class V, class B >
-std::vector< V > fpa::Base::MinimumSpanningTree< V, B >::
-GetPath( const V& a, const V& b ) const
+template< class _TVertex, class _TSuperclass >
+void fpa::Base::MinimumSpanningTree< _TVertex, _TSuperclass >::
+ClearSeeds( )
 {
-  std::vector< V > path;
-  typename TDecorated::const_iterator aIt = this->Get( ).find( a );
-  typename TDecorated::const_iterator bIt = this->Get( ).find( b );
+  this->m_Seeds.clear( );
+  this->Modified( );
+}
 
-  if( aIt == this->Get( ).end( ) || bIt == this->Get( ).end( ) )
-    return( path );
-  
-  short fa = aIt->second.second;
-  short fb = bIt->second.second;
+// -------------------------------------------------------------------------
+template< class _TVertex, class _TSuperclass >
+void fpa::Base::MinimumSpanningTree< _TVertex, _TSuperclass >::
+AddSeed( const _TVertex& seed )
+{
+  this->m_Seeds.push_back( seed );
+  this->Modified( );
+}
 
-  if( fa == fb )
+// -------------------------------------------------------------------------
+template< class _TVertex, class _TSuperclass >
+typename fpa::Base::MinimumSpanningTree< _TVertex, _TSuperclass >::
+TVertices fpa::Base::MinimumSpanningTree< _TVertex, _TSuperclass >::
+GetPath( const _TVertex& a ) const
+{
+  TVertices path;
+  _TVertex it = a;
+  _TVertex p = this->GetParent( it );
+  while( it != p )
   {
-    std::vector< V > ap, bp;
-    this->_Path( ap, a );
-    this->_Path( bp, b );
-
-    // Ignore common part
-    typename std::vector< V >::const_reverse_iterator raIt, rbIt;
-    raIt = ap.rbegin( );
-    rbIt = bp.rbegin( );
-    while( *raIt == *rbIt && raIt != ap.rend( ) && rbIt != bp.rend( ) )
-    {
-      ++raIt;
-      ++rbIt;
-
-    } // elihw
-    if( raIt != ap.rbegin( ) ) --raIt;
-    if( rbIt != bp.rbegin( ) ) --rbIt;
-
-    // Add part from a
-    typename std::vector< V >::const_iterator iaIt = ap.begin( );
-    for( ; iaIt != ap.end( ) && *iaIt != *raIt; ++iaIt )
-      path.push_back( *iaIt );
-
-    // Add part from b
-    for( ; rbIt != bp.rend( ); ++rbIt )
-      path.push_back( *rbIt );
-  }
-  else
+    path.push_front( it );
+    it = p;
+    p = this->GetParent( it );
+
+  } // elihw
+  path.push_front( it );
+  return( path );
+}
+
+// -------------------------------------------------------------------------
+template< class _TVertex, class _TSuperclass >
+typename fpa::Base::MinimumSpanningTree< _TVertex, _TSuperclass >::
+TVertices fpa::Base::MinimumSpanningTree< _TVertex, _TSuperclass >::
+GetPath( const _TVertex& a, const _TVertex& b ) const
+{
+  static const unsigned long _inf =
+    std::numeric_limits< unsigned long >::max( );
+
+  TVertices path;
+  TVertices pa = this->GetPath( a );
+  TVertices pb = this->GetPath( b );
+  if( pa.size( ) > 0 && pb.size( ) > 0 )
   {
-    // Use this->m_FrontPaths from Floyd-Warshall
-    if( this->m_FrontPaths[ fa ][ fb ] != Self::INF_VALUE )
+    // Find front identifiers
+    unsigned long ia = _inf, ib = _inf;
+    unsigned long N = this->m_Seeds.size( );
+    for( unsigned long i = 0; i < N; ++i )
     {
-      // Compute front path
-      std::vector< long > fpath;
-      fpath.push_back( fa );
-      while( fa != fb )
-      {
-        fa = this->m_FrontPaths[ fa ][ fb ];
-        fpath.push_back( fa );
+      if( this->m_Seeds[ i ] == pa.front( ) )
+        ia = i;
+      if( this->m_Seeds[ i ] == pb.front( ) )
+        ib = i;
 
-      } // elihw
+    } // rof
 
-      // Continue only if both fronts are connected
-      unsigned int N = fpath.size( );
-      if( 0 < N )
+    if( ia != ib )
+    {
+      // Use this->m_FrontPaths from Floyd-Warshall
+      if( this->m_FrontPaths[ ia ][ ib ] < _inf )
       {
-        // First path: from start vertex to first collision
-        path = this->GetPath(
-          a, this->m_Collisions[ fpath[ 0 ] ][ fpath[ 1 ] ].first
-          );
-
-        // Intermediary paths
-        for( unsigned int i = 1; i < N - 1; ++i )
+        // Compute front path
+        std::vector< long > fpath;
+        fpath.push_back( ia );
+        while( ia != ib )
         {
-          std::vector< V > ipath =
-            this->GetPath(
-              this->m_Collisions[ fpath[ i ] ][ fpath[ i - 1 ] ].first,
-              this->m_Collisions[ fpath[ i ] ][ fpath[ i + 1 ] ].first
-              );
-          path.insert( path.end( ), ipath.begin( ), ipath.end( ) );
+          ia = this->m_FrontPaths[ ia ][ ib ];
+          fpath.push_back( ia );
 
-        } // rof
+        } // elihw
 
-        // Final path: from last collision to end point
-        std::vector< V > lpath =
-          this->GetPath(
-            this->m_Collisions[ fpath[ N - 1 ] ][ fpath[ N - 2 ] ].first, b
+        // Continue only if both fronts are connected
+        unsigned int N = fpath.size( );
+        if( N > 0 )
+        {
+          // First path: from start vertex to first collision
+          path = this->GetPath(
+            a, this->m_Collisions[ fpath[ 0 ] ][ fpath[ 1 ] ].first
             );
-        path.insert( path.end( ), lpath.begin( ), lpath.end( ) );
-
-      } // fi
 
-    } // fi
-
-  } // fi
-  return( path );
-}
+          // Intermediary paths
+          for( unsigned int i = 1; i < N - 1; ++i )
+          {
+            TVertices ipath =
+              this->GetPath(
+                this->m_Collisions[ fpath[ i ] ][ fpath[ i - 1 ] ].first,
+                this->m_Collisions[ fpath[ i ] ][ fpath[ i + 1 ] ].first
+                );
+            path.insert( path.end( ), ipath.begin( ), ipath.end( ) );
 
-// -------------------------------------------------------------------------
-template< class V, class B >
-template< class I >
-std::vector< typename I::PointType > fpa::Base::MinimumSpanningTree< V, B >::
-GetPathFromImage(
-  const V& a, const V& b,
-  const I* image, unsigned int kernel
-  ) const
-{
-  typedef typename I::PointType _P;
+          } // rof
 
-  std::vector< _P > path;
-  std::vector< V > vertices = this->GetPath( a, b );
-  for( unsigned int i = 0; i < vertices.size( ); ++i )
-  {
-    _P p;
-    image->TransformIndexToPhysicalPoint( vertices[ i ], p );
-    path.push_back( p );
+          // Final path: from last collision to end point
+          TVertices lpath =
+            this->GetPath(
+              this->m_Collisions[ fpath[ N - 1 ] ][ fpath[ N - 2 ] ].first, b
+              );
+          path.insert( path.end( ), lpath.begin( ), lpath.end( ) );
 
-  } // rof
+        } // fi
 
-  // Lowpass filter
-  if( kernel > 0 )
-  {
-    int k = int( kernel ) >> 1;
-    std::vector< _P > lowpass_path;
-    for( unsigned int i = 0; i < path.size( ); ++i )
+      } // fi
+    }
+    else
     {
-      _P p;
-      p.Fill( ( typename _P::ValueType )( 0 ) );
-      unsigned int c = 0;
-      for( int j = -k; j <= k; ++j )
+      // Ignore common part: find common ancestor
+      auto aIt = pa.begin( );
+      auto bIt = pb.begin( );
+      while( *aIt == *bIt && aIt != pa.end( ) && bIt != pb.end( ) )
       {
-        int l = int( i ) + j;
-        if( l >= 0 && l < path.size( ) )
-        {
-          p += path[ l ].GetVectorFromOrigin( );
-          c++;
+        ++aIt;
+        ++bIt;
 
-        } // fi
-
-      } // rof
-      if( c > 0 )
-        for( unsigned int d = 0; d < _P::PointDimension; ++d )
-          p[ d ] /= ( typename _P::ValueType )( c );
-      lowpass_path.push_back( p );
+      } // elihw
 
-    } // rof
+      // Glue both parts
+      for( --aIt; aIt != pa.end( ); ++aIt )
+        path.push_front( *aIt );
+      for( ; bIt != pb.end( ); ++bIt )
+        path.push_back( *bIt );
 
-    path = lowpass_path;
+    } // fi
 
   } // fi
   return( path );
 }
 
 // -------------------------------------------------------------------------
-template< class V, class B >
-fpa::Base::MinimumSpanningTree< V, B >::
+template< class _TVertex, class _TSuperclass >
+fpa::Base::MinimumSpanningTree< _TVertex, _TSuperclass >::
 MinimumSpanningTree( )
   : Superclass( )
 {
 }
 
 // -------------------------------------------------------------------------
-template< class V, class B >
-fpa::Base::MinimumSpanningTree< V, B >::
+template< class _TVertex, class _TSuperclass >
+fpa::Base::MinimumSpanningTree< _TVertex, _TSuperclass >::
 ~MinimumSpanningTree( )
 {
 }
 
-// -------------------------------------------------------------------------
-template< class V, class B >
-void fpa::Base::MinimumSpanningTree< V, B >::
-_Path( std::vector< V >& path, const V& a ) const
-{
-  typename TDecorated::const_iterator dIt = this->Get( ).find( a );
-  if( dIt != this->Get( ).end( ) )
-  {
-    do
-    {
-      path.push_back( dIt->first );
-      dIt = this->Get( ).find( dIt->second.first );
-
-    } while( dIt->first != dIt->second.first && dIt != this->Get( ).end( ) );
-
-    if( dIt != this->Get( ).end( ) )
-      path.push_back( dIt->first );
-
-  } // fi
-}
-
-#endif // __FPA__BASE__MINIMUMSPANNINGTREE__HXX__
+#endif // __fpa__Base__MinimumSpanningTree__hxx__
 
 // eof - $RCSfile$