+
+bool is_orientation_superior(itk::SpatialOrientation::ValidCoordinateOrientationFlags orientation)
+{
+ itk::SpatialOrientation::CoordinateTerms sup = itk::SpatialOrientation::ITK_COORDINATE_Superior;
+ std::cout << "orientation: " << std::hex << orientation << "; superior: " << std::hex << sup << std::endl;
+ std::cout << std::dec;
+
+ 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 (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;
+ 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());
+
+ typename LabelMapToImageFilterType::Pointer map_to_label_filter = LabelMapToImageFilterType::New();
+ typename FileWriterType::Pointer writer = FileWriterType::New();
+
+ typename ImageType::IndexType index;
+ typename ImageType::RegionType region = working_input->GetLargestPossibleRegion();
+ typename ImageType::SizeType size = region.GetSize();
+ typename ImageType::SpacingType spacing = working_input->GetSpacing();
+ typename ImageType::PointType origin = working_input->GetOrigin();
+
+ int nslices = min(numberOfSlices, size[2]);
+ int start = 0, increment = 1;
+ itk::SpatialOrientationAdapter orientation;
+ typename ImageType::DirectionType dir = working_input->GetDirection();
+ if (!is_orientation_superior(orientation.FromDirectionCosines(dir))) {
+ start = size[2]-1;
+ increment = -1;
+ }
+
+ typename MaskImageType::PointType image_centre;
+ image_centre[0] = size[0]/2;
+ image_centre[1] = size[1]/2;
+ image_centre[2] = 0;