]> Creatis software - clitk.git/commitdiff
Extract vessels v0.1
authordsarrut <david.sarrut@gmail.com>
Mon, 2 May 2011 09:16:06 +0000 (11:16 +0200)
committerdsarrut <david.sarrut@gmail.com>
Mon, 2 May 2011 09:16:06 +0000 (11:16 +0200)
segmentation/clitkExtractMediastinalVessels.cxx [new file with mode: 0644]
segmentation/clitkExtractMediastinalVessels.ggo [new file with mode: 0644]
segmentation/clitkExtractMediastinalVesselsFilter.h [new file with mode: 0644]
segmentation/clitkExtractMediastinalVesselsFilter.txx [new file with mode: 0644]
segmentation/clitkExtractMediastinalVesselsGenericFilter.h [new file with mode: 0644]
segmentation/clitkExtractMediastinalVesselsGenericFilter.txx [new file with mode: 0644]

diff --git a/segmentation/clitkExtractMediastinalVessels.cxx b/segmentation/clitkExtractMediastinalVessels.cxx
new file mode 100644 (file)
index 0000000..4171f05
--- /dev/null
@@ -0,0 +1,45 @@
+/*=========================================================================
+  Program:   vv                     http://www.creatis.insa-lyon.fr/rio/vv
+
+  Authors belong to:
+  - University of LYON              http://www.universite-lyon.fr/
+  - Léon Bérard cancer center       http://oncora1.lyon.fnclcc.fr
+  - CREATIS CNRS laboratory         http://www.creatis.insa-lyon.fr
+
+  This software is distributed WITHOUT ANY WARRANTY; without even
+  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+  PURPOSE.  See the copyright notices for more information.
+
+  It is distributed under dual licence
+
+  - BSD        See included LICENSE.txt file
+  - CeCILL-B   http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+======================================================================-====*/
+
+// clitk
+#include "clitkExtractMediastinalVessels_ggo.h"
+#include "clitkExtractMediastinalVesselsGenericFilter.h"
+
+//--------------------------------------------------------------------
+int main(int argc, char * argv[])
+{
+
+  // Init command line
+  GGO(clitkExtractMediastinalVessels, args_info);
+  CLITK_INIT;
+
+  // Filter
+  typedef clitk::ExtractMediastinalVesselsGenericFilter<args_info_clitkExtractMediastinalVessels> FilterType;
+  FilterType::Pointer filter = FilterType::New();
+
+  filter->SetArgsInfo(args_info);
+  
+  try {
+    filter->Update();
+  } catch(std::runtime_error e) {
+    std::cout << e.what() << std::endl;
+  }
+
+  return EXIT_SUCCESS;
+} // This is the end, my friend
+//--------------------------------------------------------------------
diff --git a/segmentation/clitkExtractMediastinalVessels.ggo b/segmentation/clitkExtractMediastinalVessels.ggo
new file mode 100644 (file)
index 0000000..718c24e
--- /dev/null
@@ -0,0 +1,21 @@
+#File clitkExtractMediastinalVessels.ggo
+package "clitkExtractMediastinalVessels"
+version "1.0"
+purpose "Extract MediastinalVessels"
+
+option "config"                -  "Config file"                  string        no
+option "imagetypes"     -  "Display allowed image types"  flag          off
+option "verbose"        v  "Verbose"                     flag          off
+option "verboseStep"    -  "Verbose each step"           flag          off
+option "writeStep"      w  "Write image at each step"    flag          off
+option "verboseOption"  -  "Display options values"       flag          off
+option "verboseWarningOff" -  "Do not display warning"    flag          off
+option "verboseMemory"  -  "Display memory usage"         flag          off
+
+section "I/O"
+
+option "afdb"          a       "Input Anatomical Feature DB"     string        yes
+option "input"         i       "Input CT filename"               string        yes
+option "output"        o       "Output folder"                   string        yes
+
+option "threshold"     t       "Initial threshold"               double default="140" no
diff --git a/segmentation/clitkExtractMediastinalVesselsFilter.h b/segmentation/clitkExtractMediastinalVesselsFilter.h
new file mode 100644 (file)
index 0000000..665f79c
--- /dev/null
@@ -0,0 +1,139 @@
+/*=========================================================================
+  Program:   vv                     http://www.creatis.insa-lyon.fr/rio/vv
+
+  Authors belong to: 
+  - University of LYON              http://www.universite-lyon.fr/
+  - Léon Bérard cancer center       http://oncora1.lyon.fnclcc.fr
+  - CREATIS CNRS laboratory         http://www.creatis.insa-lyon.fr
+
+  This software is distributed WITHOUT ANY WARRANTY; without even
+  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+  PURPOSE.  See the copyright notices for more information.
+
+  It is distributed under dual licence
+
+  - BSD        See included LICENSE.txt file
+  - CeCILL-B   http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+  ======================================================================-====*/
+
+#ifndef CLITKEXTRACTLYMPHSTATIONSFILTER_H
+#define CLITKEXTRACTLYMPHSTATIONSFILTER_H
+
+// clitk
+#include "clitkFilterBase.h"
+#include "clitkFilterWithAnatomicalFeatureDatabaseManagement.h"
+
+namespace clitk {
+  
+  //--------------------------------------------------------------------
+  /*
+    Try to extract the some Mediastinal Vessels in a thorax CT.
+    Need a set of Anatomical Features (AFDB)
+  */
+  //--------------------------------------------------------------------
+  
+  template <class TImageType>
+  class ITK_EXPORT ExtractMediastinalVesselsFilter: 
+    public virtual clitk::FilterBase, 
+    public clitk::FilterWithAnatomicalFeatureDatabaseManagement,
+    public itk::ImageToImageFilter<TImageType, itk::Image<uchar, 3> >
+  {
+
+  public:
+    /** Standard class typedefs. */
+    typedef itk::ImageToImageFilter<TImageType, itk::Image<uchar, 3> > Superclass;
+    typedef ExtractMediastinalVesselsFilter          Self;
+    typedef itk::SmartPointer<Self>             Pointer;
+    typedef itk::SmartPointer<const Self>       ConstPointer;
+    
+    /** Method for creation through the object factory. */
+    itkNewMacro(Self);
+    
+    /** Run-time type information (and related methods). */
+    itkTypeMacro(ExtractMediastinalVesselsFilter, ImageToImageFilter);
+
+    /** Some convenient typedefs. */
+    typedef TImageType                       ImageType;
+    typedef typename ImageType::ConstPointer ImageConstPointer;
+    typedef typename ImageType::Pointer      ImagePointer;
+    typedef typename ImageType::RegionType   ImageRegionType; 
+    typedef typename ImageType::PixelType    ImagePixelType; 
+    typedef typename ImageType::SizeType     ImageSizeType; 
+    typedef typename ImageType::IndexType    ImageIndexType; 
+    typedef typename ImageType::PointType    ImagePointType; 
+        
+    typedef uchar MaskImagePixelType;
+    typedef itk::Image<MaskImagePixelType, 3>    MaskImageType;  
+    typedef typename MaskImageType::Pointer      MaskImagePointer;
+    typedef typename MaskImageType::RegionType   MaskImageRegionType; 
+    typedef typename MaskImageType::SizeType     MaskImageSizeType; 
+    typedef typename MaskImageType::IndexType    MaskImageIndexType; 
+    typedef typename MaskImageType::PointType    MaskImagePointType; 
+
+    typedef itk::Image<MaskImagePixelType, 2>    MaskSliceType;
+    typedef typename MaskSliceType::Pointer      MaskSlicePointer;
+    typedef typename MaskSliceType::PointType    MaskSlicePointType;
+
+    typedef long LabelType;
+    /** ImageDimension constants */
+    itkStaticConstMacro(ImageDimension, unsigned int, ImageType::ImageDimension);
+    FILTERBASE_INIT;
+   
+    itkGetConstMacro(BackgroundValue, MaskImagePixelType);
+    itkSetMacro(BackgroundValue, MaskImagePixelType);
+
+    itkGetConstMacro(ForegroundValue, MaskImagePixelType);
+    itkSetMacro(ForegroundValue, MaskImagePixelType);
+
+    itkGetConstMacro(TemporaryForegroundValue, MaskImagePixelType);
+    itkSetMacro(TemporaryForegroundValue, MaskImagePixelType);
+
+    itkGetConstMacro(OutputFolder, std::string);
+    itkSetMacro(OutputFolder, std::string);
+
+    itkGetConstMacro(Threshold, ImagePixelType);
+    itkSetMacro(Threshold, ImagePixelType);
+
+  protected:
+    ExtractMediastinalVesselsFilter();
+    virtual ~ExtractMediastinalVesselsFilter() {}
+    
+    virtual void GenerateOutputInformation();
+    virtual void GenerateInputRequestedRegion();
+    virtual void GenerateData();
+    
+    std::string        m_OutputFolder;
+    ImagePointer       m_Input;
+    MaskImagePointer   m_Working_Support;
+    MaskImagePointer   m_Mediastinum;
+    MaskImagePointer   m_Mask;
+    MaskImagePixelType m_BackgroundValue;
+    MaskImagePixelType m_ForegroundValue;
+    MaskImagePixelType m_TemporaryForegroundValue;
+    ImagePixelType     m_Threshold;
+
+    std::vector<MaskSlicePointer> m_slice_recon;
+    
+    void CropSupInf();
+    //void SearchBrachioCephalicArtery(int & BCA_first_slice, LabelType & BCA_first_label);
+    void TrackBifurcationFromPoint(MaskImagePointer & recon, 
+                                   std::vector<MaskSlicePointer> & slices_recon, 
+                                   MaskImagePointType BCA_p, 
+                                   LabelType newLabel, 
+                                   std::vector<MaskImagePointType> & bif);
+
+  private:
+    ExtractMediastinalVesselsFilter(const Self&); //purposely not implemented
+    void operator=(const Self&); //purposely not implemented
+    
+  }; // end class
+  //--------------------------------------------------------------------
+
+} // end namespace clitk
+//--------------------------------------------------------------------
+
+#ifndef ITK_MANUAL_INSTANTIATION
+#include "clitkExtractMediastinalVesselsFilter.txx"
+#endif
+
+#endif
diff --git a/segmentation/clitkExtractMediastinalVesselsFilter.txx b/segmentation/clitkExtractMediastinalVesselsFilter.txx
new file mode 100644 (file)
index 0000000..02576f4
--- /dev/null
@@ -0,0 +1,402 @@
+/*=========================================================================
+  Program:   vv                     http://www.creatis.insa-lyon.fr/rio/vv
+
+  Authors belong to: 
+  - University of LYON              http://www.universite-lyon.fr/
+  - Léon Bérard cancer center       http://oncora1.lyon.fnclcc.fr
+  - CREATIS CNRS laboratory         http://www.creatis.insa-lyon.fr
+
+  This software is distributed WITHOUT ANY WARRANTY; without even
+  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+  PURPOSE.  See the copyright notices for more information.
+
+  It is distributed under dual licence
+
+  - BSD        See included LICENSE.txt file
+  - CeCILL-B   http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+  ======================================================================-====*/
+
+#ifndef CLITKEXTRACTLYMPHSTATIONSFILTER_TXX
+#define CLITKEXTRACTLYMPHSTATIONSFILTER_TXX
+
+// clitk
+#include "clitkCommon.h"
+#include "clitkExtractMediastinalVesselsFilter.h"
+#include "clitkAutoCropFilter.h"
+#include "clitkSegmentationUtils.h"
+#include "clitkMorphoMathFilter.h"
+
+// itk
+#include <itkBinaryThresholdImageFilter.h>
+#include <itkGrayscaleDilateImageFilter.h>
+#include <itkMinimumMaximumImageCalculator.h>
+
+//--------------------------------------------------------------------
+template <class TImageType>
+clitk::ExtractMediastinalVesselsFilter<TImageType>::
+ExtractMediastinalVesselsFilter():
+  clitk::FilterBase(),
+  clitk::FilterWithAnatomicalFeatureDatabaseManagement(),
+  itk::ImageToImageFilter<TImageType, MaskImageType>()
+{
+  this->SetNumberOfRequiredInputs(1);
+  SetBackgroundValue(0);
+  SetForegroundValue(1);
+  SetThreshold(140);
+  SetTemporaryForegroundValue(1);
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template <class TImageType>
+void 
+clitk::ExtractMediastinalVesselsFilter<TImageType>::
+GenerateOutputInformation() { 
+  // Get inputs
+  LoadAFDB();
+  m_Input = dynamic_cast<ImageType*>(itk::ProcessObject::GetInput(0));
+  
+  // ------------------------------------------------------------------
+  // Sup-Inf crop -> Carina
+  CropSupInf();  
+
+  // ------------------------------------------------------------------
+  // Binarize
+  StartNewStep("Binarize with treshold = "+toString(GetThreshold()));
+  typedef itk::BinaryThresholdImageFilter<ImageType, MaskImageType> BinarizeFilterType; 
+  typename BinarizeFilterType::Pointer binarizeFilter = BinarizeFilterType::New();
+  binarizeFilter->SetInput(m_Input);
+  binarizeFilter->SetLowerThreshold(GetThreshold());
+  binarizeFilter->SetInsideValue(GetTemporaryForegroundValue());
+  binarizeFilter->SetOutsideValue(GetBackgroundValue());
+  binarizeFilter->Update();
+  m_Mask = binarizeFilter->GetOutput();
+  clitk::writeImage<MaskImageType>(m_Mask, "m.mhd");
+  StopCurrentStep<MaskImageType>(m_Mask);
+
+  // Keep main CCL ? 
+  m_Mask = clitk::Labelize<MaskImageType>(m_Mask, GetBackgroundValue(), false, 10);
+  m_Mask = KeepLabels<MaskImageType>(m_Mask, GetBackgroundValue(), GetTemporaryForegroundValue(), 1, 1, true);
+  clitk::writeImage<MaskImageType>(m_Mask, "m2.mhd");
+  
+  // ------------------------------------------------------------------
+  // Extract slices
+  StartNewStep("Detect vessels (slice by slice reconstruction)");
+  std::vector<MaskSlicePointer> slices_mask;
+  clitk::ExtractSlices<MaskImageType>(m_Mask, 2, slices_mask);
+  
+  DD(slices_mask.size());
+  
+  std::vector<MaskSlicePointer> debug_eroded;
+  std::vector<MaskSlicePointer> debug_labeled;
+  std::vector<MaskSlicePointer> debug_slabeled;
+  
+  int radius = 3;
+  DD(radius); // TO PUT IN OPTION
+
+  // ------------------------------------------------------------------
+  // Loop Slice by Slice -> erode find CCL and reconstruct
+  clitk::MorphoMathFilter<MaskSliceType>::Pointer f= clitk::MorphoMathFilter<MaskSliceType>::New();
+  for(uint i=0; i<slices_mask.size(); i++) {
+    // Erosion
+    f->SetInput(slices_mask[i]);
+    f->SetBackgroundValue(GetBackgroundValue());
+    f->SetForegroundValue(GetTemporaryForegroundValue());
+    f->SetRadius(radius);
+    f->SetOperationType(0); // Erode
+    f->VerboseFlagOff();
+    f->Update();
+    MaskSlicePointer eroded = f->GetOutput();
+    //    writeImage<MaskSliceType>(eroded, "er-"+toString(i)+".mhd");
+    debug_eroded.push_back(eroded);
+      
+    // CCL
+    int nb;
+    MaskSlicePointer labeled = 
+      clitk::LabelizeAndCountNumberOfObjects<MaskSliceType>(eroded, GetBackgroundValue(), false, 1, nb);
+
+    // Relabel, large CCL with large label number
+    for(int n=nb; n>0; n--) {
+      //        DD(n);
+      int li = n;
+      int lo = 2*(nb+1)-li;
+      labeled = clitk::SetBackground<MaskSliceType, MaskSliceType>(labeled, labeled, li, lo, true);
+    }
+    debug_labeled.push_back(labeled);
+
+    // Create kernel for GrayscaleDilateImageFilter
+    typedef itk::BinaryBallStructuringElement<MaskSliceType::PixelType,MaskSliceType::ImageDimension > KernelType;
+    KernelType k;
+    k.SetRadius(radius+1);
+    k.CreateStructuringElement();
+    
+    // Keep the MAX -> we prefer the opposite su change the label
+    typedef itk::GrayscaleDilateImageFilter<MaskSliceType, MaskSliceType, KernelType> FilterType;
+    FilterType::Pointer m = FilterType::New();
+    m->SetKernel(k);
+    m->SetInput(labeled);
+    // DD(m->GetAlgorithm());
+    // m->SetAlgorithm(3);
+    m->Update();
+    MaskSlicePointer s = m->GetOutput();
+
+
+    // Remove Initial BG
+    s = clitk::SetBackground<MaskSliceType, MaskSliceType>(s, slices_mask[i], 
+                                                   GetBackgroundValue(), GetBackgroundValue(), true);
+
+    m_slice_recon.push_back(s);
+  } // end loop
+  DD("end loop");
+  
+  MaskImageType::Pointer eroded = clitk::JoinSlices<MaskImageType>(debug_eroded, m_Mask, 2);
+  clitk::writeImage<MaskImageType>(eroded, "eroded.mhd");
+
+  DD("l");
+  MaskImageType::Pointer l = clitk::JoinSlices<MaskImageType>(debug_labeled, m_Mask, 2);
+  clitk::writeImage<MaskImageType>(l, "labeled.mhd");
+  
+  DD("r");
+  MaskImageType::Pointer r = clitk::JoinSlices<MaskImageType>(m_slice_recon, m_Mask, 2);
+  clitk::writeImage<MaskImageType>(r, "recon.mhd");
+  
+  // ------------------------------------------------------------------
+  // Loop Slice by Slice -> BCA not found yet
+  /*  MaskImagePointType BCA_p;
+  GetAFDB()->GetPoint3D("BrachioCephalicArteryFirstInferiorPoint", BCA_p);
+  DD(BCA_p);
+  MaskImagePointType bif1;
+  MaskImagePointType bif2;
+  TrackBifurcationFromPoint(r, BCA_p, bif1, bif2);
+  DD(bif1);
+  DD(bif2);
+  */
+  // Find max label
+  typedef itk::MinimumMaximumImageCalculator<MaskImageType> MinMaxFilterType;
+  MinMaxFilterType::Pointer ff = MinMaxFilterType::New();
+  ff->SetImage(r);
+  ff->ComputeMaximum();
+  LabelType newLabel = ff->GetMaximum()+1; 
+  DD(newLabel);
+
+  // Get all centroids of the first slice
+  std::vector<MaskSlicePointType> centroids2D;
+  clitk::ComputeCentroids<MaskSliceType>(m_slice_recon[0], GetBackgroundValue(), centroids2D);
+  DD(centroids2D.size());
+  std::vector<MaskImagePointType> bifurcations;
+  clitk::PointsUtils<MaskImageType>::Convert2DListTo3DList(centroids2D, 0, r, bifurcations);  
+  DD(bifurcations.size());
+  for(uint i=1; i<bifurcations.size()+1; i++) {
+    DD(i);
+    DD(bifurcations.size());
+    TrackBifurcationFromPoint(r, m_slice_recon, bifurcations[i], newLabel+i, bifurcations);
+    DD("end track");
+    DD(bifurcations.size());
+    MaskImageType::Pointer rr = clitk::JoinSlices<MaskImageType>(m_slice_recon, m_Mask, 2);
+    clitk::writeImage<MaskImageType>(rr, "recon"+toString(i)+".mhd");
+  }
+
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template <class TImageType>
+void 
+clitk::ExtractMediastinalVesselsFilter<TImageType>::
+GenerateInputRequestedRegion() {
+  //DD("GenerateInputRequestedRegion (nothing?)");
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template <class TImageType>
+void 
+clitk::ExtractMediastinalVesselsFilter<TImageType>::
+GenerateData() {
+  DD("GenerateData");
+  // Final Step -> graft output (if SetNthOutput => redo)
+  MaskImagePointer BrachioCephalicArtery = 
+    GetAFDB()->template GetImage<MaskImageType>("BrachioCephalicArtery");
+  this->GraftNthOutput(0, BrachioCephalicArtery);
+}
+//--------------------------------------------------------------------
+  
+
+//--------------------------------------------------------------------
+template <class TImageType>
+void 
+clitk::ExtractMediastinalVesselsFilter<TImageType>::
+CropSupInf() { 
+  StartNewStep("Inf/Sup limits (carina) and crop with mediastinum");
+  // Get Trachea and Carina
+  MaskImagePointer Trachea = GetAFDB()->template GetImage <MaskImageType>("Trachea");  
+  
+  // Get or compute Carina
+  double m_CarinaZ;
+  // Get Carina Z position
+  MaskImagePointer Carina = GetAFDB()->template GetImage<MaskImageType>("Carina");
+  
+  std::vector<MaskImagePointType> centroids;
+  clitk::ComputeCentroids<MaskImageType>(Carina, GetBackgroundValue(), centroids);
+  m_CarinaZ = centroids[1][2];
+  // DD(m_CarinaZ);
+  // add one slice to include carina ?
+  m_CarinaZ += Carina->GetSpacing()[2];
+  // We dont need Carina structure from now
+  Carina->Delete();
+  GetAFDB()->SetDouble("CarinaZ", m_CarinaZ);
+  
+  // Crop
+  m_Input = clitk::CropImageRemoveLowerThan<ImageType>(m_Input, 2, 
+                                                       m_CarinaZ, false, GetBackgroundValue());  
+
+  // Crop not post to centroid
+  double m_CarinaY = centroids[1][1];
+  DD(m_CarinaY);
+  m_Input = clitk::CropImageRemoveGreaterThan<ImageType>(m_Input, 1, // OLD ABOVE
+                                                         m_CarinaY, false, GetBackgroundValue());  
+  // Crop not ant to centroid
+  m_Input = clitk::CropImageRemoveLowerThan<ImageType>(m_Input, 1, 
+                                                         m_CarinaY-80, false, GetBackgroundValue());  
+
+  // AutoCrop with Mediastinum
+  m_Mediastinum  = GetAFDB()->template GetImage<MaskImageType>("Mediastinum");
+  // Resize like input (sup to carina)
+  m_Mediastinum = clitk::ResizeImageLike<MaskImageType>(m_Mediastinum, m_Input, GetBackgroundValue());
+  // Auto crop
+  m_Mediastinum = clitk::AutoCrop<MaskImageType>(m_Mediastinum, GetBackgroundValue());
+  // Resize input
+  m_Input = clitk::ResizeImageLike<ImageType>(m_Input, m_Mediastinum, GetBackgroundValue());
+
+  // End
+  StopCurrentStep<ImageType>(m_Input);
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template <class TImageType>
+void 
+clitk::ExtractMediastinalVesselsFilter<TImageType>::
+//SearchBrachioCephalicArtery(int & BCA_first_slice, LabelType & BCA_first_label) { 
+TrackBifurcationFromPoint(MaskImagePointer & recon, 
+                          std::vector<MaskSlicePointer> & slices_recon, 
+                          MaskImagePointType point3D, 
+                          LabelType newLabel,
+                          std::vector<MaskImagePointType> & bifurcations) {
+  StartNewStep("Search for BCA first slice and label");
+  DD(newLabel);
+
+  // Extract slices
+  //  std::vector<MaskSlicePointer> slices_recon;
+  //clitk::ExtractSlices<MaskImageType>(recon, 2, slices_recon);
+
+  // Find first slice
+  MaskImageIndexType index;
+  recon->TransformPhysicalPointToIndex(point3D, index);
+  DD(point3D);
+  DD(index);
+
+  uint i=index[2]; 
+  bool found = false;
+  LabelType previous_largest_label=recon->GetPixel(index);
+  DD(previous_largest_label);
+  do {
+    DD(i);
+    // Consider current reconstructed slice
+    MaskSlicePointer s = slices_recon[i];
+    MaskSlicePointer previous;
+    if (i==index[2]) previous = s;
+    else previous = slices_recon[i-1];
+    
+    // Get centroids of the labels in the current slice
+    typedef typename MaskSliceType::PointType SlicePointType;
+    static const unsigned int Dim = MaskSliceType::ImageDimension;
+    typedef itk::ShapeLabelObject< LabelType, Dim > LabelObjectType;
+    typedef itk::LabelMap< LabelObjectType > LabelMapType;
+    typedef itk::LabelImageToLabelMapFilter<MaskSliceType, LabelMapType> ImageToMapFilterType;
+    typename ImageToMapFilterType::Pointer imageToLabelFilter = ImageToMapFilterType::New(); 
+    typedef itk::ShapeLabelMapFilter<LabelMapType, MaskSliceType> ShapeFilterType; 
+    typename ShapeFilterType::Pointer statFilter = ShapeFilterType::New();
+    imageToLabelFilter->SetBackgroundValue(GetBackgroundValue());
+    imageToLabelFilter->SetInput(s);
+    statFilter->SetInput(imageToLabelFilter->GetOutput());
+    statFilter->Update();
+    typename LabelMapType::Pointer labelMap = statFilter->GetOutput();
+
+    // Look what centroid inside the previous largest one
+    std::vector<SlicePointType> centroids;
+    std::vector<LabelType> centroids_label;
+    for(uint c=0; c<labelMap->GetNumberOfLabelObjects(); c++) {
+      int label = labelMap->GetLabels()[c];
+      DD(label);
+      SlicePointType center = labelMap->GetLabelObject(label)->GetCentroid();
+      SlicePointType center_previous = center;
+      center_previous[2] -= m_Input->GetSpacing()[2];      
+      // Get label into previous slice
+      typename MaskSliceType::IndexType index;
+      previous->TransformPhysicalPointToIndex(center_previous, index);
+      LabelType l = previous->GetPixel(index);
+      DD(l);
+      if (l == previous_largest_label) {
+        centroids.push_back(center);
+        centroids_label.push_back(label);
+      }
+    }
+    DD(centroids.size());
+    
+    // If several centroids, we found a bifurcation
+    if (centroids.size() > 1) {
+      found = true;
+      for(uint c=0; c<centroids.size(); c++) {
+        ImagePointType bif;
+        clitk::PointsUtils<MaskImageType>::Convert2DTo3D(centroids[c], m_Mask, i, bif);
+        bifurcations.push_back(bif);
+        s = clitk::SetBackground<MaskSliceType, MaskSliceType>(s, s, centroids_label[c], newLabel+c+1, true);
+        //        slices_recon[i] = s; // (useful ?)
+      }
+      DD("FOUND");
+    }
+    
+    // if only one centroids, we change the current image with the current label 
+    if (centroids.size() == 1) {
+      s = clitk::SetBackground<MaskSliceType, MaskSliceType>(s, s, centroids_label[0], newLabel, true);
+      previous_largest_label = newLabel;
+      /*typedef itk::BinaryThresholdImageFilter<MaskSliceType, MaskSliceType> BinarizeFilterType; 
+      typename BinarizeFilterType::Pointer binarizeFilter = BinarizeFilterType::New();
+      binarizeFilter->SetInput(s);
+      binarizeFilter->SetLowerThreshold(centroids_label[0]);
+      binarizeFilter->SetUpperThreshold(centroids_label[0]+1);
+      binarizeFilter->SetInsideValue(previous_largest_label);
+      binarizeFilter->SetOutsideValue(GetBackgroundValue());
+      binarizeFilter->Update();
+      s = binarizeFilter->GetOutput();*/
+      slices_recon[i] = s; // (not useful ?)
+    }
+
+    if (centroids.size() == 0) {
+      DD("no centroid, I stop");
+      found = true;
+    }
+
+    if (i == slices_recon.size()-1) found = true;
+
+    // iterate
+    ++i;
+  } while (!found);
+
+
+  //MaskImageType::Pointer rr = clitk::JoinSlices<MaskImageType>(slices_recon, m_Mask, 2);
+  //clitk::writeImage<MaskImageType>(rr, "recon2.mhd");
+
+  // End
+  StopCurrentStep();
+}
+//--------------------------------------------------------------------
+
+
+
+#endif //#define CLITKEXTRACTMEDIASTINALVESSELSFILTER_TXX
diff --git a/segmentation/clitkExtractMediastinalVesselsGenericFilter.h b/segmentation/clitkExtractMediastinalVesselsGenericFilter.h
new file mode 100644 (file)
index 0000000..9e53d7a
--- /dev/null
@@ -0,0 +1,80 @@
+/*=========================================================================
+  Program:   vv                     http://www.creatis.insa-lyon.fr/rio/vv
+
+  Authors belong to: 
+  - University of LYON              http://www.universite-lyon.fr/
+  - Léon Bérard cancer center       http://oncora1.lyon.fnclcc.fr
+  - CREATIS CNRS laboratory         http://www.creatis.insa-lyon.fr
+
+  This software is distributed WITHOUT ANY WARRANTY; without even
+  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+  PURPOSE.  See the copyright notices for more information.
+
+  It is distributed under dual licence
+
+  - BSD        See included LICENSE.txt file
+  - CeCILL-B   http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+======================================================================-====*/
+
+#ifndef CLITKEXTRACTLYMPHSTATIONSSGENERICFILTER_H
+#define CLITKEXTRACTLYMPHSTATIONSSGENERICFILTER_H
+
+#include "clitkIO.h"
+#include "clitkImageToImageGenericFilter.h"
+#include "clitkExtractMediastinalVesselsFilter.h"
+
+//--------------------------------------------------------------------
+namespace clitk 
+{
+  
+  template<class ArgsInfoType>
+  class ITK_EXPORT ExtractMediastinalVesselsGenericFilter: 
+    public ImageToImageGenericFilter<ExtractMediastinalVesselsGenericFilter<ArgsInfoType> >
+  {
+    
+  public:
+    //--------------------------------------------------------------------
+    ExtractMediastinalVesselsGenericFilter();
+
+    //--------------------------------------------------------------------
+    typedef ImageToImageGenericFilter<ExtractMediastinalVesselsGenericFilter<ArgsInfoType> > Superclass;
+    typedef ExtractMediastinalVesselsGenericFilter Self;
+    typedef itk::SmartPointer<Self>           Pointer;
+    typedef itk::SmartPointer<const Self>     ConstPointer;
+
+    //--------------------------------------------------------------------
+    itkNewMacro(Self);  
+    itkTypeMacro(ExtractMediastinalVesselsGenericFilter, LightObject);
+
+    //--------------------------------------------------------------------
+    // Options for the GenericFilter
+    void SetArgsInfo(const ArgsInfoType & a);
+
+    //--------------------------------------------------------------------
+    // Options for the Filter
+    template<class FilterType> 
+    void SetOptionsFromArgsInfoToFilter(FilterType * f) ;
+
+    //--------------------------------------------------------------------
+    // Main function called each time the filter is updated
+    template<class ImageType>  
+    void UpdateWithInputImageType();
+
+  protected:
+    template<unsigned int Dim> void InitializeImageType();
+    ArgsInfoType mArgsInfo;
+    
+  private:
+    ExtractMediastinalVesselsGenericFilter(const Self&); //purposely not implemented
+    void operator=(const Self&); //purposely not implemented
+    
+  }; // end class
+  //--------------------------------------------------------------------
+    
+} // end namespace clitk
+
+#ifndef ITK_MANUAL_INSTANTIATION
+#include "clitkExtractMediastinalVesselsGenericFilter.txx"
+#endif
+
+#endif // #define CLITKEXTRACTLYMPHSTATIONSSGENERICFILTER_H
diff --git a/segmentation/clitkExtractMediastinalVesselsGenericFilter.txx b/segmentation/clitkExtractMediastinalVesselsGenericFilter.txx
new file mode 100644 (file)
index 0000000..b85da27
--- /dev/null
@@ -0,0 +1,109 @@
+/*=========================================================================
+  Program:   vv                     http://www.creatis.insa-lyon.fr/rio/vv
+
+  Authors belong to: 
+  - University of LYON              http://www.universite-lyon.fr/
+  - Léon Bérard cancer center       http://oncora1.lyon.fnclcc.fr
+  - CREATIS CNRS laboratory         http://www.creatis.insa-lyon.fr
+
+  This software is distributed WITHOUT ANY WARRANTY; without even
+  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+  PURPOSE.  See the copyright notices for more information.
+
+  It is distributed under dual licence
+
+  - BSD        See included LICENSE.txt file
+  - CeCILL-B   http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+  ======================================================================-====*/
+
+#ifndef CLITKEXTRACTLYMPHSTATIONSSGENERICFILTER_TXX
+#define CLITKEXTRACTLYMPHSTATIONSSGENERICFILTER_TXX
+
+#include "clitkImageCommon.h"
+  
+//--------------------------------------------------------------------
+template<class ArgsInfoType>
+clitk::ExtractMediastinalVesselsGenericFilter<ArgsInfoType>::ExtractMediastinalVesselsGenericFilter():
+  ImageToImageGenericFilter<Self>("ExtractMediastinalVessels") 
+{
+  // Default values
+  cmdline_parser_clitkExtractMediastinalVessels_init(&mArgsInfo);
+  InitializeImageType<3>(); // Only for 3D images
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template<class ArgsInfoType>
+template<unsigned int Dim>
+void clitk::ExtractMediastinalVesselsGenericFilter<ArgsInfoType>::InitializeImageType() 
+{  
+  ADD_IMAGE_TYPE(Dim, short); // Can add float later
+}
+//--------------------------------------------------------------------
+  
+
+//--------------------------------------------------------------------
+template<class ArgsInfoType>
+void clitk::ExtractMediastinalVesselsGenericFilter<ArgsInfoType>::SetArgsInfo(const ArgsInfoType & a) 
+{
+  mArgsInfo=a;
+  SetIOVerbose(mArgsInfo.verbose_flag);
+  if (mArgsInfo.imagetypes_flag) this->PrintAvailableImageTypes();
+  if (mArgsInfo.input_given) AddInputFilename(mArgsInfo.input_arg);
+  if (mArgsInfo.output_given) AddOutputFilename(mArgsInfo.output_arg);
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template<class ArgsInfoType>
+template<class FilterType>
+void 
+clitk::ExtractMediastinalVesselsGenericFilter<ArgsInfoType>::
+SetOptionsFromArgsInfoToFilter(FilterType * f)
+{
+  f->SetVerboseOptionFlag(mArgsInfo.verbose_flag);
+  f->SetVerboseStepFlag(mArgsInfo.verboseStep_flag);
+  f->SetWriteStepFlag(mArgsInfo.writeStep_flag);
+  f->SetVerboseMemoryFlag(mArgsInfo.verboseMemory_flag);
+  f->SetAFDBFilename(mArgsInfo.afdb_arg);  
+  f->SetOutputFolder(mArgsInfo.output_arg);
+
+  f->SetThreshold(mArgsInfo.threshold_arg);
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+// Update with the number of dimensions and the pixeltype
+//--------------------------------------------------------------------
+template<class ArgsInfoType>
+template<class ImageType>
+void clitk::ExtractMediastinalVesselsGenericFilter<ArgsInfoType>::UpdateWithInputImageType() 
+{ 
+  // Reading input
+  typename ImageType::Pointer input = this->template GetInput<ImageType>(0);
+
+  // Create filter
+  typedef clitk::ExtractMediastinalVesselsFilter<ImageType> FilterType;
+  typename FilterType::Pointer filter = FilterType::New();
+    
+  // Set global Options 
+  filter->SetInput(input);
+  SetOptionsFromArgsInfoToFilter<FilterType>(filter);
+
+  // Go !
+  filter->Update();
+
+  // Write/Save results
+  for(uint i=0; i<filter->GetNumberOfOutputs(); i++) {
+    typedef uchar MaskImagePixelType;
+    typedef itk::Image<MaskImagePixelType, 3> OutputImageType;
+    typename OutputImageType::Pointer output = filter->GetOutput(i);
+    this->template SetNextOutput<OutputImageType>(output); 
+  }
+}
+//--------------------------------------------------------------------
+
+#endif //#define CLITKEXTRACTLYMPHSTATIONSSGENERICFILTER_TXX