From: Leonardo Florez-Valencia Date: Thu, 19 Feb 2015 19:24:14 +0000 (-0500) Subject: Segmentation completely debugged (at least for 10 images) X-Git-Url: https://git.creatis.insa-lyon.fr/pubgit/?a=commitdiff_plain;h=c9880eac7f555da20538eb0eb95aa7e31fab3d51;p=FrontAlgorithms.git Segmentation completely debugged (at least for 10 images) --- diff --git a/appli/examples/example_ImageAlgorithmRegionGrow_MultipleThresholds.cxx b/appli/examples/example_ImageAlgorithmRegionGrow_MultipleThresholds.cxx index ff82103..deef6de 100644 --- a/appli/examples/example_ImageAlgorithmRegionGrow_MultipleThresholds.cxx +++ b/appli/examples/example_ImageAlgorithmRegionGrow_MultipleThresholds.cxx @@ -42,14 +42,14 @@ int main( int argc, char* argv[] ) { std::cerr << "Usage: " << argv[ 0 ] - << " input_image thr_0 thr_1 n_samples" << std::endl; + << " input_image thr_0 thr_1 step" << std::endl; return( 1 ); } // fi std::string input_image_fn = argv[ 1 ]; TPixel thr_0 = TPixel( std::atof( argv[ 2 ] ) ); TPixel thr_1 = TPixel( std::atof( argv[ 3 ] ) ); - unsigned int n_samples = std::atoi( argv[ 4 ] ); + unsigned int step = std::atoi( argv[ 4 ] ); // Read image TImageReader::Pointer input_image_reader = TImageReader::New( ); @@ -99,13 +99,13 @@ int main( int argc, char* argv[] ) // Configure algorithm TFrontAlgorithm::Pointer algorithm = TFrontAlgorithm::New( ); - algorithm->AddThresholds( thr_0, thr_1, n_samples ); + algorithm->AddThresholds( thr_0, thr_1, step ); algorithm->AddSeed( seed_idx, 0 ); algorithm->AddObserver( itk::AnyEvent( ), obs ); algorithm->ThrowEventsOn( ); algorithm->SetInput( input_image ); algorithm->SetNeighborhoodOrder( 1 ); - algorithm->SetDerivativeThreshold( double( 3 ) ); + algorithm->SetDifferenceThreshold( double( 3 ) ); algorithm->Update( ); // Let some interaction and close program diff --git a/lib/fpa/Image/RegionGrowWithMultipleThresholds.h b/lib/fpa/Image/RegionGrowWithMultipleThresholds.h index 5156b33..409b18b 100644 --- a/lib/fpa/Image/RegionGrowWithMultipleThresholds.h +++ b/lib/fpa/Image/RegionGrowWithMultipleThresholds.h @@ -40,8 +40,8 @@ namespace fpa itkNewMacro( Self ); itkTypeMacro( RegionGrowWithMultipleThresholds, RegionGrow ); - itkGetConstMacro( DerivativeThreshold, double ); - itkSetMacro( DerivativeThreshold, double ); + itkGetConstMacro( DifferenceThreshold, double ); + itkSetMacro( DifferenceThreshold, double ); public: void AddThreshold( const TPixel& v ); @@ -59,6 +59,8 @@ namespace fpa virtual void _AfterMainLoop( ); virtual void _AfterLoop( ); virtual bool _UpdateResult( _TNode& n ); + virtual void _Mark( const _TNode& n ); + virtual bool _CheckStopCondition( ); private: RegionGrowWithMultipleThresholds( const Self& ); // Not impl. @@ -66,8 +68,12 @@ namespace fpa protected: TThresholds m_Thresholds; - double m_DerivativeThreshold; + double m_DifferenceThreshold; THistogram m_Histogram; + unsigned long m_TotalCount; + double m_LastDiff; + bool m_StopForced; + TPixel m_StopThreshold; }; } // ecapseman diff --git a/lib/fpa/Image/RegionGrowWithMultipleThresholds.hxx b/lib/fpa/Image/RegionGrowWithMultipleThresholds.hxx index 0cb1b66..387a833 100644 --- a/lib/fpa/Image/RegionGrowWithMultipleThresholds.hxx +++ b/lib/fpa/Image/RegionGrowWithMultipleThresholds.hxx @@ -2,6 +2,7 @@ #define __FPA__IMAGE__REGIONGROWWITHMULTIPLETHRESHOLDS__HXX__ #include +#include // ------------------------------------------------------------------------- template< class I > @@ -17,9 +18,8 @@ template< class I > void fpa::Image::RegionGrowWithMultipleThresholds< I >:: AddThresholds( const TPixel& t0, const TPixel& t1, const unsigned int& s ) { - double o = ( double( t1 ) - double( t0 ) ) / ( double( s ) - double( 1 ) ); - for( unsigned int i = 0; i < s; i++ ) - this->AddThreshold( TPixel( ( double( i ) * o ) + double( t0 ) ) ); + for( TPixel t = t0; t <= t1; t += TPixel( s ) ) + this->AddThreshold( t ); } // ------------------------------------------------------------------------- @@ -27,7 +27,11 @@ template< class I > fpa::Image::RegionGrowWithMultipleThresholds< I >:: RegionGrowWithMultipleThresholds( ) : Superclass( ), - m_DerivativeThreshold( double( 3 ) ) + m_DifferenceThreshold( double( 3 ) ), + m_TotalCount( 0 ), + m_LastDiff( double( 0 ) ), + m_StopForced( false ), + m_StopThreshold( TPixel( 0 ) ) { } @@ -45,10 +49,17 @@ _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 = *tIt; - for( ++tIt; tIt != this->m_Thresholds.end( ); ++tIt ) + TPixel min_thr = std::numeric_limits< TPixel >::min( ); + for( ; tIt != this->m_Thresholds.end( ); ++tIt ) { typename TFunction::Pointer function = TFunction::New( ); function->SetInputImage( img ); @@ -58,9 +69,48 @@ _BeforeMainLoop( ) } // 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 ); + 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( ) +{ + this->Superclass::_AfterMainLoop( ); +} + // ------------------------------------------------------------------------- template< class I > void fpa::Image::RegionGrowWithMultipleThresholds< I >:: @@ -70,36 +120,111 @@ _AfterLoop( ) { TFunction* f = dynamic_cast< TFunction* >( this->m_ActualFunction->GetPointer( ) ); - // std::cout << f->GetUpperThreshold( ) << " " << this->m_Marks.size( ) << std::endl; + 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 > -void fpa::Image::RegionGrowWithMultipleThresholds< I >:: -_AfterMainLoop( ) -{ - this->Superclass::_AfterMainLoop( ); -} - // ------------------------------------------------------------------------- template< class I > bool fpa::Image::RegionGrowWithMultipleThresholds< I >:: _UpdateResult( _TNode& n ) { - #error ACA VOY -> explicar esto!!!! cambiar salida por una especie de curva de nivel + // 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( ) ); - this->GetOutput( )->SetPixel( n.Vertex, f->GetUpperThreshold( ) ); + 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 ) + { + 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$