lib/fpa/Base/MinimumSpanningTree.hxx

   1 #ifndef __FPA__BASE__MINIMUMSPANNINGTREE__HXX__
   2 #define __FPA__BASE__MINIMUMSPANNINGTREE__HXX__
   3
   4 #include <limits>
   5
   6 // -------------------------------------------------------------------------
   7 template< class V, class B >
   8 const unsigned long fpa::Base::MinimumSpanningTree< V, B >::INF_VALUE =
   9   std::numeric_limits< unsigned long >::max( ) >> 1;
  10
  11 // -------------------------------------------------------------------------
  12 template< class V, class B >
  13 void fpa::Base::MinimumSpanningTree< V, B >::
  14 SetCollisions( const TCollisions& collisions )
  15 {
  16   this->m_Collisions = collisions;
  17
  18   // Prepare fronts graph using Floyd-Warshall
  19   unsigned long nSeeds = this->m_Collisions.size( );
  20   _TMatrix dist( nSeeds, _TRow( nSeeds, Self::INF_VALUE ) );
  21   this->m_FrontPaths = dist;
  22
  23   for( unsigned long i = 0; i < nSeeds; ++i )
  24   {
  25     for( unsigned long j = 0; j < nSeeds; ++j )
  26     {
  27       if( this->m_Collisions[ i ][ j ].second )
  28       {
  29         dist[ i ][ j ] = 1;
  30         dist[ j ][ i ] = 1;
  31         this->m_FrontPaths[ i ][ j ] = j;
  32         this->m_FrontPaths[ j ][ i ] = i;
  33
  34       } // fi
  35
  36     } // rof
  37     dist[ i ][ i ] = 0;
  38     this->m_FrontPaths[ i ][ i ] = i;
  39
  40   } // rof
  41   for( unsigned long k = 0; k < nSeeds; ++k )
  42   {
  43     for( unsigned long i = 0; i < nSeeds; ++i )
  44     {
  45       for( unsigned long j = 0; j < nSeeds; ++j )
  46       {
  47         if( ( dist[ i ][ k ] + dist[ k ][ j ] ) < dist[ i ][ j ] )
  48         {
  49           dist[ i ][ j ] = dist[ i ][ k ] + dist[ k ][ j ];
  50           this->m_FrontPaths[ i ][ j ] = this->m_FrontPaths[ i ][ k ];
  51
  52         } // fi
  53
  54       } // rof
  55
  56     } // rof
  57
  58   } // rof
  59   this->Modified( );
  60 }
  61
  62 // -------------------------------------------------------------------------
  63 template< class V, class B >
  64 std::vector< V > fpa::Base::MinimumSpanningTree< V, B >::
  65 GetPath( const V& a, const V& b ) const
  66 {
  67   std::vector< V > path;
  68   typename TDecorated::const_iterator aIt = this->Get( ).find( a );
  69   typename TDecorated::const_iterator bIt = this->Get( ).find( b );
  70
  71   if( aIt == this->Get( ).end( ) || bIt == this->Get( ).end( ) )
  72     return( path );
  73
  74   short fa = aIt->second.second;
  75   short fb = bIt->second.second;
  76
  77   if( fa == fb )
  78   {
  79     std::vector< V > ap, bp;
  80     this->_Path( ap, a );
  81     this->_Path( bp, b );
  82
  83     // Ignore common part
  84     typename std::vector< V >::const_reverse_iterator raIt, rbIt;
  85     raIt = ap.rbegin( );
  86     rbIt = bp.rbegin( );
  87     while( *raIt == *rbIt && raIt != ap.rend( ) && rbIt != bp.rend( ) )
  88     {
  89       ++raIt;
  90       ++rbIt;
  91
  92     } // elihw
  93     if( raIt != ap.rbegin( ) ) --raIt;
  94     if( rbIt != bp.rbegin( ) ) --rbIt;
  95
  96     // Add part from a
  97     typename std::vector< V >::const_iterator iaIt = ap.begin( );
  98     for( ; iaIt != ap.end( ) && *iaIt != *raIt; ++iaIt )
  99       path.push_back( *iaIt );
 100
 101     // Add part from b
 102     for( ; rbIt != bp.rend( ); ++rbIt )
 103       path.push_back( *rbIt );
 104   }
 105   else
 106   {
 107     // Use this->m_FrontPaths from Floyd-Warshall
 108     if( this->m_FrontPaths[ fa ][ fb ] != Self::INF_VALUE )
 109     {
 110       // Compute front path
 111       std::vector< long > fpath;
 112       fpath.push_back( fa );
 113       while( fa != fb )
 114       {
 115         fa = this->m_FrontPaths[ fa ][ fb ];
 116         fpath.push_back( fa );
 117
 118       } // elihw
 119
 120       // Continue only if both fronts are connected
 121       unsigned int N = fpath.size( );
 122       if( 0 < N )
 123       {
 124         // First path: from start vertex to first collision
 125         path = this->GetPath(
 126           a, this->m_Collisions[ fpath[ 0 ] ][ fpath[ 1 ] ].first
 127           );
 128
 129         // Intermediary paths
 130         for( unsigned int i = 1; i < N - 1; ++i )
 131         {
 132           std::vector< V > ipath =
 133             this->GetPath(
 134               this->m_Collisions[ fpath[ i ] ][ fpath[ i - 1 ] ].first,
 135               this->m_Collisions[ fpath[ i ] ][ fpath[ i + 1 ] ].first
 136               );
 137           path.insert( path.end( ), ipath.begin( ), ipath.end( ) );
 138
 139         } // rof
 140
 141         // Final path: from last collision to end point
 142         std::vector< V > lpath =
 143           this->GetPath(
 144             this->m_Collisions[ fpath[ N - 1 ] ][ fpath[ N - 2 ] ].first, b
 145             );
 146         path.insert( path.end( ), lpath.begin( ), lpath.end( ) );
 147
 148       } // fi
 149
 150     } // fi
 151
 152   } // fi
 153   return( path );
 154 }
 155
 156 // -------------------------------------------------------------------------
 157 template< class V, class B >
 158 template< class I >
 159 std::vector< typename I::PointType > fpa::Base::MinimumSpanningTree< V, B >::
 160 GetPathFromImage(
 161   const V& a, const V& b,
 162   const I* image, unsigned int kernel
 163   ) const
 164 {
 165   typedef typename I::PointType _P;
 166
 167   std::vector< _P > path;
 168   std::vector< V > vertices = this->GetPath( a, b );
 169   for( unsigned int i = 0; i < vertices.size( ); ++i )
 170   {
 171     _P p;
 172     image->TransformIndexToPhysicalPoint( vertices[ i ], p );
 173     path.push_back( p );
 174
 175   } // rof
 176
 177   // Lowpass filter
 178   if( kernel > 0 )
 179   {
 180     int k = int( kernel ) >> 1;
 181     std::vector< _P > lowpass_path;
 182     for( unsigned int i = 0; i < path.size( ); ++i )
 183     {
 184       _P p;
 185       p.Fill( ( typename _P::ValueType )( 0 ) );
 186       unsigned int c = 0;
 187       for( int j = -k; j <= k; ++j )
 188       {
 189         int l = int( i ) + j;
 190         if( l >= 0 && l < path.size( ) )
 191         {
 192           p += path[ l ].GetVectorFromOrigin( );
 193           c++;
 194
 195         } // fi
 196
 197       } // rof
 198       if( c > 0 )
 199         for( unsigned int d = 0; d < _P::PointDimension; ++d )
 200           p[ d ] /= ( typename _P::ValueType )( c );
 201       lowpass_path.push_back( p );
 202
 203     } // rof
 204
 205     path = lowpass_path;
 206
 207   } // fi
 208   return( path );
 209 }
 210
 211 // -------------------------------------------------------------------------
 212 template< class V, class B >
 213 fpa::Base::MinimumSpanningTree< V, B >::
 214 MinimumSpanningTree( )
 215   : Superclass( )
 216 {
 217 }
 218
 219 // -------------------------------------------------------------------------
 220 template< class V, class B >
 221 fpa::Base::MinimumSpanningTree< V, B >::
 222 ~MinimumSpanningTree( )
 223 {
 224 }
 225
 226 // -------------------------------------------------------------------------
 227 template< class V, class B >
 228 void fpa::Base::MinimumSpanningTree< V, B >::
 229 _Path( std::vector< V >& path, const V& a ) const
 230 {
 231   typename TDecorated::const_iterator dIt = this->Get( ).find( a );
 232   if( dIt != this->Get( ).end( ) )
 233   {
 234     do
 235     {
 236       path.push_back( dIt->first );
 237       dIt = this->Get( ).find( dIt->second.first );
 238
 239     } while( dIt->first != dIt->second.first && dIt != this->Get( ).end( ) );
 240
 241   } // fi
 242 }
 243
 244 #endif // __FPA__BASE__MINIMUMSPANNINGTREE__HXX__
 245
 246 // eof - $RCSfile$