itk/clitkLocallyAdaptiveBinaryThresholdImageFunction.txx

   1 /*=========================================================================
   2   Program:   vv                     http://www.creatis.insa-lyon.fr/rio/vv
   3
   4   Authors belong to:
   5   - University of LYON              http://www.universite-lyon.fr/
   6   - Léon Bérard cancer center       http://www.centreleonberard.fr
   7   - CREATIS CNRS laboratory         http://www.creatis.insa-lyon.fr
   8
   9   This software is distributed WITHOUT ANY WARRANTY; without even
  10   the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
  11   PURPOSE.  See the copyright notices for more information.
  12
  13   It is distributed under dual licence
  14
  15   - BSD        See included LICENSE.txt file
  16   - CeCILL-B   http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
  17 ===========================================================================**/
  18 #ifndef __clitkLocallyAdaptiveBinaryThresholdImageFunction_txx
  19 #define __clitkLocallyAdaptiveBinaryThresholdImageFunction_txx
  20 #include "clitkLocallyAdaptiveBinaryThresholdImageFunction.h"
  21 #include "itkNumericTraits.h"
  22 #include "itkConstNeighborhoodIterator.h"
  23
  24 namespace clitk
  25 {
  26
  27   /**
  28    * Constructor
  29    */
  30   template <class TInputImage, class TCoordRep>
  31   LocallyAdaptiveBinaryThresholdImageFunction<TInputImage,TCoordRep>
  32   ::LocallyAdaptiveBinaryThresholdImageFunction()
  33   {
  34     m_Radius.Fill(1);
  35     m_LowerBorderIsGiven=true;
  36     m_UpperBorderIsGiven=true;
  37     m_MaximumSDIsGiven=true;
  38     m_Multiplier=1.0;
  39     m_MaximumSD=100.0;
  40   }
  41
  42
  43   /**
  44    *
  45    */
  46   template <class TInputImage, class TCoordRep>
  47   void
  48   LocallyAdaptiveBinaryThresholdImageFunction<TInputImage,TCoordRep>
  49   ::PrintSelf(std::ostream& os, itk::Indent indent) const
  50   {
  51     this->Superclass::PrintSelf(os,indent);
  52
  53     os << indent << "Radius: " << m_Radius << std::endl;
  54   }
  55
  56
  57   /**
  58    *
  59    */
  60   template <class TInputImage, class TCoordRep>
  61   bool
  62   LocallyAdaptiveBinaryThresholdImageFunction<TInputImage,TCoordRep>
  63   ::EvaluateAtIndex(const IndexType& index) const
  64   {
  65
  66     if( !this->GetInputImage() )
  67       {
  68         return ( false );
  69       }
  70
  71     if ( !this->IsInsideBuffer( index ) )
  72       {
  73         return ( false );
  74       }
  75
  76     // Create an N-d neighborhood kernel, using a zeroflux boundary condition
  77     itk::ConstNeighborhoodIterator<InputImageType>
  78       it(m_Radius, this->GetInputImage(), this->GetInputImage()->GetBufferedRegion());
  79
  80     // Set the iterator at the desired location
  81     it.SetLocation(index);
  82     PixelType centerValue = it.GetPixel(0);
  83     PixelType currentvalue;
  84     bool isInside=true;
  85
  86     // Walk the neighborhood for the mean and SD
  87     const unsigned int size = it.Size();
  88     typename   itk::NumericTraits<PixelType>::RealType mean=0;
  89     typename   itk::NumericTraits<PixelType>::RealType sd=0;
  90     for (unsigned int i = 1; i < size; ++i)
  91       {
  92         currentvalue=it.GetPixel(i);
  93         mean+=currentvalue;
  94         sd+=currentvalue*currentvalue;
  95       }
  96     mean/=( typename itk::NumericTraits<PixelType>::RealType) size-1.;
  97     sd= sqrt( (sd /(typename itk::NumericTraits<PixelType>::RealType)size-1.) - (mean*mean) );
  98
  99     // Verify fixed borders
 100     if (this->GetLower() > centerValue || centerValue > this->GetUpper())
 101       isInside = false;
 102
 103     // Verify lower adaptive borders
 104     if( (m_LowerBorderIsGiven)  &&  (centerValue <  (    mean - (typename  itk::NumericTraits<PixelType>::RealType) ( m_Multiplier*sd) )  )  )
 105       isInside = false;
 106
 107     // Verify upper adaptive border
 108     if ( (m_UpperBorderIsGiven)  &&  (centerValue >  (    mean + (typename  itk::NumericTraits<PixelType>::RealType) ( m_Multiplier*sd) )  )  )
 109       isInside = false;
 110
 111     // Verify SD
 112     if ( (m_MaximumSDIsGiven)  &&  ( sd> m_MaximumSD) )
 113       isInside=false;
 114
 115 //     DD(centerValue);
 116 //     DD(mean);
 117 //     DD(sd);
 118 //     DD(isInside);
 119
 120     return ( isInside );
 121   }
 122
 123
 124 } // end namespace itk
 125
 126 #endif