Trachea finding algorithm

- version used in LOLA11
- preserved original version
- new "type" parameter to choose between algorithms
- algorithm-specific params

diff --git a/segmentation/clitkExtractLung.ggo b/segmentation/clitkExtractLung.ggo
index 98ddbd80d3dc29f1dcb53355d663a6dc7e1a4638..0cae0f89ad0a926ba70e365d0e0e779e56b06435 100644 (file)
--- a/segmentation/clitkExtractLung.ggo
+++ b/segmentation/clitkExtractLung.ggo
@@ -31,13 +31,17 @@ option "minSize"     -      "Minimum component size in voxels"        int           no      default="1
 
 section "Step 2 : find trachea"
 
-option "skipslices"                  -  "Number of slices to skip before searching seed" int no default="0"
-option "upperThresholdForTrachea"    - "Initial upper threshold for trachea"  double   no  default="-900"
-option "multiplierForTrachea"       -  "Multiplier for the region growing"    double   no  default="5"
-option "thresholdStepSizeForTrachea" - "Threshold step size"                  int      no  default="64"
-option "seed"                        - "Index of the trachea seed point (in pixel, not in mm)"      int        no  multiple
-option "doNotCheckTracheaVolume"     -  "If set, do not check the trachea volume" flag off
-option "verboseRG"                   -  "Verbose RegionGrowing"   flag off
+option "type"                        -  "trachea finding algorithm (0: original; 1: LOLA11)" int no default="0"
+option "skipslices"                  -  "0: Number of slices to skip before searching seed" int no default="0"
+option "upperThresholdForTrachea"    - "0: Initial upper threshold for trachea"  double        no  default="-900"
+option "multiplierForTrachea"       -  "0: Multiplier for the region growing"    double        no  default="5"
+option "thresholdStepSizeForTrachea" - "0: Threshold step size"                       int      no  default="64"
+option "seed"                        - "0,1: Index of the trachea seed point (in pixel, not in mm)"      int   no  multiple
+option "doNotCheckTracheaVolume"     -  "0,1: If set, do not check the trachea volume" flag off
+option "verboseRG"                   -  "0,1: Verbose RegionGrowing"   flag off
+option "maxElongation"               -  "1: Maximum allowed elongation of candidate regions for the trachea" float no default="0.5"
+option "numSlices"                   -  "1: Number of slices to search for trachea" int no default="50"
+option "seedPreProcessingThreshold"   -  "1: Threshold for the pre-processing step of the algorithm" int no default="-400"
 
 section "Step 3 : auto extract lung"
 

