+++ /dev/null
-#ifndef __FPA__IMAGE__REGIONGROWWITHMULTIPLETHRESHOLDS__HXX__
-#define __FPA__IMAGE__REGIONGROWWITHMULTIPLETHRESHOLDS__HXX__
-
-#include <limits>
-#include <itkBinaryThresholdImageFilter.h>
-#include <itkConstNeighborhoodIterator.h>
-
-// -------------------------------------------------------------------------
-template< class I >
-void fpa::Image::RegionGrowWithMultipleThresholds< I >::
-AddThreshold( const TPixel& v )
-{
- this->m_Thresholds.insert( v );
- this->Modified( );
-}
-
-// -------------------------------------------------------------------------
-template< class I >
-void fpa::Image::RegionGrowWithMultipleThresholds< I >::
-AddThresholds( const TPixel& t0, const TPixel& t1, const unsigned int& s )
-{
- for( TPixel t = t0; t <= t1; t += TPixel( s ) )
- this->AddThreshold( t );
-}
-
-// -------------------------------------------------------------------------
-template< class I >
-fpa::Image::RegionGrowWithMultipleThresholds< I >::
-RegionGrowWithMultipleThresholds( )
- : Superclass( ),
- m_InsideValue( TPixel( 1 ) ),
- m_OutsideValue( TPixel( 0 ) ),
- m_DifferenceThreshold( double( 3 ) ),
- m_TotalCount( 0 ),
- m_LastDiff( double( 0 ) ),
- m_StopForced( false ),
- m_StopThreshold( TPixel( 0 ) )
-{
-}
-
-// -------------------------------------------------------------------------
-template< class I >
-fpa::Image::RegionGrowWithMultipleThresholds< I >::
-~RegionGrowWithMultipleThresholds( )
-{
-}
-
-// -------------------------------------------------------------------------
-template< class I >
-void fpa::Image::RegionGrowWithMultipleThresholds< I >::
-_BeforeMainLoop( )
-{
- const I* img = this->GetInput( );
-
- // Clear all states
- this->ClearMembershipFunctions( );
- this->m_Histogram.clear( );
- this->m_TotalCount = img->GetLargestPossibleRegion( ).GetNumberOfPixels( );
- this->m_LastDiff = double( -1 );
- this->m_StopForced = false;
-
- // Initialize thresholding functors
- typename TThresholds::const_iterator tIt = this->m_Thresholds.begin( );
- TPixel min_thr = std::numeric_limits< TPixel >::min( );
- for( ; tIt != this->m_Thresholds.end( ); ++tIt )
- {
- typename TFunction::Pointer function = TFunction::New( );
- function->SetInputImage( img );
- function->SetLowerThreshold( min_thr );
- function->SetUpperThreshold( *tIt );
- this->AddMembershipFunction( function );
-
- } // rof
-
- // Correct seeds
- typename I::SizeType radius;
- radius.Fill( 3 );
- itk::ConstNeighborhoodIterator< I > it(
- radius, img, img->GetRequestedRegion( )
- );
- for( unsigned int s = 0; s < this->m_Seeds.size( ); ++s )
- {
- _TNode seed = this->m_Seeds[ s ];
- it.SetLocation( seed.Vertex );
-
- typename I::SizeType size = it.GetSize( );
- unsigned int nN = 1;
- for( unsigned int d = 0; d < I::ImageDimension; ++d )
- nN *= size[ d ];
- TPixel min_value = img->GetPixel( seed.Vertex );
- for( unsigned int i = 0; i < nN; ++i )
- {
- if( it.GetPixel( i ) < min_value )
- {
- seed.Vertex = it.GetIndex( i );
- seed.Parent = seed.Vertex;
- min_value = it.GetPixel( i );
-
- } // fi
-
- } // rof
- this->m_Seeds[ s ] = seed;
-
- } // rof
-
- // Continue all initializations
- this->Superclass::_BeforeMainLoop( );
-}
-
-// -------------------------------------------------------------------------
-template< class I >
-void fpa::Image::RegionGrowWithMultipleThresholds< I >::
-_AfterMainLoop( )
-{
- typedef itk::BinaryThresholdImageFilter< I, I > _TBinFilter;
-
- // Binarize, inplace, the grown region
- if( this->m_Histogram.size( ) > 1 )
- {
- typename _TBinFilter::Pointer bin = _TBinFilter::New( );
- bin->SetInput( this->GetOutput( ) );
- bin->SetInsideValue( this->m_InsideValue );
- bin->SetOutsideValue( this->m_OutsideValue );
- bin->InPlaceOn( );
- if( this->m_StopForced )
- bin->SetUpperThreshold( this->m_StopThreshold );
- else
- bin->SetUpperThreshold( this->m_Histogram.rbegin( )->first );
- bin->GraftOutput( this->GetOutput( ) );
- bin->Update( );
- this->GraftOutput( bin->GetOutput( ) );
-
- } // fi
-
- this->Superclass::_AfterMainLoop( );
-}
-
-// -------------------------------------------------------------------------
-template< class I >
-void fpa::Image::RegionGrowWithMultipleThresholds< I >::
-_AfterLoop( )
-{
- if( this->m_ActualFunction != this->m_Functions.end( ) )
- {
- TFunction* f =
- dynamic_cast< TFunction* >( this->m_ActualFunction->GetPointer( ) );
- this->m_Histogram[ f->GetUpperThreshold( ) ] = this->m_Marks.size( );
-
- /* TODO: code to print the signal for demo purposes
- std::cout
- << f->GetUpperThreshold( ) << " "
- << this->m_Marks.size( )
- << std::endl;
- */
-
- // Get the two last complete count, if any, and compute last
- // finite difference
- if( this->m_Histogram.size( ) > 1 )
- {
- typename THistogram::const_reverse_iterator hIt =
- this->m_Histogram.rbegin( );
- typename THistogram::const_reverse_iterator gIt = hIt;
- gIt++;
-
- this->m_LastDiff = double( hIt->second ) - double( gIt->second );
- this->m_LastDiff /= double( hIt->first ) - double( gIt->first );
- this->m_LastDiff *= m_DifferenceThreshold;
- }
- else
- this->m_LastDiff = double( -1 );
-
- } // fi
- this->Superclass::_AfterLoop( );
-}
-
-// -------------------------------------------------------------------------
-template< class I >
-bool fpa::Image::RegionGrowWithMultipleThresholds< I >::
-_UpdateResult( _TNode& n )
-{
- // Here, the output is changed to the threshold used to compute the
- // growing condition of the actual vertex n
- bool ret = this->Superclass::_UpdateResult( n );
- if( this->m_ActualFunction != this->m_Functions.end( ) )
- {
- TFunction* f =
- dynamic_cast< TFunction* >( this->m_ActualFunction->GetPointer( ) );
- if( f != NULL )
- this->GetOutput( )->SetPixel( n.Vertex, f->GetUpperThreshold( ) );
-
- } // fi
- return( ret );
-}
-
-// -------------------------------------------------------------------------
-template< class I >
-void fpa::Image::RegionGrowWithMultipleThresholds< I >::
-_Mark( const _TNode& n )
-{
- this->Superclass::_Mark( n );
-
- // Check if the histogram's support is enough
- if( this->m_Histogram.size( ) < 2 )
- return;
- typename THistogram::const_reverse_iterator hIt =
- this->m_Histogram.rbegin( );
-
- // Get the actual function
- TFunction* f =
- dynamic_cast< TFunction* >( this->m_ActualFunction->GetPointer( ) );
- if( f == NULL )
- return;
-
- // Previous pixel count
- unsigned long prev_aCount = hIt->second;
- double prev_p = double( prev_aCount ) / double( this->m_TotalCount );
-
- // Actual pixel count
- unsigned long aCount = this->m_Marks.size( );
- double p = double( aCount ) / double( this->m_TotalCount );
-
- // Loop over, at least, 0.1% of the total number of pixels before
- // performing stop analysis
- if( double( 1e-3 ) < p && double( 1e-3 ) < prev_p )
- {
- // Does the difference is worthy to be analyzed?
- if( aCount > hIt->second )
- {
- // Compute finite difference and compare it to the previous one
- double diff = double( aCount ) - double( hIt->second );
- diff /= double( f->GetUpperThreshold( ) ) - double( hIt->first );
- if( diff > this->m_LastDiff )
- {
- /* TODO: comment this for demo purposes
- */
- this->m_StopForced = true;
- this->m_StopThreshold = hIt->first;
-
- } // fi
-
- } // fi
-
- } // fi
-}
-
-// -------------------------------------------------------------------------
-template< class I >
-bool fpa::Image::RegionGrowWithMultipleThresholds< I >::
-_CheckStopCondition( )
-{
- return( this->m_StopForced );
-}
-
-#endif // __FPA__IMAGE__REGIONGROWWITHMULTIPLETHRESHOLDS__HXX__
-
-// eof - $RCSfile$