X-Git-Url: https://git.creatis.insa-lyon.fr/pubgit/?a=blobdiff_plain;f=libs%2Ffpa%2FImage%2FMoriRegionGrow.hxx;h=5efdf65179c415022262bef81b57b23e19311e93;hb=1bde955a3f32330dcbbd2c190de75a8f77865de9;hp=737bf1f79d15acb515c1d9fb6f94a2f959c80bb3;hpb=49336779a8037ecd49fc0c4f6038c02636eb538b;p=FrontAlgorithms.git

diff --git a/libs/fpa/Image/MoriRegionGrow.hxx b/libs/fpa/Image/MoriRegionGrow.hxx
index 737bf1f..5efdf65 100644
--- a/libs/fpa/Image/MoriRegionGrow.hxx
+++ b/libs/fpa/Image/MoriRegionGrow.hxx
@@ -1,160 +1,252 @@
+// =========================================================================
+// @author Leonardo Florez Valencia
+// @email florez-l@javeriana.edu.co
+// =========================================================================
+
 #ifndef __fpa__Image__MoriRegionGrow__hxx__
 #define __fpa__Image__MoriRegionGrow__hxx__
 
+#include <queue>
+#include <itkProgressReporter.h>
+
 // -------------------------------------------------------------------------
-template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
+template< class _TInputImage, class _TOutputImage, class _TAuxiliaryPixel >
 typename
-fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
-TAuxImage*
-fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
+fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxiliaryPixel >::
+TAuxiliaryImage*
+fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxiliaryPixel >::
 GetAuxiliaryImage( )
 {
-  return( this->m_Helper->GetOutput( ) );
+  return(
+    dynamic_cast< TAuxiliaryImage* >(
+      this->itk::ProcessObject::GetOutput( 1 )
+      )
+    );
 }
 
 // -------------------------------------------------------------------------
-template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
+template< class _TInputImage, class _TOutputImage, class _TAuxiliaryPixel >
 const typename
-fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
-TAuxImage*
-fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
+fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxiliaryPixel >::
+TAuxiliaryImage*
+fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxiliaryPixel >::
 GetAuxiliaryImage( ) const
 {
-  return( this->m_Helper->GetOutput( ) );
+  return(
+    dynamic_cast< const TAuxiliaryImage* >(
+      this->itk::ProcessObject::GetOutput( 1 )
+      )
+    );
 }
 
 // -------------------------------------------------------------------------
-template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
-typename
-fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
-TInputPixel
-fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
-GetLower( ) const
+template< class _TInputImage, class _TOutputImage, class _TAuxiliaryPixel >
+void
+fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxiliaryPixel >::
+SetThresholdRange(
+  const TInputPixel& lower, const TInputPixel& upper,
+  const TInputPixel& delta
+  )
 {
-  return( this->m_Helper->GetLower( ) );
+  this->SetLowerThreshold( lower );
+  this->SetUpperThreshold( upper );
+  this->SetDeltaThreshold( delta );
 }
 
 // -------------------------------------------------------------------------
-template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
-typename
-fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
-TInputPixel
-fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
-GetUpper( ) const
-{
-  return( this->m_Helper->GetUpper( ) );
-}
-
-// -------------------------------------------------------------------------
-template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
-typename
-fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
-TInputPixel
-fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
-GetStep( ) const
-{
-  return( this->m_Helper->GetStep( ) );
-}
-
-// -------------------------------------------------------------------------
-template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
-typename
-fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
-TOutputPixel
-fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
-GetInsideValue( ) const
-{
-  return( this->m_Threshold->GetInsideValue( ) );
-}
-
-// -------------------------------------------------------------------------
-template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
-typename
-fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
-TOutputPixel
-fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
-GetOutsideValue( ) const
-{
-  return( this->m_Threshold->GetInsideValue( ) );
-}
-
-// -------------------------------------------------------------------------
-template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
-void fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
-SetLower( const TInputPixel& v )
-{
-  this->m_Helper->SetLower( v );
-  this->Modified( );
-}
-
-// -------------------------------------------------------------------------
-template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
-void fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
-SetUpper( const TInputPixel& v )
-{
-  this->m_Helper->SetUpper( v );
-  this->Modified( );
-}
-
-// -------------------------------------------------------------------------
-template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
-void fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
-SetStep( const TInputPixel& v )
-{
-  this->m_Helper->SetStep( v );
-  this->Modified( );
-}
-
-// -------------------------------------------------------------------------
-template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
-void fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
-SetInsideValue( const TOutputPixel& v )
+template< class _TInputImage, class _TOutputImage, class _TAuxiliaryPixel >
+fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxiliaryPixel >::
+MoriRegionGrow( )
+  : Superclass( )
 {
-  this->m_Threshold->SetInsideValue( v );
-  this->Modified( );
+  this->m_InsideValue = TOutputPixel( 1 );
+  this->m_OutsideValue = TOutputPixel( 0 );
+  this->m_UpperThreshold = std::numeric_limits< TInputPixel >::max( );
+  if( std::numeric_limits< TInputPixel >::is_integer )
+    this->m_LowerThreshold = std::numeric_limits< TInputPixel >::min( );
+  else
+    this->m_LowerThreshold = -this->m_UpperThreshold;
+  this->m_DeltaThreshold = TInputPixel( 1 );
+
+  this->SetNumberOfRequiredOutputs( 2 );
+  this->SetNthOutput( 0, TOutputImage::New( ) );
+  this->SetNthOutput( 1, TAuxiliaryImage::New( ) );
+
+  this->m_ThresholdFilter = TThresholdFilter::New( );
 }
 
 // -------------------------------------------------------------------------
