segmentation/clitkFillMaskFilter.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
  19 #ifndef CLITKFILLMASKFILTER_TXX
  20 #define CLITKFILLMASKFILTER_TXX
  21
  22 //--------------------------------------------------------------------
  23 template <class ImageType>
  24 clitk::FillMaskFilter<ImageType>::FillMaskFilter():
  25   itk::ImageToImageFilter<ImageType, ImageType>()
  26 {
  27   // Default global options
  28   this->SetNumberOfRequiredInputs(1);
  29   ResetDirections();
  30 }
  31 //--------------------------------------------------------------------
  32
  33
  34 //--------------------------------------------------------------------
  35 template <class ImageType>
  36 void
  37 clitk::FillMaskFilter<ImageType>::
  38 AddDirection(int dir)
  39 {
  40   m_Directions.push_back(dir);
  41 }
  42 //--------------------------------------------------------------------
  43
  44
  45 //--------------------------------------------------------------------
  46 template <class ImageType>
  47 void
  48 clitk::FillMaskFilter<ImageType>::
  49 ResetDirections()
  50 {
  51   m_Directions.clear();
  52 }
  53 //--------------------------------------------------------------------
  54
  55
  56 //--------------------------------------------------------------------
  57 template <class ImageType>
  58 template<class ArgsInfoType>
  59 void
  60 clitk::FillMaskFilter<ImageType>::
  61 SetOptionsFromArgsInfo(ArgsInfoType & mArgsInfo)
  62 {
  63   ResetDirections();
  64   for(uint i=0; i<mArgsInfo.dir_given; i++) {
  65     AddDirection(mArgsInfo.dir_arg[i]);
  66   }
  67 }
  68 //--------------------------------------------------------------------
  69
  70
  71 //--------------------------------------------------------------------
  72 template <class ImageType>
  73 template<class ArgsInfoType>
  74 void
  75 clitk::FillMaskFilter<ImageType>::
  76 SetOptionsToArgsInfo(ArgsInfoType & mArgsInfo)
  77 {
  78   mArgsInfo.dir_arg = new int[m_Directions.size()];
  79   mArgsInfo.dir_given = m_Directions.size();
  80   for(uint i=0; i<mArgsInfo.dir_given; i++) {
  81     mArgsInfo.dir_arg[i] = m_Directions[i];
  82   }
  83 }
  84 //--------------------------------------------------------------------
  85
  86
  87 //--------------------------------------------------------------------
  88 template <class ImageType>
  89 void
  90 clitk::FillMaskFilter<ImageType>::
  91 GenerateData()
  92 {
  93   InputImagePointer input  = dynamic_cast<ImageType*>(itk::ProcessObject::GetInput(0));
  94
  95   // Cast
  96   typedef int InternalPixelType;
  97   typedef itk::Image<InternalPixelType, ImageType::ImageDimension> InternalImageType;
  98   typedef itk::CastImageFilter<ImageType,InternalImageType> InputCastImageFilterType;
  99   typename InputCastImageFilterType::Pointer inputCaster = InputCastImageFilterType::New();
 100   inputCaster->SetInput(input);
 101   inputCaster->Update();
 102
 103   //----------------------------------------
 104   // Loop over directions
 105   //----------------------------------------
 106   typename InternalImageType::Pointer output = inputCaster->GetOutput();
 107   for (unsigned int i=0; i<m_Directions.size();i++) {
 108     //----------------------------------------
 109     // Fill the holes of a mask in 2D
 110     //----------------------------------------
 111
 112     // We define the region to be extracted.
 113     typedef  itk::Image<InternalPixelType, ImageType::ImageDimension-1> ImageSliceType;
 114     typedef  itk::Image<InternalPixelType, ImageType::ImageDimension-1> MaskSliceType;
 115     typename InternalImageType::RegionType region3D = input->GetLargestPossibleRegion();
 116     typename InternalImageType::RegionType::SizeType size3D = region3D.GetSize();
 117     typename InternalImageType::RegionType::SizeType size2D = size3D;
 118     size2D[m_Directions[i]]=0;
 119     typename InternalImageType::IndexType start2D;
 120     start2D.Fill(0);
 121     typename InternalImageType::RegionType desiredRegion;
 122     desiredRegion.SetSize( size2D );
 123     desiredRegion.SetIndex( start2D );
 124
 125     // Extract and Join
 126     typedef itk::ExtractImageFilter<InternalImageType, ImageSliceType> ExtractImageFilterType;
 127     typedef itk::JoinSeriesImageFilter<ImageSliceType, InternalImageType> JoinSeriesFilterType;
 128     typename JoinSeriesFilterType::Pointer joinFilter=JoinSeriesFilterType::New();
 129     joinFilter->SetSpacing(input->GetSpacing()[m_Directions[i]]);
 130
 131     //----------------------------------------
 132     // Run over the sliceIndexs
 133     // ----------------------------------------
 134     for(unsigned int sliceIndex=0; sliceIndex <size3D[m_Directions[i]]; sliceIndex++)
 135       {
 136         //----------------------------------------
 137         // Extract mask sliceIndex
 138         //----------------------------------------
 139         typename ExtractImageFilterType::Pointer extractFilter=ExtractImageFilterType::New();
 140         extractFilter->SetInput(output);
 141         start2D[m_Directions[i]]=sliceIndex;
 142         desiredRegion.SetIndex( start2D );
 143         extractFilter->SetExtractionRegion( desiredRegion );
 144 #if ITK_VERSION_MAJOR == 4
 145         extractFilter->SetDirectionCollapseToSubmatrix();
 146 #endif
 147         extractFilter->Update( );
 148         typename ImageSliceType::Pointer slice= extractFilter->GetOutput();
 149
 150         // Binarize the image (Before: OBJECT!=0, rest=0, After: object=1, rest=0 )
 151         typedef itk::BinaryThresholdImageFilter<ImageSliceType,ImageSliceType> BinarizeFilterType;
 152         typename BinarizeFilterType::Pointer binarizeFilter=BinarizeFilterType::New();
 153         binarizeFilter->SetInput(slice);
 154         binarizeFilter->SetUpperThreshold(0);
 155         binarizeFilter->SetOutsideValue(0);
 156         binarizeFilter->SetInsideValue(1);
 157         // writeImage<ImageSliceType>(binarizeFilter->GetOutput(),"/home/jef/tmp/input.mhd");
 158
 159         // Perform connected labelling on the slice (body+air=0 )
 160         typedef itk::ConnectedComponentImageFilter<ImageSliceType, ImageSliceType> ConnectFilterType;
 161         typename ConnectFilterType::Pointer connectFilter=ConnectFilterType::New();
 162         connectFilter->SetInput(binarizeFilter->GetOutput());
 163         connectFilter->SetBackgroundValue(0);
 164         connectFilter->SetFullyConnected(false);
 165         //connectFilter->Update();
 166         //writeImage<ImageSliceType>(connectFilter->GetOutput(),"/home/jef/tmp/connect.mhd");
 167
 168         // Sort the labels
 169         typedef itk::RelabelComponentImageFilter<ImageSliceType, ImageSliceType> RelabelFilterType;
 170         typename RelabelFilterType::Pointer relabelFilter=RelabelFilterType::New();
 171         relabelFilter->SetInput(connectFilter->GetOutput());
 172         //relabelFilter->Update();
 173         //writeImage<ImageSliceType>(relabelFilter->GetOutput(),"/home/jef/tmp/label.mhd");
 174
 175         // Keep the first
 176         typedef itk::ThresholdImageFilter<ImageSliceType> ThresholdFilterType;
 177         typename ThresholdFilterType::Pointer thresholdFilter=ThresholdFilterType::New();
 178         thresholdFilter->SetInput(relabelFilter->GetOutput());
 179         thresholdFilter->SetUpper(1);
 180         thresholdFilter->SetOutsideValue(0);
 181         // thresholdFilter->Update();
 182         // writeImage<ImageSliceType>(thresholdFilter->GetOutput(),"/home/jef/tmp/bin.mhd");
 183
 184         // Invert the labels (lung 1, rest 0)
 185         typename BinarizeFilterType::Pointer switchFilter=BinarizeFilterType::New();
 186         switchFilter->SetInput(thresholdFilter->GetOutput());
 187         switchFilter->SetUpperThreshold(0);
 188         switchFilter->SetOutsideValue(0);
 189         switchFilter->SetInsideValue(1);
 190         switchFilter->Update();
 191         //writeImage<ImageSliceType>(switchFilter->GetOutput(),"/home/jef/tmp/inv_bin.mhd");
 192
 193         //Join
 194         joinFilter->SetInput( sliceIndex, switchFilter->GetOutput());
 195       }
 196
 197     // Join to a 3D image
 198     joinFilter->Update();
 199
 200     // Permute the axes to reset to orientation
 201     typedef itk::PermuteAxesImageFilter<InternalImageType> PermuteFilterType;
 202     typename PermuteFilterType::Pointer permuteFilter=PermuteFilterType::New();
 203     permuteFilter->SetInput(joinFilter->GetOutput());
 204     typename PermuteFilterType::PermuteOrderArrayType order;
 205     order[m_Directions[i]]=2;
 206     if( m_Directions[i]==2)
 207       {
 208         order[0]=0;
 209         order[1]=1;
 210       }
 211     else if ( m_Directions[i]==1)
 212       {
 213         order[0]=0;
 214         order[2]=1;
 215       }
 216     else if (m_Directions[i]==0)
 217       {
 218         order[1]=0;
 219         order[2]=1;
 220       }
 221     permuteFilter->SetOrder(order);
 222     permuteFilter->Update();
 223     output = permuteFilter->GetOutput();
 224
 225     // Set the image direction to the input one
 226     output->SetDirection(input->GetDirection());
 227     output->SetOrigin(input->GetOrigin());
 228   }
 229
 230   // Cast
 231   typedef itk::CastImageFilter<InternalImageType,ImageType> OutputCastImageFilterType;
 232   typename OutputCastImageFilterType::Pointer outputCaster =OutputCastImageFilterType::New();
 233   outputCaster->SetInput(output);
 234   outputCaster->Update();
 235
 236   // Output
 237   this->GraftOutput(outputCaster->GetOutput());
 238 }
 239 //--------------------------------------------------------------------
 240
 241 #endif //#define CLITKFILLMASKFILTER_TXX