tools/clitkImageGradientMagnitudeGenericFilter.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 clitkImageGradientMagnitudeGenericFilter_txx
  19 #define clitkImageGradientMagnitudeGenericFilter_txx
  20
  21 /* =================================================
  22  * @file   clitkImageGradientMagnitudeGenericFilter.txx
  23  * @author Jef Vandemeulebroucke <jef@creatis.insa-lyon.fr>
  24  * @date   29 june 2009
  25  *
  26  * @brief
  27  *
  28  ===================================================*/
  29
  30 // itk include
  31 #include "itkGradientMagnitudeImageFilter.h"
  32 #include "itkGradientMagnitudeRecursiveGaussianImageFilter.h"
  33 #include "itkLabelStatisticsImageFilter.h"
  34 #include "itkMaskImageFilter.h"
  35 #include "itkMaskNegatedImageFilter.h"
  36 #include <clitkCommon.h>
  37
  38 namespace clitk
  39 {
  40
  41     //--------------------------------------------------------------------
  42     template<class args_info_type>
  43     ImageGradientMagnitudeGenericFilter<args_info_type>::ImageGradientMagnitudeGenericFilter():
  44         ImageToImageGenericFilter<Self>("ImageGradientMagnitude")
  45     {
  46         InitializeImageType<2>();
  47         InitializeImageType<3>();
  48     }
  49     //--------------------------------------------------------------------
  50
  51
  52     //--------------------------------------------------------------------
  53     template<class args_info_type>
  54     template<unsigned int Dim>
  55     void ImageGradientMagnitudeGenericFilter<args_info_type>::InitializeImageType()
  56     {
  57         ADD_DEFAULT_IMAGE_TYPES(Dim);
  58     }
  59     //--------------------------------------------------------------------
  60
  61
  62     //--------------------------------------------------------------------
  63     template<class args_info_type>
  64     void ImageGradientMagnitudeGenericFilter<args_info_type>::SetArgsInfo(const args_info_type & a)
  65     {
  66         mArgsInfo=a;
  67         this->SetIOVerbose(mArgsInfo.verbose_flag);
  68         if (mArgsInfo.imagetypes_flag) this->PrintAvailableImageTypes();
  69
  70         if (mArgsInfo.input_given) {
  71             this->SetInputFilename(mArgsInfo.input_arg);
  72         }
  73         if (mArgsInfo.output_given) {
  74             this->SetOutputFilename(mArgsInfo.output_arg);
  75         }
  76         //
  77         if (mArgsInfo.normalize_flag) {
  78             this->m_NormalizeOutput = 1;
  79         }
  80         else {
  81             this->m_NormalizeOutput = 0;
  82         }
  83         if (mArgsInfo.mask_given) {
  84           this->AddInputFilename(mArgsInfo.mask_arg);
  85         }
  86     }
  87     //--------------------------------------------------------------------
  88
  89     //--------------------------------------------------------------------
  90     // Update with the number of dimensions and the pixeltype
  91     //--------------------------------------------------------------------
  92     template<class args_info_type>
  93     template<class InputImageType>
  94     void
  95     ImageGradientMagnitudeGenericFilter<args_info_type>::UpdateWithInputImageType()
  96     {
  97         // Main filter
  98         typedef typename InputImageType::PixelType InputPixelType;
  99         typedef itk::Image<float, InputImageType::ImageDimension> OutputImageType;
 100         typedef itk::Image<unsigned char, OutputImageType::ImageDimension> MaskImageType;
 101
 102         // Reading input
 103         typename InputImageType::Pointer input = this->template GetInput<InputImageType>(0);
 104
 105         typename MaskImageType::Pointer mask = ITK_NULLPTR;
 106         if(mArgsInfo.mask_given) {
 107             mask = this->template GetInput<MaskImageType>(1);
 108         }
 109         else {
 110             mask = MaskImageType::New();
 111             mask->SetRegions(input->GetLargestPossibleRegion());
 112             mask->SetOrigin(input->GetOrigin());
 113             mask->SetSpacing(input->GetSpacing());
 114             mask->Allocate();
 115             mask->FillBuffer(1);
 116         }
 117
 118
 119         // Create output image
 120         typename OutputImageType::Pointer outputImage = OutputImageType::New();
 121         outputImage->SetRegions(input->GetLargestPossibleRegion());
 122         outputImage->SetOrigin(input->GetOrigin());
 123         outputImage->SetSpacing(input->GetSpacing());
 124         outputImage->Allocate();
 125         outputImage->FillBuffer(0.0);
 126         // Set output iterator
 127         typedef itk::ImageRegionIterator<OutputImageType> IteratorOutputType;
 128         IteratorOutputType ito = IteratorOutputType(outputImage, outputImage->GetLargestPossibleRegion());
 129
 130         // Filter
 131         typename OutputImageType::Pointer outputGradientFilter;
 132         if (mArgsInfo.gaussian_filter_flag == 0) {
 133         typedef itk::GradientMagnitudeImageFilter<InputImageType, OutputImageType> GradientMagnitudeImageFilterType;
 134         typename GradientMagnitudeImageFilterType::Pointer gradientFilter=GradientMagnitudeImageFilterType::New();
 135         gradientFilter->SetInput(input);
 136         gradientFilter->Update();
 137         outputGradientFilter = gradientFilter->GetOutput();
 138         }
 139         else {
 140             typedef itk::GradientMagnitudeRecursiveGaussianImageFilter<InputImageType, OutputImageType> GradientMagnitudeImageFilterType;
 141             typename GradientMagnitudeImageFilterType::Pointer gradientFilter=GradientMagnitudeImageFilterType::New();
 142             gradientFilter->SetInput(input);
 143             gradientFilter->Update();
 144             //std::cout<<"sigma value="<<gradientFilter->GetSigma()<<std::endl;
 145             outputGradientFilter = gradientFilter->GetOutput();
 146         }
 147         // Set iterator
 148         typedef itk::ImageRegionIterator<OutputImageType> IteratorType;
 149         IteratorType it(outputGradientFilter, outputGradientFilter->GetLargestPossibleRegion());
 150
 151         // Set mask iterator
 152         typedef itk::ImageRegionIterator<MaskImageType> IteratorMaskType;
 153         IteratorMaskType itm(mask, mask->GetLargestPossibleRegion());
 154
 155         //typedef itk::MinimumMaximumImageCalculator <OutputImageType> ImageCalculatorFilterType;
 156         //typename ImageCalculatorFilterType::Pointer imageCalculatorFilter = ImageCalculatorFilterType::New();
 157         //imageCalculatorFilter->SetImage(outputGradientFilter);
 158         //imageCalculatorFilter->Compute();
 159         typedef itk::LabelStatisticsImageFilter< OutputImageType, MaskImageType > LabelStatisticsImageFilterType;
 160         typename LabelStatisticsImageFilterType::Pointer labelStatisticsImageFilter = LabelStatisticsImageFilterType::New();
 161         labelStatisticsImageFilter->SetLabelInput( mask );
 162         labelStatisticsImageFilter->SetInput(outputGradientFilter);
 163         labelStatisticsImageFilter->Update();
 164
 165         //std::cout << "Number of labels: " << labelStatisticsImageFilter->GetNumberOfLabels() << std::endl;
 166
 167         float minImg = labelStatisticsImageFilter->GetMinimum(1);
 168         //std::cout << "minImg: " << minImg << std::endl;
 169         float maxImg = labelStatisticsImageFilter->GetMaximum(1);
 170         //std::cout << "maxImg: " << maxImg << std::endl;
 171
 172         it.GoToBegin();
 173         ito.GoToBegin();
 174         itm.GoToBegin();
 175
 176         while (!ito.IsAtEnd()) {
 177             if(m_NormalizeOutput && itm.Get() == 1) {
 178                 ito.Set(((float) it.Get() - minImg)/(maxImg-minImg));
 179             }
 180             if (m_NormalizeOutput == 0 && itm.Get() == 1) {
 181                 ito.Set((float) it.Get());
 182             }
 183             ++it;
 184             ++ito;
 185             ++itm;
 186         }
 187
 188         //typename OutputImageType::Pointer outputImage = outputGradientFilter;
 189         this->template SetNextOutput<OutputImageType>(outputImage);
 190     }
 191     //--------------------------------------------------------------------
 192
 193
 194 }//end clitk
 195
 196 #endif //#define clitkImageGradientMagnitudeGenericFilter_txx