-template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
-void fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
-SetOutsideValue( const TOutputPixel& v )
+template< class _TInputImage, class _TOutputImage, class _TAuxiliaryPixel >
+fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxiliaryPixel >::
+~MoriRegionGrow( )
 {
-  this->m_Threshold->SetOutsideValue( v );
-  this->Modified( );
 }
 
 // -------------------------------------------------------------------------
-template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
-fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
-MoriRegionGrow( )
-  : Superclass( )
+template< class _TInputImage, class _TOutputImage, class _TAuxiliaryPixel >
+void
+fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxiliaryPixel >::
+GenerateInputRequestedRegion( )
 {
-  this->m_Helper = THelper::New( );
-  this->m_Threshold = TThreshold::New( );
+  this->Superclass::GenerateInputRequestedRegion( );
+  if( this->GetInput( ) )
+  {
+    TInputImage* input =
+      const_cast< TInputImage* >( this->GetInput( ) );
+    input->SetRequestedRegionToLargestPossibleRegion( );
+
+  } // fi
 }
 
 // -------------------------------------------------------------------------
-template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
-fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
-~MoriRegionGrow( )
+template< class _TInputImage, class _TOutputImage, class _TAuxiliaryPixel >
+void
+fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxiliaryPixel >::
+EnlargeOutputRequestedRegion( itk::DataObject* output )
 {
+  this->Superclass::EnlargeOutputRequestedRegion( output );
+  output->SetRequestedRegionToLargestPossibleRegion( );
 }
 
 // -------------------------------------------------------------------------
