1 /*=========================================================================
2 Program: vv http://www.creatis.insa-lyon.fr/rio/vv
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
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.
13 It is distributed under dual licence
15 - BSD See included LICENSE.txt file
16 - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
17 ===========================================================================**/
19 #ifndef CLITKEXTRACTBONESSFILTER_TXX
20 #define CLITKEXTRACTBONESSFILTER_TXX
23 #include "clitkImageCommon.h"
24 #include "clitkSetBackgroundImageFilter.h"
25 #include "clitkSegmentationUtils.h"
26 #include "clitkAutoCropFilter.h"
27 #include "clitkFillMaskFilter.h"
30 #include "itkBinaryThresholdImageFilter.h"
31 #include "itkConnectedComponentImageFilter.h"
32 #include "itkRelabelComponentImageFilter.h"
33 #include "itkNeighborhoodConnectedImageFilter.h"
34 #include "itkCurvatureAnisotropicDiffusionImageFilter.h"
36 //--------------------------------------------------------------------
37 template <class TInputImageType>
38 clitk::ExtractBonesFilter<TInputImageType>::
41 clitk::FilterWithAnatomicalFeatureDatabaseManagement(),
42 itk::ImageToImageFilter<TInputImageType, MaskImageType>()
44 // Default global options
45 this->SetNumberOfRequiredInputs(1);
46 SetBackgroundValue(0); // Must be zero
47 SetForegroundValue(1);
48 SetInitialSmoothing(false);
50 SetSmoothingConductanceParameter(3.0);
51 SetSmoothingNumberOfIterations(5);
52 SetSmoothingTimeStep(0.0625);
53 SetSmoothingUseImageSpacing(false);
55 SetMinimalComponentSize(100);
56 SetUpperThreshold1(1500);
57 SetLowerThreshold1(100);
58 SetFullConnectivity(false);
60 SetUpperThreshold2(1500);
61 SetLowerThreshold2(10);
69 //--------------------------------------------------------------------
72 //--------------------------------------------------------------------
73 template <class TInputImageType>
75 clitk::ExtractBonesFilter<TInputImageType>::
76 SetInput(const TInputImageType * image)
78 this->SetNthInput(0, const_cast<TInputImageType *>(image));
80 //--------------------------------------------------------------------
84 //--------------------------------------------------------------------
85 template <class TInputImageType>
87 clitk::ExtractBonesFilter<TInputImageType>::
88 GenerateOutputInformation() {
91 InputImagePointer input = dynamic_cast<TInputImageType*>(itk::ProcessObject::GetInput(0));
92 Superclass::GenerateOutputInformation();
93 MaskImagePointer outputImage = this->GetOutput(0);
94 outputImage->SetRegions(input->GetLargestPossibleRegion());
100 typedef itk::BinaryThresholdImageFilter<InputImageType, InternalImageType> InputBinarizeFilterType;
101 typedef itk::BinaryThresholdImageFilter<InternalImageType, InternalImageType> BinarizeFilterType;
102 typedef itk::ConnectedComponentImageFilter<InternalImageType, InternalImageType> ConnectFilterType;
103 typedef itk::RelabelComponentImageFilter<InternalImageType, InternalImageType> RelabelFilterType;
104 typedef clitk::SetBackgroundImageFilter<InternalImageType,InternalImageType, InternalImageType> SetBackgroundFilterType;
105 typedef itk::CastImageFilter<InternalImageType,MaskImageType> CastImageFilterType;
106 typedef itk::ImageFileWriter<MaskImageType> WriterType;
108 //---------------------------------
109 // Smoothing [Optional]
110 //---------------------------------
111 if (GetInitialSmoothing()) {
112 StartNewStep("Initial Smoothing");
113 typedef itk::CurvatureAnisotropicDiffusionImageFilter<InputImageType, InputImageType> FilterType;
114 typename FilterType::Pointer df = FilterType::New();
115 df->SetConductanceParameter(GetSmoothingConductanceParameter());
116 df->SetNumberOfIterations(GetSmoothingNumberOfIterations());
117 df->SetTimeStep(GetSmoothingTimeStep());
118 df->SetUseImageSpacing(GetSmoothingUseImageSpacing());
121 filtered_input = df->GetOutput();
122 StopCurrentStep<InputImageType>(filtered_input);
125 filtered_input = input;
128 //--------------------------------------------------------------------
129 //--------------------------------------------------------------------
130 StartNewStep("Initial Labeling");
131 typename InternalImageType::Pointer firstLabelImage;
133 //---------------------------------
134 // Binarize the image
135 //---------------------------------
136 typename InputBinarizeFilterType::Pointer binarizeFilter=InputBinarizeFilterType::New();
137 binarizeFilter->SetInput(filtered_input);
138 binarizeFilter->SetLowerThreshold(GetLowerThreshold1());
139 binarizeFilter->SetUpperThreshold(GetUpperThreshold1());
140 binarizeFilter->SetInsideValue(this->GetForegroundValue());
141 binarizeFilter->SetOutsideValue(this->GetBackgroundValue());
143 //---------------------------------
144 // Label the connected components
145 //---------------------------------
146 typename ConnectFilterType::Pointer connectFilter=ConnectFilterType::New();
147 connectFilter->SetInput(binarizeFilter->GetOutput());
148 connectFilter->SetBackgroundValue(this->GetBackgroundValue());
149 connectFilter->SetFullyConnected(GetFullConnectivity());
151 //---------------------------------
152 // Sort the labels according to size
153 //---------------------------------
154 typename RelabelFilterType::Pointer relabelFilter=RelabelFilterType::New();
155 relabelFilter->SetInput(connectFilter->GetOutput());
156 relabelFilter->SetMinimumObjectSize(GetMinimalComponentSize());
158 //---------------------------------
160 //---------------------------------
161 typename BinarizeFilterType::Pointer binarizeFilter2=BinarizeFilterType::New();
162 binarizeFilter2->SetInput(relabelFilter->GetOutput());
163 binarizeFilter2->SetLowerThreshold(1);
164 binarizeFilter2->SetUpperThreshold(1);
165 binarizeFilter2->SetInsideValue(this->GetForegroundValue());
166 binarizeFilter2->SetOutsideValue(this->GetBackgroundValue());
167 binarizeFilter2->Update();
169 firstLabelImage = binarizeFilter2->GetOutput();
170 StopCurrentStep<InternalImageType>(firstLabelImage);
172 //--------------------------------------------------------------------
173 //--------------------------------------------------------------------
174 StartNewStep("Neighborhood connected filter");
175 typename InternalImageType::Pointer secondLabelImage;
177 //---------------------------------
178 //Neighborhood connected RG
179 //---------------------------------
180 typedef itk::NeighborhoodConnectedImageFilter<InputImageType, InternalImageType>
181 NeighborhoodConnectedImageFilterType;
182 typename NeighborhoodConnectedImageFilterType::Pointer neighborhoodConnectedImageFilter=
183 NeighborhoodConnectedImageFilterType::New();
186 neighborhoodConnectedImageFilter->SetLower(GetLowerThreshold2());
187 neighborhoodConnectedImageFilter->SetUpper(GetUpperThreshold2());
188 neighborhoodConnectedImageFilter->SetReplaceValue(this->GetForegroundValue());
189 neighborhoodConnectedImageFilter->SetRadius(GetRadius2());
190 neighborhoodConnectedImageFilter->SetInput(filtered_input);
192 // Seeds from label image
193 typedef itk::ImageRegionIteratorWithIndex<InternalImageType> IteratorType;
194 IteratorType it(firstLabelImage, firstLabelImage->GetLargestPossibleRegion());
195 typename InputImageType::IndexType index;
196 unsigned int counter=0;
197 while (!it.IsAtEnd())
199 if (it.Get()==this->GetForegroundValue())
203 neighborhoodConnectedImageFilter->AddSeed(index);
206 while (!it.IsAtEnd() && i< (unsigned int) GetSampleRate2())
215 neighborhoodConnectedImageFilter->Update();
216 secondLabelImage = neighborhoodConnectedImageFilter->GetOutput();
218 //--------------------------------------------------------------------
219 //--------------------------------------------------------------------
220 StartNewStep("Combine the images");
221 typedef clitk::SetBackgroundImageFilter<InternalImageType, InternalImageType, InternalImageType>
222 SetBackgroundImageFilterType;
223 typename SetBackgroundImageFilterType::Pointer setBackgroundFilter=SetBackgroundImageFilterType::New();
224 setBackgroundFilter->SetInput(firstLabelImage);
225 setBackgroundFilter->SetInput2(secondLabelImage);
226 setBackgroundFilter->SetMaskValue(this->GetForegroundValue());
227 setBackgroundFilter->SetOutsideValue(this->GetForegroundValue());
228 setBackgroundFilter->Update();
230 output = setBackgroundFilter->GetOutput();
232 //--------------------------------------------------------------------
233 //--------------------------------------------------------------------
235 if (GetFillHoles()) {
236 StartNewStep("Fill Holes");
237 typedef clitk::FillMaskFilter<InternalImageType> FillMaskFilterType;
238 typename FillMaskFilterType::Pointer fillMaskFilter = FillMaskFilterType::New();
239 fillMaskFilter->SetInput(output);
240 fillMaskFilter->AddDirection(2);
241 fillMaskFilter->Update();
242 output = fillMaskFilter->GetOutput();
243 StopCurrentStep<InternalImageType>(output);
246 //--------------------------------------------------------------------
247 //--------------------------------------------------------------------
250 StartNewStep("AutoCrop");
251 typedef clitk::AutoCropFilter<InternalImageType> CropFilterType;
252 typename CropFilterType::Pointer cropFilter = CropFilterType::New();
253 cropFilter->SetInput(output);
254 cropFilter->SetBackgroundValue(GetBackgroundValue());
255 cropFilter->Update();
256 output = cropFilter->GetOutput();
257 StopCurrentStep<InternalImageType>(output);
258 outputImage->SetRegions(output->GetLargestPossibleRegion());
262 //--------------------------------------------------------------------
265 //--------------------------------------------------------------------
266 template <class TInputImageType>
268 clitk::ExtractBonesFilter<TInputImageType>::
271 //--------------------------------------------------------------------
272 //--------------------------------------------------------------------
274 typedef itk::CastImageFilter<InternalImageType, MaskImageType> CastImageFilterType;
275 typename CastImageFilterType::Pointer caster= CastImageFilterType::New();
276 caster->SetInput(output);
278 this->GraftOutput(caster->GetOutput());
280 // Store image filenames into AFDB
281 GetAFDB()->SetImageFilename("Bones", this->GetOutputBonesFilename());
285 //--------------------------------------------------------------------
288 #endif //#define CLITKBOOLEANOPERATORLABELIMAGEFILTER_TXX