#define __FPA__IMAGE__REGIONGROWWITHMULTIPLETHRESHOLDS__HXX__
#include <limits>
-#include <fpa/Image/Functors/RegionGrowThresholdFunction.h>
+#include <itkBinaryThresholdImageFilter.h>
+#include <itkConstNeighborhoodIterator.h>
// -------------------------------------------------------------------------
template< class I >
void fpa::Image::RegionGrowWithMultipleThresholds< I >::
AddThreshold( const TPixel& v )
{
- typedef
- fpa::Image::Functors::RegionGrowThresholdFunction< I >
- TFunction;
- typename TFunction::Pointer function = TFunction::New( );
-
- if( this->GetNumberOfMembershipFunctions( ) > 0 )
- {
- }
- else
- function->SetLowerThreshold( std::numeric_limits< TPixel >::min( ) );
- std::cout
- << typeid( TPixel ).name( ) << " <<<----->>> "
- << function->GetLowerThreshold( )
- << std::endl;
- std::exit( 1 );
- function->SetUpperThreshold( v );
- this->AddMembershipFunction( function );
-
- /* TODO
- this->m_Histogram[ v ] = 0;
-
- TFunction* function =
- dynamic_cast< TFunction* >( this->GetMembershipFunction( ) );
- if( function != NULL )
- {
- function->SetLowerThreshold( this->m_Histogram.begin( )->first );
-
- } // fi
- this->Modified( );
- */
+ this->m_Thresholds.insert( v );
+ this->Modified( );
}
// -------------------------------------------------------------------------
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 );
}
// -------------------------------------------------------------------------
fpa::Image::RegionGrowWithMultipleThresholds< I >::
RegionGrowWithMultipleThresholds( )
: Superclass( ),
- m_DerivativeThreshold( double( 3 ) )
+ 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 ) )
{
- /* TODO
- typedef
- fpa::Image::Functors::RegionGrowThresholdFunction< I >
- TFunction;
- typename TFunction::Pointer function = TFunction::New( );
- this->SetMembershipFunction( function );
- */
}
// -------------------------------------------------------------------------
// -------------------------------------------------------------------------
template< class I >
-bool fpa::Image::RegionGrowWithMultipleThresholds< I >::
-_UpdateResult( _TNode& n )
+void fpa::Image::RegionGrowWithMultipleThresholds< I >::
+_BeforeMainLoop( )
{
- bool ret = this->Superclass::_UpdateResult( n );
+ const I* img = this->GetInput( );
- /* TODO
- if( ret )
- {
- TPixel v = TPixel( this->_Cost( n.Vertex, n.Vertex ) );
-
- typename THistogram::reverse_iterator hIt = this->m_Histogram.rbegin( );
- while( hIt != this->m_Histogram.rend( ) )
- {
- if( v <= hIt->first )
- {
- hIt->second += 1;
- hIt++;
- }
- else
- hIt = this->m_Histogram.rend( );
-
- } // elihw
- this->GetOutput( )->SetPixel( n.Vertex, v );
-
- } // fi
- */
- return( ret );
+ // 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 >::
-_BeforeLoop( )
+_AfterMainLoop( )
{
- std::cout << "**1" << std::endl;
- const I* img = this->GetInput( );
- std::cout << "**2" << std::endl;
- typename TFunctions::iterator fIt = this->m_Functions.begin( );
- for( ; fIt != this->m_Functions.end( ); ++fIt )
- {
- TMembershipFunction* f =
- dynamic_cast< TMembershipFunction* >( fIt->GetPointer( ) );
- if( f != NULL )
- f->SetInputImage( this->GetInput( ) );
+ typedef itk::BinaryThresholdImageFilter< I, I > _TBinFilter;
- } // rof
- this->Superclass::_BeforeLoop( );
+ // 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( );
}
// -------------------------------------------------------------------------
void fpa::Image::RegionGrowWithMultipleThresholds< I >::
_AfterLoop( )
{
- /*
- typename THistogram::iterator prevIt = this->m_Histogram.begin( );
- typename THistogram::iterator currIt = prevIt; currIt++;
- typename THistogram::iterator nextIt = currIt; nextIt++;
- double prev_d1 = double( 0 );
- for( ; nextIt != this->m_Histogram.end( ); ++prevIt, ++currIt, ++nextIt )
+ 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 )
{
- double d1 = double( nextIt->second ) - double( prevIt->second );
- d1 /= double( nextIt->first ) - double( prevIt->first );
-
- std::cout
- << currIt->first << " "
- << currIt->second << " "
- << d1 << " "
- << ( d1 - prev_d1 ) << std::endl;
+ 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( );
+}
- prev_d1 = d1;
+// -------------------------------------------------------------------------
+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( ) );
- } // rof
- */
+ } // fi
+ return( ret );
+}
- /*
- double prev = double( 0 );
- while( hIt != this->m_Histogram.end( ) )
+// -------------------------------------------------------------------------
+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 )
{
- double curr = double( hIt->second );
- double deri = curr - prev;
-
- if( hIt != this->m_Histogram.begin( ) )
- std::cout << hIt->first << " " << curr << " " << deri << std::endl;
-
- prev = curr;
- hIt++;
+ // 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
+}
- } // rof
- */
- this->Superclass::_AfterLoop( );
+// -------------------------------------------------------------------------
+template< class I >
+bool fpa::Image::RegionGrowWithMultipleThresholds< I >::
+_CheckStopCondition( )
+{
+ return( this->m_StopForced );
}
#endif // __FPA__IMAGE__REGIONGROWWITHMULTIPLETHRESHOLDS__HXX__
git://git.creatis.insa-lyon.fr / FrontAlgorithms.git / blobdiff