-template< class _TInputImage, class _TOutputImage, class _TAuxPixel >
-void fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxPixel >::
+template< class _TInputImage, class _TOutputImage, class _TAuxiliaryPixel >
+void fpa::Image::MoriRegionGrow< _TInputImage, _TOutputImage, _TAuxiliaryPixel >::
 GenerateData( )
 {
-  this->m_Helper->ClearSeeds( );
-  this->m_Helper->AddSeed( this->m_Seed, 0 );
-  this->m_Helper->SetInput( this->GetInput( ) );
-  this->m_Helper->Update( );
-
-  this->m_Threshold->SetInput( this->m_Helper->GetOutput( ) );
-  this->m_Threshold->SetLowerThreshold( std::numeric_limits< _TAuxPixel >::min( ) );
-  this->m_Threshold->SetUpperThreshold( this->m_Helper->GetOptimumThreshold( ) );
-  this->m_Threshold->Update( );
-  this->GetOutput( )->Graft( this->m_Threshold->GetOutput( ) );
+  const TInputImage* input = this->GetInput( );
+  TAuxiliaryImage* auxiliary = this->GetAuxiliaryImage( );
+  TRegion region = input->GetRequestedRegion( );
+
+  // Prepare progress counter
+  double prog_size = double( this->m_UpperThreshold - this->m_LowerThreshold );
+  prog_size /= double( this->m_DeltaThreshold );
+  itk::ProgressReporter progress( this, 0, prog_size, 100, 0, 0.7 );
+
+  // Configure outputs
+  auxiliary->SetBufferedRegion( region );
+  auxiliary->Allocate( );
+  auxiliary->FillBuffer( 0 );
+
+  // Init marks
+  typedef itk::Image< bool, TInputImage::ImageDimension > _TMarks;
+  typename _TMarks::Pointer marks = _TMarks::New( );
+  marks->CopyInformation( input );
+  marks->SetRequestedRegion( region );
+  marks->SetBufferedRegion( input->GetBufferedRegion( ) );
+  marks->Allocate( );
+  marks->FillBuffer( false );
+
+  // Prepare queues
+  std::queue< TIndex > queues[ 2 ];
+  queues[ 0 ].push( this->m_Seed );
+  unsigned int cur_queue = 0;
+  unsigned int aux_queue = 1;
+
+  // Incremental strategy
+  TInputPixel upper = this->m_LowerThreshold + this->m_DeltaThreshold;
+  this->m_Curve.clear( );
+  unsigned long count = 0;
+  while( upper <= this->m_UpperThreshold )
+  {
+    // Main growing strategy
+    while( queues[ cur_queue ].size( ) > 0 )
+    {
+      // Get next candidate
+      TIndex node = queues[ cur_queue ].front( );
+      queues[ cur_queue ].pop( );
+      if( marks->GetPixel( node ) )
+        continue;
+      marks->SetPixel( node, true );
+
+      // Apply inclusion predicate
+      TInputPixel value = input->GetPixel( node );
+      bool in = ( ( this->m_LowerThreshold < value ) && ( value < upper ) );
+      if( !in )
+      {
+        if( value < this->m_UpperThreshold )
+          queues[ aux_queue ].push( node );
+        marks->SetPixel( node, false );
+        continue;
+
+      } // fi
+
+      // Ok, pixel lays inside region
+      auxiliary->SetPixel( node, this->m_Curve.size( ) + 1 );
+      count++;
+
+      // Add neighborhood
+      for( unsigned int d = 0; d < TInputImage::ImageDimension; ++d )
+      {
+        TIndex neigh = node;
+        for( int i = -1; i <= 1; i += 2 )
+        {
+          neigh[ d ] = node[ d ] + i;
+          if( region.IsInside( neigh ) )
+            queues[ cur_queue ].push( neigh );
+
+        } // rof
+
+      } // rof
+
+    } // elihw
+
+    // Update curve
+    if( this->m_Curve.size( ) > 0 )
+    {
+      if( this->m_Curve.back( ).YValue < count )
+        this->m_Curve.push_back( TCurveData( upper, count ) );
+      if( this->m_Curve.size( ) > 2 )
+      {
+        long j = this->m_Curve.size( ) - 2;
+        double yp = double( this->m_Curve[ j + 1 ].YValue );
+        double yn = double( this->m_Curve[ j - 1 ].YValue );
+        double xp = double( this->m_Curve[ j + 1 ].XValue );
+        double xn = double( this->m_Curve[ j - 1 ].XValue );
+        this->m_Curve[ j ].Diff1 = ( yp - yn ) / ( xp - xn );
+
+      } // fi
+    }
+    else
+      this->m_Curve.push_back( TCurveData( upper, count ) );
+
+    // Update queue
+    cur_queue = aux_queue;
+    aux_queue = ( cur_queue + 1 ) % 2;
+
+    // Update threshold
+    upper += this->m_DeltaThreshold;
+    progress.CompletedPixel( );
+
+  } // elihw
+
+  // Compute optimum threshold
+  double dmax = -std::numeric_limits< double >::max( );
+  int jmax = 0;
+  for( int j = 1; j < this->m_Curve.size( ) - 1; ++j )
+  {
+    double dp = this->m_Curve[ j + 1 ].Diff1;
+    double dn = this->m_Curve[ j - 1 ].Diff1;
+    double xp = double( this->m_Curve[ j + 1 ].XValue );
+    double xn = double( this->m_Curve[ j - 1 ].XValue );
+    double d2 = ( dp - dn ) / ( xp - xn );
+    if( d2 > dmax )
+    {
+      dmax = d2;
+      jmax = j;
+
+    } // fi
+
+  } // rof
+  this->m_OptimumThreshold = this->m_Curve[ jmax ].XValue;
+
+  // Final threshold
+  this->m_ThresholdFilter->SetInput( auxiliary );
+  this->m_ThresholdFilter->SetInsideValue( this->m_InsideValue );
+  this->m_ThresholdFilter->SetOutsideValue( this->m_OutsideValue );
+  this->m_ThresholdFilter->SetLowerThreshold( 1 );
+  this->m_ThresholdFilter->SetUpperThreshold( jmax );
+  this->m_ThresholdFilter->Update( );
+  this->GetOutput( )->Graft( this->m_ThresholdFilter->GetOutput( ) );
 }
 
 #endif // __fpa__Image__MoriRegionGrow__hxx__