]> Creatis software - clitk.git/blob - segmentation/clitkExtractAirwayTreeInfoFilter.h
now dump anatomical info (carina position) into a DB
[clitk.git] / segmentation / clitkExtractAirwayTreeInfoFilter.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 CLITKEXTRACTAIRWAYTREEINFOSFILTER_H
20 #define CLITKEXTRACTAIRWAYTREEINFOSFILTER_H
21
22 // clitk 
23 #include "clitkFilterBase.h"
24 #include "clitkDecomposeAndReconstructImageFilter.h"
25 #include "clitkExplosionControlledThresholdConnectedImageFilter.h"
26 #include "clitkSegmentationUtils.h"
27
28 // itk
29 #include "itkStatisticsImageFilter.h"
30
31 namespace clitk {
32   
33   //--------------------------------------------------------------------
34   /*
35     From a trachea binary image, compute the skeleton and track the
36     path to find the carena.
37   */
38   //--------------------------------------------------------------------
39   
40
41   //--------------------------------------------------------------------
42   template<class IndexType, class PixelType>
43   class Bifurcation
44   {
45   public:
46     Bifurcation(IndexType _index, PixelType _l, PixelType _l1, PixelType _l2) {
47       index = _index;
48       _l = l;
49       _l1 = l1;
50       _l2 = l2;
51     }
52     IndexType index;
53     PixelType l;
54     PixelType l1;
55     PixelType l2;
56   };
57   //--------------------------------------------------------------------
58
59
60   //--------------------------------------------------------------------
61   template <class TImageType>
62   class ITK_EXPORT ExtractAirwayTreeInfoFilter: 
63     public clitk::FilterBase, 
64     public itk::ImageToImageFilter<TImageType, TImageType> 
65   {
66     
67   public:
68     /** Standard class typedefs. */
69     typedef itk::ImageToImageFilter<TImageType, TImageType> Superclass;
70     typedef ExtractAirwayTreeInfoFilter         Self;
71     typedef itk::SmartPointer<Self>             Pointer;
72     typedef itk::SmartPointer<const Self>       ConstPointer;
73     
74     /** Method for creation through the object factory. */
75     itkNewMacro(Self);  
76     
77     /** Run-time type information (and related methods). */
78     itkTypeMacro(ExtractAirwayTreeInfoFilter, ImageToImageFilter);
79     FILTERBASE_INIT;
80
81     /** Some convenient typedefs */
82     typedef TImageType                       ImageType;
83     typedef typename ImageType::ConstPointer ImageConstPointer;
84     typedef typename ImageType::Pointer      ImagePointer;
85     typedef typename ImageType::RegionType   ImageRegionType; 
86     typedef typename ImageType::PixelType    ImagePixelType; 
87     typedef typename ImageType::SizeType     ImageSizeType; 
88     typedef typename ImageType::IndexType    ImageIndexType; 
89     typedef typename ImageType::PointType    ImagePointType; 
90         
91     itkStaticConstMacro(ImageDimension, unsigned int, ImageType::ImageDimension);
92     typedef int InternalPixelType;
93     typedef itk::Image<InternalPixelType, ImageType::ImageDimension> InternalImageType;
94     typedef typename InternalImageType::Pointer                      InternalImagePointer;
95     typedef typename InternalImageType::IndexType                    InternalIndexType;
96     typedef LabelizeParameters<InternalPixelType>                    LabelParamType;
97     
98     /** Connect inputs */
99     void SetInput(const ImageType * image);
100
101     // Set all options at a time
102     template<class ArgsInfoType>
103       void SetArgsInfo(ArgsInfoType arg);
104
105     // Background / Foreground
106     itkGetConstMacro(BackgroundValue, ImagePixelType);
107     itkGetConstMacro(ForegroundValue, ImagePixelType);
108     
109     // Get results
110     itkGetConstMacro(FirstTracheaPoint, ImagePointType);
111     itkGetConstMacro(CarinaPoint, ImagePointType);
112
113   protected:
114     ExtractAirwayTreeInfoFilter();
115     virtual ~ExtractAirwayTreeInfoFilter() {}
116
117     // Main members
118     ImageConstPointer input;
119     ImagePointer skeleton;
120     ImagePointer working_input;
121
122     // Global options
123     itkSetMacro(BackgroundValue, ImagePixelType);
124     itkSetMacro(ForegroundValue, ImagePixelType);
125     ImagePixelType m_BackgroundValue;
126     ImagePixelType m_ForegroundValue;
127
128     // Results
129     ImagePointType m_FirstTracheaPoint;
130     ImagePointType m_CarinaPoint;
131     
132     virtual void GenerateOutputInformation();
133     virtual void GenerateData();
134
135     typedef Bifurcation<ImageIndexType,ImagePixelType> BifurcationType;
136     void TrackFromThisIndex(std::vector<BifurcationType> & listOfBifurcations, 
137                             ImagePointer skeleton, 
138                             ImageIndexType index,
139                             ImagePixelType label);
140   private:
141     ExtractAirwayTreeInfoFilter(const Self&); //purposely not implemented
142     void operator=(const Self&); //purposely not implemented
143     
144   }; // end class
145   //--------------------------------------------------------------------
146
147 } // end namespace clitk
148 //--------------------------------------------------------------------
149
150 #ifndef ITK_MANUAL_INSTANTIATION
151 #include "clitkExtractAirwayTreeInfoFilter.txx"
152 #endif
153
154 #endif