1 #ifndef __fpa__Image__EndPointsFilter__hxx__
2 #define __fpa__Image__EndPointsFilter__hxx__
7 #include <itkImageRegionConstIteratorWithIndex.h>
9 // -------------------------------------------------------------------------
10 template< class _TDistanceMap, class _TCostMap >
11 const typename fpa::Image::EndPointsFilter< _TDistanceMap, _TCostMap >::
12 TIndices& fpa::Image::EndPointsFilter< _TDistanceMap, _TCostMap >::
13 GetBifurcations( ) const
15 return( this->m_Bifurcations );
18 // -------------------------------------------------------------------------
19 template< class _TDistanceMap, class _TCostMap >
20 const typename fpa::Image::EndPointsFilter< _TDistanceMap, _TCostMap >::
21 TIndices& fpa::Image::EndPointsFilter< _TDistanceMap, _TCostMap >::
24 return( this->m_EndPoints );
27 // -------------------------------------------------------------------------
28 template< class _TDistanceMap, class _TCostMap >
29 void fpa::Image::EndPointsFilter< _TDistanceMap, _TCostMap >::
32 typedef itk::ImageRegionConstIteratorWithIndex< _TDistanceMap > _TDistMapIt;
33 typedef itk::ImageRegionConstIteratorWithIndex< _TCostMap > _TCostMapIt;
34 typedef std::multimap< double, TIndex, std::greater< double > > _TQueue;
35 typedef typename _TQueue::value_type _TQueueValue;
38 typename _TDistanceMap::RegionType region =
39 this->m_DistanceMap->GetRequestedRegion( );
43 _TDistMapIt dIt( this->m_DistanceMap, region );
44 _TCostMapIt cIt( this->m_CostMap, region );
47 for( ; !dIt.IsAtEnd( ) && !cIt.IsAtEnd( ); ++dIt, ++cIt )
49 double d = double( dIt.Get( ) );
52 double v = double( cIt.Get( ) ) / d;
53 queue.insert( _TQueueValue( v, dIt.GetIndex( ) ) );
59 // BFS from maximum queue
61 while( queue.size( ) > 0 )
64 auto nIt = queue.begin( );
67 if( marks.find( n.second ) != marks.end( ) )
71 marks.insert( n.second );
72 this->m_EndPoints.insert( n.second );
75 typename TMST::TPath::Pointer path;
76 this->m_MST->GetPath( path, n.second );
77 for( unsigned long i = 0; i < path->GetSize( ); ++i )
79 typename TMST::TPath::TContinuousIndex cidx = path->GetVertex( i );
80 typename _TCostMap::PointType cnt;
81 this->m_CostMap->TransformContinuousIndexToPhysicalPoint( cidx, cnt );
83 this->m_CostMap->TransformPhysicalPointToIndex( cnt, idx );
84 double d = double( this->m_DistanceMap->GetPixel( idx ) );
88 for( unsigned int d = 0; d < _TCostMap::ImageDimension; ++d )
100 while( queue.size( ) > 0 )
102 marks.insert( n.second );
103 this->m_EndPoints.insert( n.second );
104 auto path = this->m_MST->GetPath( n.second );
105 std::cout << path.size( ) << std::endl;
106 for( auto pIt = path.begin( ); pIt != path.end( ); ++pIt )
108 double d = double( this->m_DistanceMap->GetPixel( *pIt ) );
109 d = std::sqrt( std::fabs( d ) );
110 typename _TCostMap::PointType center;
113 std::queue< TIndex > q;
116 while( q.size( ) > 0 )
118 TIndex idx = q.front( );
120 if( m.find( idx ) != m.end( ) )
124 for( unsigned int x = 0; x < _TCostMap::ImageDimension; ++x )
126 for( int y = -1; y <= 1; y += 2 )
130 typename _TCostMap::PointType c;
131 this->m_CostMap->TransformIndexToPhysicalPoint( idx2, c );
132 if( this->m_CostMap->GetRequestedRegion( ).IsInside( idx2 ) )
134 if( center.EuclideanDistanceTo( c ) <= ( d * 1.5 ) )
153 // -------------------------------------------------------------------------
154 template< class _TDistanceMap, class _TCostMap >
155 fpa::Image::EndPointsFilter< _TDistanceMap, _TCostMap >::
161 // -------------------------------------------------------------------------
162 template< class _TDistanceMap, class _TCostMap >
163 fpa::Image::EndPointsFilter< _TDistanceMap, _TCostMap >::
168 #endif // __fpa__Image__EndPointsFilter__hxx__