+ // output = working_mask;
+ //this->GetOutput(0)->SetRegions(output->GetLargestPossibleRegion());
+
+ // this->GetOutput(0)->SetRegions(working_mask->GetLargestPossibleRegion());
+
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template <class TImageType>
+void
+clitk::ExtractLungFilter<TImageType>::
+GenerateInputRequestedRegion() {
+ // DD("GenerateInputRequestedRegion (nothing?)");
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template <class ImageType>
+void
+clitk::ExtractLungFilter<ImageType>::
+GenerateData()
+{
+ // Set the output
+ // this->GraftOutput(output); // not SetNthOutput
+ this->GraftOutput(working_mask); // not SetNthOutput
+ // Store image filenames into AFDB
+ GetAFDB()->SetImageFilename("Lungs", this->GetOutputLungFilename());
+ GetAFDB()->SetImageFilename("Trachea", this->GetOutputTracheaFilename());
+ WriteAFDB();
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template <class ImageType>
+bool
+clitk::ExtractLungFilter<ImageType>::
+SearchForTracheaSeed(int skip)
+{
+ if (m_Seeds.size() == 0) { // try to find seed (if not zero, it is given by user)
+ // Restart the search until a seed is found, skipping more and more slices
+ bool stop = false;
+ while (!stop) {
+ // Search seed (parameters = UpperThresholdForTrachea)
+ static const unsigned int Dim = ImageType::ImageDimension;
+ typename InternalImageType::RegionType sliceRegion = working_input->GetLargestPossibleRegion();
+ typename InternalImageType::SizeType sliceRegionSize = sliceRegion.GetSize();
+ typename InternalImageType::IndexType sliceRegionIndex = sliceRegion.GetIndex();
+ sliceRegionIndex[Dim-1]=sliceRegionSize[Dim-1]-skip-5;
+ sliceRegionSize[Dim-1]=5;
+ sliceRegion.SetSize(sliceRegionSize);
+ sliceRegion.SetIndex(sliceRegionIndex);
+ typedef itk::ImageRegionConstIterator<ImageType> IteratorType;
+ IteratorType it(working_input, sliceRegion);
+ it.GoToBegin();
+ while (!it.IsAtEnd()) {
+ if(it.Get() < GetUpperThresholdForTrachea() ) {
+ AddSeed(it.GetIndex());
+ // DD(it.GetIndex());
+ }
+ ++it;
+ }
+
+ // if we do not found : restart
+ stop = (m_Seeds.size() != 0);
+ if (!stop) {
+ if (GetVerboseStepFlag()) {
+ std::cout << "\t No seed found this time. I skip some slices and restart." << std::endl;
+ }
+ if (skip > 0.5 * working_input->GetLargestPossibleRegion().GetSize()[2]) {
+ // we want to skip more than a half of the image, it is probably a bug
+ std::cerr << "Number of slices to skip to find trachea to high = " << skip << std::endl;
+ stop = true;
+ }
+ skip += 5;
+ }
+ else {
+ // DD(m_Seeds[0]);
+ // DD(m_Seeds.size());
+ }
+ }
+ }
+ return (m_Seeds.size() != 0);
+}
+//--------------------------------------------------------------------
+
+
+bool is_orientation_superior(itk::SpatialOrientation::ValidCoordinateOrientationFlags orientation)
+{
+ itk::SpatialOrientation::CoordinateTerms sup = itk::SpatialOrientation::ITK_COORDINATE_Superior;
+ bool primary = (orientation & 0x0000ff) == sup;
+ bool secondary = ((orientation & 0x00ff00) >> 8) == sup;
+ bool tertiary = ((orientation & 0xff0000) >> 16) == sup;
+ return primary || secondary || tertiary;
+}
+
+//--------------------------------------------------------------------
+template <class ImageType>
+bool
+clitk::ExtractLungFilter<ImageType>::
+SearchForTracheaSeed2(int numberOfSlices)
+{
+ if (m_Seeds.size() == 0) { // try to find seed (if not zero, it is given by user)
+ if (GetVerboseRegionGrowingFlag())
+ 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;
+ typedef itk::BinaryThresholdImageFilter<ImageType, MaskImageType> ThresholdFilterType;
+ typedef itk::BinaryBallStructuringElement<MaskPixelType, MaskImageType::ImageDimension> KernelType;
+ typedef itk::BinaryMorphologicalClosingImageFilter<MaskImageType, MaskImageType, KernelType> ClosingFilterType;
+ typedef itk::BinaryMorphologicalOpeningImageFilter<MaskImageType, MaskImageType, KernelType> OpeningFilterType;
+ typedef itk::ExtractImageFilter<MaskImageType, MaskImageType2D> SlicerFilterType;
+ typedef itk::ConnectedComponentImageFilter<MaskImageType2D, MaskImageType2D> LabelFilterType;
+ typedef itk::ShapeLabelObject<MaskPixelType, MaskImageType2D::ImageDimension> ShapeLabelType;
+ typedef itk::LabelMap<ShapeLabelType> LabelImageType;
+ typedef itk::LabelImageToShapeLabelMapFilter<MaskImageType2D, LabelImageType> ImageToLabelMapFilterType;
+ typedef itk::LabelMapToLabelImageFilter<LabelImageType, MaskImageType2D> LabelMapToImageFilterType;
+ typedef itk::ImageFileWriter<MaskImageType2D> FileWriterType;
+
+ // threshold to isolate airawys and lungs
+ typename ThresholdFilterType::Pointer threshold = ThresholdFilterType::New();
+ threshold->SetLowerThreshold(-2000);
+ threshold->SetUpperThreshold(GetSeedPreProcessingThreshold());
+ threshold->SetInput(working_input);
+ threshold->Update();
+
+ KernelType kernel_closing, kernel_opening;
+
+ // remove small noise
+ typename OpeningFilterType::Pointer opening = OpeningFilterType::New();
+ kernel_opening.SetRadius(1);
+ opening->SetKernel(kernel_opening);
+ opening->SetInput(threshold->GetOutput());
+ opening->Update();
+
+ typename SlicerFilterType::Pointer slicer = SlicerFilterType::New();
+ slicer->SetInput(opening->GetOutput());
+
+ // label result
+ typename LabelFilterType::Pointer label_filter = LabelFilterType::New();
+ label_filter->SetInput(slicer->GetOutput());
+
+ // extract shape information from labels
+ typename ImageToLabelMapFilterType::Pointer label_to_map_filter = ImageToLabelMapFilterType::New();
+ label_to_map_filter->SetInput(label_filter->GetOutput());