i cpPlugins_Paths.h
+t fpa/Image/Functors/SimpleNeighborhood
+t fpa/Image/Functors/SimpleDijkstraCost
t fpa/Image/MinimumSpanningTree
t fpa/Base/MinimumSpanningTree
+c fpa::Image::Functors::SimpleNeighborhood< itk::ImageBase< #process_dims# > >
+c fpa::Image::Functors::SimpleDijkstraCost< itk::Image< #pixels#, #process_dims# >, #reals# >
c fpa::Base::MinimumSpanningTree< itk::Index< #process_dims# >, cpExtensions::DataStructures::PolyLineParametricPath< #process_dims# >, itk::Image< itk::Offset< #process_dims# >, #process_dims# > >
c fpa::Image::MinimumSpanningTree< #process_dims# >
} // elihw
vertices.push_back( it );
- path = _TPath::New( );
+ if( path.IsNull( ) )
+ path = _TPath::New( );
for( auto v = vertices.begin( ); v != vertices.end( ); ++v )
path->AddVertex( *v );
}
typename _TPath::Pointer& path, const _TVertex& a, const _TVertex& b
) const
{
-/*
static const unsigned long _inf =
std::numeric_limits< unsigned long >::max( );
- _TPath path;
- _TPath pa = this->GetPath( a );
- _TPath pb = this->GetPath( b );
- if( pa.size( ) > 0 && pb.size( ) > 0 )
+ if( path.IsNull( ) )
+ path = _TPath::New( );
+ typename _TPath::Pointer pa, pb;
+ this->GetPath( pa, a );
+ this->GetPath( pb, b );
+ if( pa->GetSize( ) > 0 && pb->GetSize( ) > 0 )
{
// Find front identifiers
unsigned long ia = _inf, ib = _inf;
unsigned long N = this->m_Seeds.size( );
for( unsigned long i = 0; i < N; ++i )
{
- if( this->m_Seeds[ i ] == pa.front( ) )
+ if( this->m_Seeds[ i ] == pa->GetVertex( 0 ) )
ia = i;
- if( this->m_Seeds[ i ] == pb.front( ) )
+ if( this->m_Seeds[ i ] == pb->GetVertex( 0 ) )
ib = i;
} // rof
+ // Check if there is a front-jump between given seeds
if( ia != ib )
{
- // Use this->m_FrontPaths from Floyd-Warshall
- if( this->m_FrontPaths[ ia ][ ib ] < _inf )
- {
- // Compute front path
- std::vector< long > fpath;
- fpath.push_back( ia );
- while( ia != ib )
- {
- ia = this->m_FrontPaths[ ia ][ ib ];
- fpath.push_back( ia );
-
- } // elihw
-
- // 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
- );
-
- // 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( ) );
-
- } // rof
-
- // 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( ) );
-
- } // fi
-
- } // fi
+ // Compute front path
+ /* TODO
+ std::vector< long > fpath;
+ fpath.push_back( ia );
+ while( ia != ib )
+ {
+ ia = this->m_FrontPaths[ ia ][ ib ];
+ fpath.push_back( ia );
+
+ } // elihw
+
+ // Continue only if both fronts are connected
+ unsigned int N = fpath.size( );
+ if( N > 0 )
+ {
+ // First path: from start vertex to first collision
+ this->GetPath(
+ path, a,
+ this->m_Collisions[ fpath[ 0 ] ][ fpath[ 1 ] ].first
+ );
+
+ // Intermediary paths
+ for( unsigned int i = 1; i < N - 1; ++i )
+ {
+ TPathTVertices 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( ) );
+
+ } // rof
+
+ // 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( ) );
+
+ } // fi
+
+ } // fi
+ */
}
else
{
// Ignore common part: find common ancestor
- auto aIt = pa.begin( );
- auto bIt = pb.begin( );
- while( *aIt == *bIt && aIt != pa.end( ) && bIt != pb.end( ) )
+ long aIt = pa->GetSize( ) - 1;
+ long bIt = pb->GetSize( ) - 1;
+ bool cont = true;
+ while( aIt >= 0 && bIt >= 0 && cont )
{
- ++aIt;
- ++bIt;
+ cont = ( pa->GetVertex( aIt ) == pb->GetVertex( bIt ) );
+ aIt--;
+ bIt--;
} // elihw
+ aIt++;
+ bIt++;
// Glue both parts
- for( --aIt; aIt != pa.end( ); ++aIt )
- path.push_front( *aIt );
- for( ; bIt != pb.end( ); ++bIt )
- path.push_back( *bIt );
+ for( long cIt = 0; cIt <= aIt; ++cIt )
+ path->AddVertex( pa->GetVertex( cIt ) );
+ for( ; bIt >= 0; --bIt )
+ path->AddVertex( pb->GetVertex( bIt ) );
} // fi
} // fi
- return( path );
-*/
}
// -------------------------------------------------------------------------
// Send Piotr's code to Anna
#include <itkImage.h>
+#include <fpa/Image/Functors/SimpleNeighborhood.h>
// -------------------------------------------------------------------------
template< class _TInputImage, class _TOutputImage >
Algorithm( )
: Superclass( )
{
+ typedef itk::ImageBase< _TInputImage::ImageDimension > _TImageBase;
+ typedef fpa::Image::Functors::SimpleNeighborhood< _TImageBase > _TNeigh;
typedef itk::Image< TFrontId, _TInputImage::ImageDimension > _TMarks;
+
this->m_MarksIdx = this->GetNumberOfRequiredOutputs( );
this->itk::ProcessObject::SetNumberOfRequiredOutputs( this->m_MarksIdx + 1 );
typename _TMarks::Pointer marks = _TMarks::New( );
this->SetNthOutput( this->m_MarksIdx, marks );
+
+ typename _TNeigh::Pointer neigh = _TNeigh::New( );
+ this->SetNeighborhoodFunction( neigh );
}
// -------------------------------------------------------------------------
this->AllocateOutputs( );
TNeighborhoodFunction* neighFunc =
- dynamic_cast< TNeighborhoodFunction* >( this->GetNeighborhoodFunction( ) );
+ dynamic_cast< TNeighborhoodFunction* >(
+ this->GetNeighborhoodFunction( )
+ );
if( neighFunc == NULL )
itkExceptionMacro( << "NeighborhoodFunction not well defined." );
neighFunc->SetImage( this->GetInput( ) );
#ifndef __fpa__Image__Dijkstra__hxx__
#define __fpa__Image__Dijkstra__hxx__
+#include <fpa/Image/Functors/SimpleDijkstraCost.h>
+
// -------------------------------------------------------------------------
template< class _TInputImage, class _TOutputImage >
fpa::Image::Dijkstra< _TInputImage, _TOutputImage >::
Dijkstra( )
: Superclass( )
{
+ typedef
+ fpa::Image::Functors::SimpleDijkstraCost< _TInputImage, TOutput >
+ _TCost;
+
+ typename _TCost::Pointer cost = _TCost::New( );
+ this->SetCostFunction( cost );
}
// -------------------------------------------------------------------------
typedef itk::ImageRegionConstIteratorWithIndex< _TCostMap > _TCostMapIt;
typedef std::multimap< double, TIndex, std::greater< double > > _TQueue;
typedef typename _TQueue::value_type _TQueueValue;
+ typedef itk::Image< unsigned char, _TCostMap::ImageDimension > _TMarks;
// Some values
typename _TDistanceMap::RegionType region =
this->m_DistanceMap->GetRequestedRegion( );
+ typename _TMarks::Pointer marks = _TMarks::New( );
+ marks->SetLargestPossibleRegion( this->m_DistanceMap->GetLargestPossibleRegion( ) );
+ marks->SetRequestedRegion( this->m_DistanceMap->GetRequestedRegion( ) );
+ marks->SetBufferedRegion( this->m_DistanceMap->GetBufferedRegion( ) );
+ marks->SetSpacing( this->m_DistanceMap->GetSpacing( ) );
+ marks->SetOrigin( this->m_DistanceMap->GetOrigin( ) );
+ marks->SetDirection( this->m_DistanceMap->GetDirection( ) );
+ marks->Allocate( );
+ marks->FillBuffer( 0 );
// Create queue
_TQueue queue;
} // rof
// BFS from maximum queue
- TIndices marks;
while( queue.size( ) > 0 )
{
// Get node
auto nIt = queue.begin( );
auto n = *nIt;
queue.erase( nIt );
- if( marks.find( n.second ) != marks.end( ) )
+
+ unsigned char m = marks->GetPixel( n.second );
+ if( ( m & 0x01 ) == 0x01 )
continue;
// Mark it
- marks.insert( n.second );
+ m |= 0x01;
+ marks->SetPixel( n.second, m );
this->m_EndPoints.insert( n.second );
// Get path
this->m_MST->GetPath( path, n.second );
for( unsigned long i = 0; i < path->GetSize( ); ++i )
{
- typename TMST::TPath::TContinuousIndex cidx = path->GetVertex( i );
+ TIndex idx = path->GetVertex( path->GetSize( ) - 1 - i );
typename _TCostMap::PointType cnt;
- this->m_CostMap->TransformContinuousIndexToPhysicalPoint( cidx, cnt );
- TIndex idx;
- this->m_CostMap->TransformPhysicalPointToIndex( cnt, idx );
+ this->m_CostMap->TransformIndexToPhysicalPoint( idx, cnt );
double d = double( this->m_DistanceMap->GetPixel( idx ) );
-
- /* TODO
- TIndex idx;
- for( unsigned int d = 0; d < _TCostMap::ImageDimension; ++d )
- idx[ d ] = cidx[ d ];
- */
+ d *= double( 2 );
+
+ // Mark sphere
+ std::queue< TIndex > q;
+ q.push( idx );
+ while( q.size( ) > 0 )
+ {
+ TIndex v = q.front( );
+ q.pop( );
+ unsigned char m = marks->GetPixel( v );
+ if( ( m & 0x02 ) == 0x02 )
+ continue;
+ m |= 0x03;
+ marks->SetPixel( v, m );
+
+ for( unsigned int x = 0; x < _TCostMap::ImageDimension; ++x )
+ {
+ for( int y = -1; y <= 1; y += 2 )
+ {
+ TIndex w = v;
+ w[ x ] += y;
+ if( region.IsInside( w ) )
+ {
+ typename _TCostMap::PointType p;
+ this->m_CostMap->TransformIndexToPhysicalPoint( w, p );
+ if( cnt.EuclideanDistanceTo( p ) <= d )
+ q.push( w );
+
+ } // fi
+
+ } // rof
+
+ } // rof
+
+ } // elihw
} // rof
- // TODO: temporary
- queue.clear( );
-
} // elihw
-
- /* TODO
- while( queue.size( ) > 0 )
- {
- marks.insert( n.second );
- this->m_EndPoints.insert( n.second );
- auto path = this->m_MST->GetPath( n.second );
- std::cout << path.size( ) << std::endl;
- for( auto pIt = path.begin( ); pIt != path.end( ); ++pIt )
- {
- double d = double( this->m_DistanceMap->GetPixel( *pIt ) );
- d = std::sqrt( std::fabs( d ) );
- typename _TCostMap::PointType center;
-
-
- std::queue< TIndex > q;
- TIndices m;
- q.push( *pIt );
- while( q.size( ) > 0 )
- {
- TIndex idx = q.front( );
- q.pop( );
- if( m.find( idx ) != m.end( ) )
- continue;
- m.insert( idx );
- marks.insert( idx );
- for( unsigned int x = 0; x < _TCostMap::ImageDimension; ++x )
- {
- for( int y = -1; y <= 1; y += 2 )
- {
- TIndex idx2 = idx;
- idx2[ x ] += y;
- typename _TCostMap::PointType c;
- this->m_CostMap->TransformIndexToPhysicalPoint( idx2, c );
- if( this->m_CostMap->GetRequestedRegion( ).IsInside( idx2 ) )
- {
- if( center.EuclideanDistanceTo( c ) <= ( d * 1.5 ) )
- q.push( idx2 );
-
- } // fi
-
- } // rof
-
- } // rof
-
- } // elihw
-
- } // rof
-
- } // fi
-
- } // elihw
- */
}
// -------------------------------------------------------------------------
// Configure filter
auto filter = this->_CreateITK< _TFilter >( );
filter->SetInput( image );
- filter->SetNeighborhoodFunction( neig );
- filter->SetCostFunction( cost );
- filter->SetCostConversionFunction( conv );
+ if( neig != NULL )
+ filter->SetNeighborhoodFunction( neig );
+ if( cost != NULL )
+ filter->SetCostFunction( cost );
+ if( conv != NULL )
+ filter->SetCostConversionFunction( conv );
filter->SetStopAtOneFront( this->m_Parameters.GetBool( "StopAtOneFront" ) );
// Assign seeds
#include <plugins/Plugins/SimpleImageDijkstraCost.h>
#include <cpPlugins/DataObjects/Image.h>
-#include <fpa/Image/Functors/SimpleDijkstraCost.h>
-#include <fpa/Image/Functors/SimpleDijkstraCost.hxx>
+#include <fpa_DataObjects.h>
// -------------------------------------------------------------------------
fpaPlugins::SimpleImageDijkstraCost::
#include <plugins/Plugins/SimpleImageNeighborhood.h>
#include <cpPlugins/DataObjects/Image.h>
-#include <fpa/Image/Functors/SimpleNeighborhood.h>
-#include <fpa/Image/Functors/SimpleNeighborhood.hxx>
+#include <fpa_DataObjects.h>
// -------------------------------------------------------------------------
fpaPlugins::SimpleImageNeighborhood::