diff --git a/segmentation/clitkExtractLungFilter.h b/segmentation/clitkExtractLungFilter.h
index f55023f35b68e0382d04574bfae5c1612a41bf70..ead4d988f36597401d3684ada7e69da62ee68851 100644 (file)
--- a/segmentation/clitkExtractLungFilter.h
+++ b/segmentation/clitkExtractLungFilter.h
@@ -142,6 +142,9 @@ namespace clitk {
 
     void SetLabelizeParameters1(LabelParamType * a) { m_LabelizeParameters1 = a; }
     itkGetConstMacro(LabelizeParameters1, LabelParamType*);
+    
+    itkSetMacro(TracheaSeedAlgorithm, int);
+    itkGetConstMacro(TracheaSeedAlgorithm, int);
 
     // Step 2 options FindTrachea
     itkSetMacro(UpperThresholdForTrachea, InputImagePixelType);
@@ -153,6 +156,14 @@ namespace clitk {
     itkSetMacro(ThresholdStepSizeForTrachea, InputImagePixelType);
     itkGetConstMacro(ThresholdStepSizeForTrachea, InputImagePixelType);
 
+    // options FindTrachea2
+    itkSetMacro(NumSlices, int);
+    itkGetConstMacro(NumSlices, int);
+    itkSetMacro(MaxElongation, double);
+    itkGetConstMacro(MaxElongation, double);
+    itkSetMacro(SeedPreProcessingThreshold, int);
+    itkGetConstMacro(SeedPreProcessingThreshold, int);
+
     void AddSeed(InternalIndexType s);
     std::vector<InternalIndexType> & GetSeeds() { return  m_Seeds; }
 
@@ -200,7 +211,7 @@ namespace clitk {
     itkSetMacro(AutoCrop, bool);
     itkGetConstMacro(AutoCrop, bool);
     itkBooleanMacro(AutoCrop);
-
+    
   protected:
     ExtractLungFilter();
     virtual ~ExtractLungFilter() {}
@@ -231,6 +242,7 @@ namespace clitk {
     LabelParamType* m_LabelizeParameters1;
 
     // Step 2
+    int m_TracheaSeedAlgorithm;
     InputImagePixelType m_UpperThresholdForTrachea;
     InputImagePixelType m_ThresholdStepSizeForTrachea;
     double m_MultiplierForTrachea;
@@ -238,6 +250,9 @@ namespace clitk {
     int m_NumberOfSlicesToSkipBeforeSearchingSeed;
     bool m_TracheaVolumeMustBeCheckedFlag;
     bool m_VerboseRegionGrowingFlag;
+    int m_NumSlices;
+    double m_MaxElongation;
+    int m_SeedPreProcessingThreshold;
 
     // Step 3
     int m_NumberOfHistogramBins;

diff --git a/segmentation/clitkExtractLungFilter.txx b/segmentation/clitkExtractLungFilter.txx
index 7dd27c7a3b8e1df027be41065b17e9de4be2e493..b7b28665baff5006bc2de75cfc158fd8be9821fe 100644 (file)
--- a/segmentation/clitkExtractLungFilter.txx
+++ b/segmentation/clitkExtractLungFilter.txx
@@ -79,11 +79,16 @@ ExtractLungFilter():
   SetLabelizeParameters1(p1);
 
   // Step 2 default values
-  SetUpperThresholdForTrachea(-900);
+  SetTracheaSeedAlgorithm(0);
+  SetUpperThresholdForTrachea(-700);
   SetMultiplierForTrachea(5);
   SetThresholdStepSizeForTrachea(64);
   SetNumberOfSlicesToSkipBeforeSearchingSeed(0);
   TracheaVolumeMustBeCheckedFlagOn();
+  SetNumSlices(50);
+  SetMaxElongation(0.5);
+  SetSeedPreProcessingThreshold(-400);
+ 
   
   // Step 3 default values
   SetNumberOfHistogramBins(500);
@@ -591,6 +596,9 @@ clitk::ExtractLungFilter<ImageType>::
 SearchForTracheaSeed2(int numberOfSlices)
 {
   if (m_Seeds.size() == 0) { // try to find seed (if not zero, it is given by user)    
+    if (GetVerboseFlag())
+      std::cout << "SearchForTracheaSeed2(" << numberOfSlices << ", " << GetMaxElongation() << ")" << std::endl;
+    
     typedef unsigned char MaskPixelType;
     typedef itk::Image<MaskPixelType, ImageType::ImageDimension> MaskImageType;
     typedef itk::Image<typename MaskImageType::PixelType, 2> MaskImageType2D;
@@ -609,7 +617,7 @@ SearchForTracheaSeed2(int numberOfSlices)
     // threshold to isolate airawys and lungs
     typename ThresholdFilterType::Pointer threshold = ThresholdFilterType::New();
     threshold->SetLowerThreshold(-2000);
-    threshold->SetUpperThreshold(-400);
+    threshold->SetUpperThreshold(GetSeedPreProcessingThreshold());
     threshold->SetInput(working_input);
     threshold->Update();
     
@@ -698,7 +706,7 @@ SearchForTracheaSeed2(int numberOfSlices)
       typename LabelImageType::LabelObjectContainerType shapes_map = label_map->GetLabelObjectContainer();
       typename LabelImageType::LabelObjectContainerType::const_iterator s;
       typename ShapeLabelType::Pointer shape, max_e_shape;
-      double max_elongation = 0.5;
+      double max_elongation = GetMaxElongation();
       double max_size = size[0]*size[1]/128;
       double max_e = 0;
       int nshapes = 0;
@@ -738,7 +746,7 @@ SearchForTracheaSeed2(int numberOfSlices)
         p2[2] = prev_e_centre[2];
         
         double mag = (p2 - p1).GetNorm();
-        if (GetVerboseStepFlag()) {
+        if (GetVerboseFlag()) {
           cout.precision(3);
           cout << index[2] << ": ";
           cout << "region(" << max_e_centre[0] << ", " << max_e_centre[1] << ", " << max_e_centre[2] << "); ";
@@ -770,7 +778,7 @@ SearchForTracheaSeed2(int numberOfSlices)
       }
     }
     
-    if (GetVerboseStepFlag()) 
+    if (GetVerboseFlag()) 
       std::cout << "seed at: " << trachea_centre << std::endl;
     m_Seeds.push_back(trachea_centre);
   }
@@ -859,12 +867,18 @@ SearchForTrachea()
   // compute trachea volume
   // if volume not plausible  -> skip more slices and restart 
 
+  bool has_seed;
   bool stop = false;
   double volume = 0.0;
   int skip = GetNumberOfSlicesToSkipBeforeSearchingSeed();
   while (!stop) {
     stop = true;
-    if (SearchForTracheaSeed2(50)) {
+    if (GetTracheaSeedAlgorithm() == 0)
+      has_seed = SearchForTracheaSeed(skip);
+    else
+      has_seed = SearchForTracheaSeed2(GetNumSlices());
+    
+    if (has_seed) {
       TracheaRegionGrowing();
       volume = ComputeTracheaVolume()/1000; // assume mm3, so divide by 1000 to get cc
       if (GetWriteStepFlag()) {

diff --git a/segmentation/clitkExtractLungGenericFilter.txx b/segmentation/clitkExtractLungGenericFilter.txx
index b40dc28fc49de81bd0a4784fba74158b7dd333f6..b76b365e12ef34dca1af9e8fc50c4e8090b5e1f9 100644 (file)
--- a/segmentation/clitkExtractLungGenericFilter.txx
+++ b/segmentation/clitkExtractLungGenericFilter.txx
@@ -89,6 +89,10 @@ SetOptionsFromArgsInfoToFilter(FilterType * f)
   f->SetUpperThresholdForTrachea(mArgsInfo.upperThresholdForTrachea_arg);
   f->SetMultiplierForTrachea(mArgsInfo.multiplierForTrachea_arg);
   f->SetThresholdStepSizeForTrachea(mArgsInfo.thresholdStepSizeForTrachea_arg);
+  f->SetTracheaSeedAlgorithm(mArgsInfo.type_arg);
+  f->SetNumSlices(mArgsInfo.numSlices_arg);
+  f->SetMaxElongation(mArgsInfo.maxElongation_arg);
+  f->SetSeedPreProcessingThreshold(mArgsInfo.seedPreProcessingThreshold_arg);
 
   typename FilterType::InputImageIndexType s;
   if (mArgsInfo.seed_given) {