segmentation/clitkExtractLungFilter.h

   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://oncora1.lyon.fnclcc.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 CLITKEXTRACTLUNGSFILTER_H
  20 #define CLITKEXTRACTLUNGSFILTER_H
  21
  22 // clitk
  23 #include "clitkFilterBase.h"
  24 #include "clitkDecomposeAndReconstructImageFilter.h"
  25 #include "clitkExplosionControlledThresholdConnectedImageFilter.h"
  26 #include "clitkSegmentationUtils.h"
  27 #include "clitkFilterWithAnatomicalFeatureDatabaseManagement.h"
  28 #include "tree.hh"
  29
  30 // itk
  31 #include "itkStatisticsImageFilter.h"
  32 #include "itkTreeContainer.h"
  33
  34 namespace clitk {
  35
  36   //--------------------------------------------------------------------
  37   /*
  38     Try to extract the Lung part of a thorax CT. Inspired by
  39     Rikxoort2009, Section IIA, MedPhys.
  40
  41     - First, all air besides lungs and thrachea is removed, by
  42     removing the second largest label of the firstLabelImage and
  43     setting the remainder to 0HU . This modified input is optimally
  44     thresholded (Otsu1979).
  45
  46     - Trachea and bronchi are grown from seeds in the top of the image
  47     by explosion controlled region growing, slightly dilated and
  48     removed from the second label image.
  49
  50     - Left and right lung are separated (if necessary) by erosion and
  51     reconstructed by conditional dilation.
  52
  53     - TRACHEA is available at the end
  54
  55     TODO ********** Remaining holes can be       filled afterwards (clitkFillMask).
  56
  57   */
  58   //--------------------------------------------------------------------
  59
  60
  61   //--------------------------------------------------------------------
  62
  63 class Bifurcation
  64 {
  65 public:
  66   typedef itk::Index<3> IndexType;
  67   typedef itk::Point<double, 3> PointType;
  68   typedef double PixelType;
  69   Bifurcation(IndexType _index, PixelType _l, PixelType _l1, PixelType _l2) {
  70     index = _index;
  71     _l = l;
  72     _l1 = l1;
  73     _l2 = l2;
  74   }
  75   IndexType index;
  76   PointType point;
  77   PixelType l;
  78   PixelType l1;
  79   PixelType l2;
  80   typedef itk::Index<3> NodeType;
  81   typedef tree<NodeType> TreeType;
  82   typedef TreeType::iterator TreeIterator;
  83   TreeIterator treeIter;
  84 };
  85   //--------------------------------------------------------------------
  86
  87
  88   //--------------------------------------------------------------------
  89   template <class TImageType, class TMaskImageType>
  90   class ITK_EXPORT ExtractLungFilter:
  91     public virtual clitk::FilterBase,
  92     public clitk::FilterWithAnatomicalFeatureDatabaseManagement,
  93     public itk::ImageToImageFilter<TImageType, TMaskImageType>
  94   {
  95
  96   public:
  97     /** Standard class typedefs. */
  98     typedef itk::ImageToImageFilter<TImageType, TMaskImageType> Superclass;
  99     typedef ExtractLungFilter              Self;
 100     typedef itk::SmartPointer<Self>        Pointer;
 101     typedef itk::SmartPointer<const Self>  ConstPointer;
 102
 103     /** Method for creation through the object factory. */
 104     itkNewMacro(Self);
 105
 106     /** Run-time type information (and related methods). */
 107     itkTypeMacro(ExtractLungFilter, ImageToImageFilter);
 108     FILTERBASE_INIT;
 109
 110     /** Some convenient typedefs */
 111     typedef TImageType                       ImageType;
 112     typedef typename ImageType::ConstPointer InputImageConstPointer;
 113     typedef typename ImageType::Pointer      InputImagePointer;
 114     typedef typename ImageType::RegionType   InputImageRegionType;
 115     typedef typename ImageType::PixelType    InputImagePixelType;
 116     typedef typename ImageType::SizeType     InputImageSizeType;
 117     typedef typename ImageType::IndexType    InputImageIndexType;
 118     typedef typename ImageType::PointType    InputImagePointType;
 119
 120     typedef TMaskImageType                       MaskImageType;
 121     typedef typename MaskImageType::ConstPointer MaskImageConstPointer;
 122     typedef typename MaskImageType::Pointer      MaskImagePointer;
 123     typedef typename MaskImageType::RegionType   MaskImageRegionType;
 124     typedef typename MaskImageType::PixelType    MaskImagePixelType;
 125     typedef typename MaskImageType::SizeType     MaskImageSizeType;
 126     typedef typename MaskImageType::IndexType    MaskImageIndexType;
 127     typedef typename MaskImageType::PointType    MaskImagePointType;
 128
 129     itkStaticConstMacro(ImageDimension, unsigned int, ImageType::ImageDimension);
 130     typedef int InternalPixelType;
 131     typedef itk::Image<InternalPixelType, ImageType::ImageDimension> InternalImageType;
 132     typedef typename InternalImageType::Pointer                      InternalImagePointer;
 133     typedef typename InternalImageType::IndexType                    InternalIndexType;
 134     typedef LabelizeParameters<InternalPixelType>                    LabelParamType;
 135
 136     typedef Bifurcation BifurcationType;
 137     typedef MaskImageIndexType NodeType;
 138     typedef tree<NodeType> TreeType;
 139     typedef typename TreeType::iterator TreeIterator;
 140
 141     /** Connect inputs */
 142     void SetInput(const ImageType * image);
 143     void SetInputPatientMask(MaskImageType * mask, MaskImagePixelType BG);
 144     itkSetMacro(PatientMaskBackgroundValue, MaskImagePixelType);
 145     itkGetConstMacro(PatientMaskBackgroundValue, MaskImagePixelType);
 146     GGO_DefineOption(patientBG, SetPatientMaskBackgroundValue, MaskImagePixelType);
 147
 148     // Set all options at a time
 149     template<class ArgsInfoType>
 150       void SetArgsInfo(ArgsInfoType arg);
 151
 152     // Get output (only availabe after update !)
 153     typename MaskImageType::Pointer GetTracheaImage() { return trachea; }
 154
 155     // Background / Foreground
 156     itkGetConstMacro(BackgroundValue, MaskImagePixelType);
 157     itkGetConstMacro(ForegroundValue, MaskImagePixelType);
 158
 159     // For common segmentation processes
 160     itkSetMacro(MinimalComponentSize, int);
 161     itkGetConstMacro(MinimalComponentSize, int);
 162     GGO_DefineOption(minSize, SetMinimalComponentSize, int);
 163
 164     // Step 1 options RemoveAir
 165     itkSetMacro(UpperThreshold, InputImagePixelType);
 166     itkGetConstMacro(UpperThreshold, InputImagePixelType);
 167     GGO_DefineOption(upper, SetUpperThreshold, InputImagePixelType);
 168
 169     itkSetMacro(NumberOfSlicesToSkipBeforeSearchingSeed, int);
 170     itkGetConstMacro(NumberOfSlicesToSkipBeforeSearchingSeed, int);
 171     GGO_DefineOption(skipslices, SetNumberOfSlicesToSkipBeforeSearchingSeed, int);
 172
 173     itkSetMacro(LowerThreshold, InputImagePixelType);
 174     itkGetConstMacro(LowerThreshold, InputImagePixelType);
 175     itkSetMacro(UseLowerThreshold, bool);
 176     itkGetConstMacro(UseLowerThreshold, bool);
 177     itkBooleanMacro(UseLowerThreshold);
 178     GGO_DefineOption_WithTest(lower, SetLowerThreshold, InputImagePixelType, UseLowerThreshold);
 179
 180     void SetLabelizeParameters1(LabelParamType * a) { m_LabelizeParameters1 = a; }
 181     itkGetConstMacro(LabelizeParameters1, LabelParamType*);
 182     GGO_DefineOption_LabelParam(1, SetLabelizeParameters1, LabelParamType);
 183
 184     // Step 2 options FindTrachea
 185     itkSetMacro(UpperThresholdForTrachea, InputImagePixelType);
 186     itkGetConstMacro(UpperThresholdForTrachea, InputImagePixelType);
 187     GGO_DefineOption(upperThresholdForTrachea, SetUpperThresholdForTrachea, InputImagePixelType);
 188
 189     itkSetMacro(MultiplierForTrachea, double);
 190     itkGetConstMacro(MultiplierForTrachea, double);
 191     GGO_DefineOption(multiplierForTrachea, SetMultiplierForTrachea, double);
 192
 193     itkSetMacro(ThresholdStepSizeForTrachea, InputImagePixelType);
 194     itkGetConstMacro(ThresholdStepSizeForTrachea, InputImagePixelType);
 195     GGO_DefineOption(thresholdStepSizeForTrachea, SetThresholdStepSizeForTrachea, InputImagePixelType);
 196
 197     void AddSeed(InternalIndexType s);
 198     std::vector<InternalIndexType> & GetSeeds() { return  m_Seeds; }
 199     GGO_DefineOption_Vector(seed, AddSeed, InternalIndexType, ImageType::ImageDimension, true);
 200
 201     // Step 3 options ExtractLung
 202     itkSetMacro(NumberOfHistogramBins, int);
 203     itkGetConstMacro(NumberOfHistogramBins, int);
 204     GGO_DefineOption(bins, SetNumberOfHistogramBins, int);
 205
 206     void SetLabelizeParameters2(LabelParamType* a) { m_LabelizeParameters2 = a; }
 207     itkGetConstMacro(LabelizeParameters2, LabelParamType*);
 208     GGO_DefineOption_LabelParam(2, SetLabelizeParameters2, LabelParamType);
 209
 210     // Step 4 options RemoveTrachea
 211     itkSetMacro(RadiusForTrachea, int);
 212     itkGetConstMacro(RadiusForTrachea, int);
 213     GGO_DefineOption(radius, SetRadiusForTrachea, int);
 214
 215     void SetLabelizeParameters3(LabelParamType * a) { m_LabelizeParameters3 = a; }
 216     itkGetConstMacro(LabelizeParameters3, LabelParamType*);
 217     GGO_DefineOption_LabelParam(3, SetLabelizeParameters3, LabelParamType);
 218
 219     // Step 5 options LungSeparation
 220     //     itkSetMacro(FinalOpenClose, bool);
 221     //     itkGetConstMacro(FinalOpenClose, bool);
 222     //     itkBooleanMacro(FinalOpenClose);
 223
 224     // Bronchial bifurcations
 225     itkSetMacro(FindBronchialBifurcations, bool);
 226     itkGetConstMacro(FindBronchialBifurcations, bool);
 227     itkBooleanMacro(FindBronchialBifurcations);
 228
 229   protected:
 230     ExtractLungFilter();
 231     virtual ~ExtractLungFilter() {}
 232
 233     // Main members
 234     InputImageConstPointer input;
 235     MaskImageConstPointer patient;
 236     InputImagePointer working_input;
 237     typename InternalImageType::Pointer working_image;
 238     typename InternalImageType::Pointer trachea_tmp;
 239     MaskImagePointer trachea;
 240     MaskImagePointer output;
 241     unsigned int m_MaxSeedNumber;
 242
 243     // Global options
 244     itkSetMacro(BackgroundValue, MaskImagePixelType);
 245     itkSetMacro(ForegroundValue, MaskImagePixelType);
 246     MaskImagePixelType m_PatientMaskBackgroundValue;
 247     MaskImagePixelType m_BackgroundValue;
 248     MaskImagePixelType m_ForegroundValue;
 249     int m_MinimalComponentSize;
 250
 251     // Step 1
 252     InputImagePixelType m_UpperThreshold;
 253     InputImagePixelType m_LowerThreshold;
 254     bool m_UseLowerThreshold;
 255     LabelParamType* m_LabelizeParameters1;
 256
 257     // Step 2
 258     InputImagePixelType m_UpperThresholdForTrachea;
 259     InputImagePixelType m_ThresholdStepSizeForTrachea;
 260     double m_MultiplierForTrachea;
 261     std::vector<InternalIndexType> m_Seeds;
 262     int m_NumberOfSlicesToSkipBeforeSearchingSeed;
 263
 264     // Step 3
 265     int m_NumberOfHistogramBins;
 266     LabelParamType* m_LabelizeParameters2;
 267
 268     // Step 4
 269     int m_RadiusForTrachea;
 270     LabelParamType* m_LabelizeParameters3;
 271
 272     // Step 5
 273     //     bool m_FinalOpenClose;
 274     bool m_FindBronchialBifurcations;
 275
 276     virtual void GenerateOutputInformation();
 277     virtual void GenerateData();
 278
 279     TreeType m_SkeletonTree;
 280
 281     void TrackFromThisIndex(std::vector<BifurcationType> & listOfBifurcations,
 282                             MaskImagePointer skeleton,
 283                             MaskImageIndexType index,
 284                             MaskImagePixelType label,
 285                             TreeIterator currentNode);
 286
 287
 288     bool SearchForTracheaSeed(int skip);
 289     void SearchForTrachea();
 290     void TracheaRegionGrowing();
 291     double ComputeTracheaVolume();
 292
 293   private:
 294     ExtractLungFilter(const Self&); //purposely not implemented
 295     void operator=(const Self&); //purposely not implemented
 296
 297   }; // end class
 298   //--------------------------------------------------------------------
 299
 300 } // end namespace clitk
 301 //--------------------------------------------------------------------
 302
 303 #ifndef ITK_MANUAL_INSTANTIATION
 304 #include "clitkExtractLungFilter.txx"
 305 #endif
 306
 307 #endif