Authors belong to:
- University of LYON http://www.universite-lyon.fr/
- - Léon Bérard cancer center http://oncora1.lyon.fnclcc.fr
+ - Léon Bérard cancer center http://www.centreleonberard.fr
- CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
This software is distributed WITHOUT ANY WARRANTY; without even
- BSD See included LICENSE.txt file
- CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
- ======================================================================-====*/
+ ===========================================================================**/
#ifndef CLITKEXTRACTMEDIASTINUMFILTER_TXX
#define CLITKEXTRACTMEDIASTINUMFILTER_TXX
#include "clitkCommon.h"
#include "clitkExtractMediastinumFilter.h"
#include "clitkAddRelativePositionConstraintToLabelImageFilter.h"
+#include "clitkSliceBySliceRelativePositionFilter.h"
#include "clitkSegmentationUtils.h"
+#include "clitkExtractAirwaysTreeInfoFilter.h"
+#include "clitkCropLikeImageFilter.h"
-// itk
+// std
#include <deque>
+
+// itk
#include "itkStatisticsLabelMapFilter.h"
#include "itkLabelImageToStatisticsLabelMapFilter.h"
#include "itkRegionOfInterestImageFilter.h"
#include "itkBinaryThresholdImageFilter.h"
+#include "itkScalarImageKmeansImageFilter.h"
// itk ENST
#include "RelativePositionPropImageFilter.h"
//--------------------------------------------------------------------
-template <class TImageType>
-clitk::ExtractMediastinumFilter<TImageType>::
+template <class ImageType>
+clitk::ExtractMediastinumFilter<ImageType>::
ExtractMediastinumFilter():
- clitk::FilterBase(),
- itk::ImageToImageFilter<TImageType, TImageType>()
+ clitk::StructuresExtractionFilter<ImageType>()
{
- this->SetNumberOfRequiredInputs(3);
-
+ this->SetNumberOfRequiredInputs(1);
SetBackgroundValuePatient(0);
SetBackgroundValueLung(0);
SetBackgroundValueBones(0);
SetForegroundValueLeftLung(1);
SetForegroundValueRightLung(2);
- SetBackgroundValue(0);
- SetForegroundValue(1);
-
- SetIntermediateSpacing(6);
- SetFuzzyThreshold1(0.6);
- SetFuzzyThreshold2(0.7);
+ SetDistanceMaxToAnteriorPartOfTheVertebralBody(10);
+ SetOutputMediastinumFilename("mediastinum.mha");
+ UseBonesOff();
}
//--------------------------------------------------------------------
//--------------------------------------------------------------------
-template <class TImageType>
+template <class ImageType>
void
-clitk::ExtractMediastinumFilter<TImageType>::
-SetInputPatientLabelImage(const TImageType * image, ImagePixelType bg) {
- this->SetNthInput(0, const_cast<TImageType *>(image));
+clitk::ExtractMediastinumFilter<ImageType>::
+SetInputPatientLabelImage(const MaskImageType * image, MaskImagePixelType bg)
+{
+ this->SetNthInput(0, const_cast<MaskImageType *>(image));
m_BackgroundValuePatient = bg;
}
//--------------------------------------------------------------------
//--------------------------------------------------------------------
-template <class TImageType>
+template <class ImageType>
void
-clitk::ExtractMediastinumFilter<TImageType>::
-SetInputLungLabelImage(const TImageType * image, ImagePixelType bg,
- ImagePixelType fgLeft, ImagePixelType fgRight) {
- this->SetNthInput(1, const_cast<TImageType *>(image));
+clitk::ExtractMediastinumFilter<ImageType>::
+SetInputLungLabelImage(const MaskImageType * image, MaskImagePixelType bg,
+ MaskImagePixelType fgLeft, MaskImagePixelType fgRight)
+{
+ this->SetNthInput(1, const_cast<MaskImageType *>(image));
m_BackgroundValueLung = bg;
m_ForegroundValueLeftLung = fgLeft;
m_ForegroundValueRightLung = fgRight;
//--------------------------------------------------------------------
-template <class TImageType>
+template <class ImageType>
void
-clitk::ExtractMediastinumFilter<TImageType>::
-SetInputBonesLabelImage(const TImageType * image, ImagePixelType bg) {
- this->SetNthInput(2, const_cast<TImageType *>(image));
+clitk::ExtractMediastinumFilter<ImageType>::
+SetInputBonesLabelImage(const MaskImageType * image, MaskImagePixelType bg)
+{
+ this->SetNthInput(2, const_cast<MaskImageType *>(image));
m_BackgroundValueBones = bg;
}
//--------------------------------------------------------------------
//--------------------------------------------------------------------
-template <class TImageType>
-template<class ArgsInfoType>
+template <class ImageType>
void
-clitk::ExtractMediastinumFilter<TImageType>::
-SetArgsInfo(ArgsInfoType mArgsInfo)
+clitk::ExtractMediastinumFilter<ImageType>::
+SetInputTracheaLabelImage(const MaskImageType * image, MaskImagePixelType bg)
{
- SetVerboseOption_GGO(mArgsInfo);
- SetVerboseStep_GGO(mArgsInfo);
- SetWriteStep_GGO(mArgsInfo);
- SetVerboseWarningOff_GGO(mArgsInfo);
+ this->SetNthInput(3, const_cast<MaskImageType *>(image));
+ m_BackgroundValueTrachea = bg;
+}
+//--------------------------------------------------------------------
- SetBackgroundValuePatient_GGO(mArgsInfo);
- SetBackgroundValueLung_GGO(mArgsInfo);
- SetBackgroundValueBones_GGO(mArgsInfo);
+//--------------------------------------------------------------------
+template <class ImageType>
+void
+clitk::ExtractMediastinumFilter<ImageType>::
+SetFuzzyThreshold(std::string tag, double value)
+{
+ m_FuzzyThreshold[tag] = value;
+}
+//--------------------------------------------------------------------
- SetForegroundValueLeftLung_GGO(mArgsInfo);
- SetForegroundValueRightLung_GGO(mArgsInfo);
- SetIntermediateSpacing_GGO(mArgsInfo);
- SetFuzzyThreshold1_GGO(mArgsInfo);
- SetFuzzyThreshold2_GGO(mArgsInfo);
+//--------------------------------------------------------------------
+template <class ImageType>
+double
+clitk::ExtractMediastinumFilter<ImageType>::
+GetFuzzyThreshold(std::string tag)
+{
+ if (m_FuzzyThreshold.find(tag) == m_FuzzyThreshold.end()) {
+ clitkExceptionMacro("Could not find options "+tag+" in the list of FuzzyThresholds.");
+ return 0.0;
+ }
+
+ return m_FuzzyThreshold[tag];
}
//--------------------------------------------------------------------
//--------------------------------------------------------------------
-template <class TImageType>
+template <class ImageType>
void
-clitk::ExtractMediastinumFilter<TImageType>::
-GenerateOutputInformation() {
- ImagePointer input = dynamic_cast<TImageType*>(itk::ProcessObject::GetInput(0));
- ImagePointer outputImage = this->GetOutput(0);
- outputImage->SetRegions(outputImage->GetLargestPossibleRegion());
+clitk::ExtractMediastinumFilter<ImageType>::
+GenerateInputRequestedRegion()
+{
+ // DD("GenerateInputRequestedRegion");
+ // Do not call default
+ // Superclass::GenerateInputRequestedRegion();
+ // DD("End GenerateInputRequestedRegion");
}
//--------------------------------------------------------------------
//--------------------------------------------------------------------
-template <class TImageType>
+template <class ImageType>
void
-clitk::ExtractMediastinumFilter<TImageType>::
-GenerateInputRequestedRegion() {
- // Call default
- Superclass::GenerateInputRequestedRegion();
- // Get input pointers
- ImagePointer patient = dynamic_cast<TImageType*>(itk::ProcessObject::GetInput(0));
- ImagePointer lung = dynamic_cast<TImageType*>(itk::ProcessObject::GetInput(1));
- ImagePointer bones = dynamic_cast<TImageType*>(itk::ProcessObject::GetInput(2));
-
- patient->SetRequestedRegion(patient->GetLargestPossibleRegion());
- lung->SetRequestedRegion(lung->GetLargestPossibleRegion());
- bones->SetRequestedRegion(bones->GetLargestPossibleRegion());
+clitk::ExtractMediastinumFilter<ImageType>::
+SetInput(const ImageType * image)
+{
+ this->SetNthInput(0, const_cast<ImageType *>(image));
}
//--------------------------------------------------------------------
+
//--------------------------------------------------------------------
-template <class TImageType>
+template <class ImageType>
void
-clitk::ExtractMediastinumFilter<TImageType>::
-GenerateData() {
- // Get input pointers
- ImageConstPointer patient = dynamic_cast<const TImageType*>(itk::ProcessObject::GetInput(0));
- ImageConstPointer lung = dynamic_cast<const TImageType*>(itk::ProcessObject::GetInput(1));
- ImageConstPointer bones = dynamic_cast<const TImageType*>(itk::ProcessObject::GetInput(2));
-
- // Get output pointer
- ImagePointer output;
+clitk::ExtractMediastinumFilter<ImageType>::
+GenerateOutputInformation() {
+ // Do not call default
+ // Superclass::GenerateOutputInformation();
- // Step 1: patient minus lungs, minus bones
- StartNewStep("Patient contours minus lungs and minus bones");
- typedef clitk::BooleanOperatorLabelImageFilter<ImageType> BoolFilterType;
+ //--------------------------------------------------------------------
+ // Get input pointers
+ this->LoadAFDB();
+ ImageConstPointer input = dynamic_cast<const ImageType*>(itk::ProcessObject::GetInput(0));
+ MaskImagePointer patient, lung, bones, trachea;
+ patient = this->GetAFDB()->template GetImage <MaskImageType>("Patient");
+ lung = this->GetAFDB()->template GetImage <MaskImageType>("Lungs");
+ if (GetUseBones()) {
+ bones = this->GetAFDB()->template GetImage <MaskImageType>("Bones");
+ }
+ trachea = this->GetAFDB()->template GetImage <MaskImageType>("Trachea");
+
+ //this->VerboseImageSizeFlagOn();
+
+ //--------------------------------------------------------------------
+ // Step : Crop support (patient) to lung extend in RL
+ this->StartNewStep("[Mediastinum] Crop support like lungs along LR");
+ typedef clitk::CropLikeImageFilter<MaskImageType> CropFilterType;
+ typename CropFilterType::Pointer cropFilter = CropFilterType::New();
+ cropFilter->SetInput(patient);
+ cropFilter->SetCropLikeImage(lung, 0);// Indicate that we only crop in X (Left-Right) axe
+ cropFilter->Update();
+ output = cropFilter->GetOutput();
+ this->template StopCurrentStep<MaskImageType>(output);
+
+ //--------------------------------------------------------------------
+ // Step : Crop support (previous) to bones extend in AP
+ if (GetUseBones()) {
+ this->StartNewStep("[Mediastinum] Crop support like bones along AP");
+ cropFilter = CropFilterType::New();
+ cropFilter->SetInput(output);
+ cropFilter->SetCropLikeImage(bones, 1);// Indicate that we only crop in Y (Ant-Post) axe
+ cropFilter->Update();
+ output = cropFilter->GetOutput();
+ this->template StopCurrentStep<MaskImageType>(output);
+ }
+
+ //--------------------------------------------------------------------
+ // Step : patient minus lungs, minus bones, minus trachea
+ this->StartNewStep("[Mediastinum] Patient contours minus lungs, trachea [and bones]");
+ typedef clitk::BooleanOperatorLabelImageFilter<MaskImageType> BoolFilterType;
typename BoolFilterType::Pointer boolFilter = BoolFilterType::New();
boolFilter->InPlaceOn();
- boolFilter->SetInput1(patient);
+ boolFilter->SetInput1(output);
boolFilter->SetInput2(lung);
boolFilter->SetOperationType(BoolFilterType::AndNot);
boolFilter->Update();
+ if (GetUseBones()) {
+ boolFilter->SetInput1(boolFilter->GetOutput());
+ boolFilter->SetInput2(bones);
+ boolFilter->SetOperationType(BoolFilterType::AndNot);
+ boolFilter->Update();
+ }
boolFilter->SetInput1(boolFilter->GetOutput());
- boolFilter->SetInput2(bones);
+ boolFilter->SetInput2(trachea);
boolFilter->SetOperationType(BoolFilterType::AndNot);
boolFilter->Update();
output = boolFilter->GetOutput();
- output = clitk::AutoCrop<ImageType>(output, GetBackgroundValue());
- ////autoCropFilter->GetOutput(); typedef clitk::AutoCropFilter<ImageType> CropFilterType;
- //typename CropFilterType::Pointer cropFilter = CropFilterType::New();
- //cropFilter->SetInput(output);
- //cropFilter->Update();
- //output = cropFilter->GetOutput();
-
- this->template StopCurrentStep<TImageType>(output);
-
- // Step 2: LR limits from lung (need separate lung ?)
- StartNewStep("Left/Right limits with lungs");
- typedef clitk::AddRelativePositionConstraintToLabelImageFilter<TImageType> RelPosFilterType;
- typename RelPosFilterType::Pointer relPosFilter = RelPosFilterType::New();
- relPosFilter->SetCurrentStepBaseId(this->GetCurrentStepId());
- relPosFilter->VerboseStepOff();
- relPosFilter->WriteStepOff();
- relPosFilter->SetInput(output);
- relPosFilter->SetInputObject(lung);
- relPosFilter->SetOrientationType(RelPosFilterType::LeftTo);
- relPosFilter->SetIntermediateSpacing(GetIntermediateSpacing());
- relPosFilter->SetFuzzyThreshold(GetFuzzyThreshold1());
- relPosFilter->Update();
- output = relPosFilter->GetOutput();
-
- relPosFilter->SetInput(output);
- relPosFilter->SetCurrentStepBaseId(this->GetCurrentStepId());
- relPosFilter->VerboseStepOff();
- relPosFilter->WriteStepOff();
- relPosFilter->SetInputObject(lung);
- relPosFilter->SetOrientationType(RelPosFilterType::RightTo);
- relPosFilter->SetIntermediateSpacing(GetIntermediateSpacing());
- relPosFilter->SetFuzzyThreshold(GetFuzzyThreshold1());
- relPosFilter->Update();
- output = relPosFilter->GetOutput();
- this->template StopCurrentStep<TImageType>(output);
-
- // Step 3: AP limits from bones
- StartNewStep("Ant/Post limits with bones");
- relPosFilter->SetCurrentStepNumber(0);
- relPosFilter->ResetPipeline();// = RelPosFilterType::New();
- relPosFilter->SetCurrentStepBaseId(this->GetCurrentStepId());
- relPosFilter->VerboseStepOff();
- relPosFilter->WriteStepOff();
- relPosFilter->SetInput(output);
- relPosFilter->SetInputObject(bones);
- relPosFilter->SetOrientationType(RelPosFilterType::AntTo);
- relPosFilter->SetIntermediateSpacing(GetIntermediateSpacing());
- relPosFilter->SetFuzzyThreshold(GetFuzzyThreshold2());
- relPosFilter->Update();
-
- relPosFilter->SetInput(relPosFilter->GetOutput());
- relPosFilter->SetCurrentStepBaseId(this->GetCurrentStepId());
- relPosFilter->VerboseStepOff();
- relPosFilter->WriteStepOff();
- relPosFilter->SetInputObject(bones);
- relPosFilter->SetOrientationType(RelPosFilterType::PostTo);
- relPosFilter->SetIntermediateSpacing(GetIntermediateSpacing());
- relPosFilter->SetFuzzyThreshold(GetFuzzyThreshold2());
- relPosFilter->Update();
- output = relPosFilter->GetOutput();
- this->template StopCurrentStep<TImageType>(output);
-
- // Get CCL
- output = clitk::Labelize<TImageType>(output, GetBackgroundValue(), true, 100);
- // output = RemoveLabels<TImageType>(output, BG, param->GetLabelsToRemove());
- output = clitk::KeepLabels<TImageType>(output, GetBackgroundValue(),
- GetForegroundValue(), 1, 1, 0);
-
- output = clitk::AutoCrop<ImageType>(output, GetBackgroundValue());
- // cropFilter = CropFilterType::New();
- //cropFilter->SetInput(output);
- //cropFilter->Update();
- //output = cropFilter->GetOutput();
-
- // Final Step -> set output
- this->SetNthOutput(0, output);
+ // Auto crop to gain support area
+ output = clitk::AutoCrop<MaskImageType>(output, this->GetBackgroundValue());
+ this->template StopCurrentStep<MaskImageType>(output);
+
+ //--------------------------------------------------------------------
+ // Step : LR limits from lung (need separate lung ?)
+ // Get separate lung images to get only the right and left lung
+ // (because RelativePositionPropImageFilter only consider fg=1);
+ // (label must be '1' because right is greater than left). (WE DO
+ // NOT NEED TO SEPARATE ? )
+ this->StartNewStep("[Mediastinum] Left/Right limits with lungs");
+
+ // The following cannot be "inplace" because mask is the same than input ...
+ MaskImagePointer right_lung =
+ clitk::SetBackground<MaskImageType, MaskImageType>(lung, lung, 2, 0, false);
+ MaskImagePointer left_lung =
+ clitk::SetBackground<MaskImageType, MaskImageType>(lung, lung, 1, 0, false);
+ left_lung = clitk::SetBackground<MaskImageType, MaskImageType>(left_lung, left_lung, 2, 1, false);
+ right_lung = clitk::ResizeImageLike<MaskImageType>(right_lung, output, this->GetBackgroundValue());
+ left_lung = clitk::ResizeImageLike<MaskImageType>(left_lung, output, this->GetBackgroundValue());
+ this->GetAFDB()->template SetImage<MaskImageType>("RightLung", "seg/RightLung.mha",
+ right_lung, true);
+ this->GetAFDB()->template SetImage<MaskImageType>("LeftLung", "seg/LeftLung.mha",
+ left_lung, true);
+ this->GetAFDB()->Write();
+ this->template StopCurrentStep<MaskImageType>(output);
+
+ //--------------------------------------------------------------------
+ // Step : AP limits from bones
+ // Separate the bones in the ant-post middle
+ MaskImagePointer bones_ant;
+ MaskImagePointer bones_post;
+ MaskImagePointType middle_AntPost_position;
+ middle_AntPost_position.Fill(NumericTraits<MaskImagePointType::ValueType>::Zero);
+ if (GetUseBones()) {
+ this->StartNewStep("[Mediastinum] Ant/Post limits with bones");
+
+ // To define ant and post part of the bones with a single horizontal line
+ MaskImageIndexType index_bones_middle;
+
+ // Method1: cut in the middle (not optimal)
+ /*
+ middle_AntPost_position[0] = middle_AntPost_position[2] = 0;
+ middle_AntPost_position[1] = bones->GetOrigin()[1]+
+ (bones->GetLargestPossibleRegion().GetSize()[1]*bones->GetSpacing()[1])/2.0;
+ DD(middle_AntPost_position);
+ bones->TransformPhysicalPointToIndex(middle_AntPost_position, index_bones_middle);
+ */
+
+ // Method2: Use VertebralBody, take most ant point
+ MaskImagePointer VertebralBody = this->GetAFDB()->template GetImage<MaskImageType>("VertebralBody");
+ FindExtremaPointInAGivenDirection<MaskImageType>(VertebralBody, this->GetBackgroundValue(),
+ 1, true, middle_AntPost_position);
+ bones->TransformPhysicalPointToIndex(middle_AntPost_position, index_bones_middle);
+
+ // Split bone image first into two parts (ant and post), and crop
+ // lateraly to get vertebral
+ typedef itk::RegionOfInterestImageFilter<MaskImageType, MaskImageType> ROIFilterType;
+ // typedef itk::ExtractImageFilter<MaskImageType, MaskImageType> ROIFilterType;
+ typename ROIFilterType::Pointer roiFilter = ROIFilterType::New();
+ MaskImageRegionType region = bones->GetLargestPossibleRegion();
+ MaskImageSizeType size = region.GetSize();
+ MaskImageIndexType index = region.GetIndex();
+ // ANT part
+ // crop LR to keep 1/4 center part
+ // index[0] = size[0]/4+size[0]/8;
+ //size[0] = size[0]/4;
+ // crop AP to keep first (ant) part
+ size[1] = index_bones_middle[1]; //size[1]/2.0;
+ region.SetSize(size);
+ region.SetIndex(index);
+ roiFilter->SetInput(bones);
+ roiFilter->SetRegionOfInterest(region);
+ roiFilter->ReleaseDataFlagOff();
+ roiFilter->Update();
+ bones_ant = roiFilter->GetOutput();
+ // writeImage<MaskImageType>(bones_ant, "b_ant.mhd");
+ // POST part
+ roiFilter = ROIFilterType::New();
+ index[1] = bones->GetLargestPossibleRegion().GetIndex()[1] + size[1]-1;
+ size[1] = bones->GetLargestPossibleRegion().GetSize()[1] - size[1];
+ region.SetIndex(index);
+ region.SetSize(size);
+ roiFilter->SetInput(bones);
+ roiFilter->SetRegionOfInterest(region);
+ roiFilter->ReleaseDataFlagOff();
+ roiFilter->Update();
+ bones_post = roiFilter->GetOutput();
+ // writeImage<MaskImageType>(bones_post, "b_post.mhd");
+
+ // Now, insert this image in the AFDB ==> (needed because used in the RelativePosition config file)
+ this->GetAFDB()->template SetImage<MaskImageType>("Bones_Post", "seg/Bones_Post.mha",
+ bones_post, true);
+ this->GetAFDB()->template SetImage<MaskImageType>("Bones_Ant", "seg/Bones_Ant.mha",
+ bones_ant, true);
+ this->GetAFDB()->Write();
+
+ this->template StopCurrentStep<MaskImageType>(output);
+ }
+
+ //--------------------------------------------------------------------
+ // Remove VertebralBody part
+ this->StartNewStep("[Mediastinum] Remove VertebralBody");
+ MaskImagePointer VertebralBody = this->GetAFDB()->template GetImage<MaskImageType>("VertebralBody");
+ clitk::AndNot<MaskImageType>(output, VertebralBody, this->GetBackgroundValue());
+ this->template StopCurrentStep<MaskImageType>(output);
+
+ //--------------------------------------------------------------------
+ // Generic RelativePosition processes
+ output = this->ApplyRelativePositionList("Mediastinum", output);
+
+
+ //--------------------------------------------------------------------
+ // FIXME --> do not put this limits here !
+ /*
+ // Step : SI limits It is better to do this limit *AFTER* the
+ // RelativePosition to avoid some issue due to superior boundaries.
+ this->StartNewStep("[Mediastinum] Keep inferior to CricoidCartilag");
+ // load Cricoid, get centroid, cut above (or below), lower bound
+ MaskImagePointType p;
+ try {
+ MaskImagePointer CricoidCartilag = this->GetAFDB()->template GetImage <MaskImageType>("CricoidCartilag");
+ p[0] = p[1] = p[2] = 0.0; // to avoid warning
+ clitk::FindExtremaPointInAGivenDirection<MaskImageType>(CricoidCartilag,
+ this->GetBackgroundValue(), 2, true, p);
+ } catch (clitk::ExceptionObject e) {
+ //DD("CricoidCartilag image not found, try CricoidCartilagZ");
+ this->GetAFDB()->GetPoint3D("CricoidCartilagPoint", p);
+ }
+ output = clitk::CropImageRemoveGreaterThan<MaskImageType>(output, 2, p[2], true, this->GetBackgroundValue());
+ this->template StopCurrentStep<MaskImageType>(output);
+ */
+
+ //--------------------------------------------------------------------
+ // Step: Get CCL
+ this->StartNewStep("[Mediastinum] Keep main connected component");
+ output = clitk::Labelize<MaskImageType>(output, this->GetBackgroundValue(), false, 500);
+ // output = RemoveLabels<MaskImageType>(output, BG, param->GetLabelsToRemove());
+ output = clitk::KeepLabels<MaskImageType>(output, this->GetBackgroundValue(),
+ this->GetForegroundValue(), 1, 1, 0);
+ this->template StopCurrentStep<MaskImageType>(output);
+
+ //--------------------------------------------------------------------
+ // Step: Remove post part from VertebralBody
+ this->StartNewStep("[Mediastinum] Remove post part according to VertebralBody");
+ RemovePostPartOfVertebralBody();
+ this->template StopCurrentStep<MaskImageType>(output);
+
+ //--------------------------------------------------------------------
+ // Step: Slice by Slice CCL
+ this->StartNewStep("[Mediastinum] Slice by Slice keep only one component");
+ typedef clitk::ExtractSliceFilter<MaskImageType> ExtractSliceFilterType;
+ // typename ExtractSliceFilterType::Pointer
+ ExtractSliceFilterType::Pointer extractSliceFilter = ExtractSliceFilterType::New();
+ extractSliceFilter->SetInput(output);
+ extractSliceFilter->SetDirection(2);
+ extractSliceFilter->Update();
+ typedef typename ExtractSliceFilterType::SliceType SliceType;
+ std::vector<typename SliceType::Pointer> mSlices;
+ extractSliceFilter->GetOutputSlices(mSlices);
+ for(unsigned int i=0; i<mSlices.size(); i++) {
+ mSlices[i] = Labelize<SliceType>(mSlices[i], 0, true, 100);
+ mSlices[i] = KeepLabels<SliceType>(mSlices[i], 0, 1, 1, 1, true);
+ }
+ typedef itk::JoinSeriesImageFilter<SliceType, MaskImageType> JoinSeriesFilterType;
+ typename JoinSeriesFilterType::Pointer joinFilter = JoinSeriesFilterType::New();
+ joinFilter->SetOrigin(output->GetOrigin()[2]);
+ joinFilter->SetSpacing(output->GetSpacing()[2]);
+ for(unsigned int i=0; i<mSlices.size(); i++) {
+ joinFilter->PushBackInput(mSlices[i]);
+ }
+ joinFilter->Update();
+ output = joinFilter->GetOutput();
+ this->template StopCurrentStep<MaskImageType>(output);
+
+ //--------------------------------------------------------------------
+ // Step 10 : AutoCrop
+ this->StartNewStep("[Mediastinum] AutoCrop");
+ output = clitk::AutoCrop<MaskImageType>(output, this->GetBackgroundValue());
+ this->template StopCurrentStep<MaskImageType>(output);
+
+ // End, set the real size
+ this->GetOutput(0)->SetRegions(output->GetLargestPossibleRegion());
+ this->GetOutput(0)->SetLargestPossibleRegion(output->GetLargestPossibleRegion());
+ this->GetOutput(0)->SetRequestedRegion(output->GetLargestPossibleRegion());
+ this->GetOutput(0)->SetBufferedRegion(output->GetLargestPossibleRegion());
}
//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template <class ImageType>
+void
+clitk::ExtractMediastinumFilter<ImageType>::
+GenerateData()
+{
+ this->GraftOutput(output);
+ // Store image filenames into AFDB
+ this->GetAFDB()->SetImageFilename("Mediastinum", this->GetOutputMediastinumFilename());
+ this->WriteAFDB();
+}
+//--------------------------------------------------------------------
+
+//--------------------------------------------------------------------
+template <class ImageType>
+void
+clitk::ExtractMediastinumFilter<ImageType>::
+RemovePostPartOfVertebralBody()
+{
+
+ /*
+ Posteriorly, Station 8 abuts the descending aorta and the anterior
+ aspect of the vertebral body until an imaginary horizontal line
+ running 1 cm posterior to the anterior border of the vertebral
+ body (Fig. 3C).
+
+ => We use this definition for all the mediastinum
+
+ Find most Ant point in VertebralBody. Consider the horizontal line
+ which is 'DistanceMaxToAnteriorPartOfTheVertebralBody' away from
+ the most ant point.
+ */
+
+ // Get VertebralBody mask image
+ MaskImagePointer VertebralBody =
+ this->GetAFDB()->template GetImage <MaskImageType>("VertebralBody");
+
+ // Consider vertebral body slice by slice
+ std::vector<MaskSlicePointer> vertebralSlices;
+ clitk::ExtractSlices<MaskImageType>(VertebralBody, 2, vertebralSlices);
+
+ // For each slice, compute the most anterior point
+ std::map<int, MaskSlicePointType> vertebralAntPositionBySlice;
+ for(uint i=0; i<vertebralSlices.size(); i++) {
+ MaskSlicePointType p;
+ bool found = clitk::FindExtremaPointInAGivenDirection<MaskSliceType>(vertebralSlices[i],
+ this->GetBackgroundValue(),
+ 1, true, p);
+ if (found) {
+ vertebralAntPositionBySlice[i] = p;
+ }
+ else {
+ // It should not happen ! But sometimes, a contour is missing or
+ // the VertebralBody is not delineated enough inferiorly ... in
+ // those cases, we consider the first found slice.
+ // std::cerr << "No foreground pixels in this VertebralBody slices !?? I try with the previous/next slice" << std::endl;
+ // [ Possible alternative -> consider previous limit ]
+ }
+ }
+
+ // Convert 2D points in slice into 3D points
+ std::vector<MaskImagePointType> vertebralAntPositions;
+ clitk::PointsUtils<MaskImageType>::Convert2DMapTo3DList(vertebralAntPositionBySlice,
+ VertebralBody,
+ vertebralAntPositions);
+
+ // DEBUG : write list of points
+ clitk::WriteListOfLandmarks<MaskImageType>(vertebralAntPositions,
+ "VertebralBodyMostAnteriorPoints.txt");
+
+ // Cut support posteriorly 1cm the most anterior point of the
+ // VertebralBody. Do nothing below and above the VertebralBody. To
+ // do that compute several region, slice by slice and fill.
+ MaskImageRegionType region;
+ MaskImageSizeType size;
+ MaskImageIndexType start;
+ size[2] = 1; // a single slice
+ start[0] = output->GetLargestPossibleRegion().GetIndex()[0];
+ size[0] = output->GetLargestPossibleRegion().GetSize()[0];
+ for(uint i=0; i<vertebralAntPositions.size(); i++) {
+ typedef typename itk::ImageRegionIterator<MaskImageType> IteratorType;
+ IteratorType iter =
+ IteratorType(output, output->GetLargestPossibleRegion());
+ MaskImageIndexType index;
+ // Consider some cm posterior to most anterior positions (usually
+ // 1 cm).
+ vertebralAntPositions[i][1] += GetDistanceMaxToAnteriorPartOfTheVertebralBody();
+ // Get index of this point
+ output->TransformPhysicalPointToIndex(vertebralAntPositions[i], index);
+ // Compute region (a single slice)
+ start[2] = index[2];
+ start[1] = output->GetLargestPossibleRegion().GetIndex()[1]+index[1];
+ size[1] = output->GetLargestPossibleRegion().GetSize()[1]-start[1];
+ region.SetSize(size);
+ region.SetIndex(start);
+ // Fill region
+ if (output->GetLargestPossibleRegion().IsInside(start)) {
+ itk::ImageRegionIterator<MaskImageType> it(output, region);
+ it.GoToBegin();
+ while (!it.IsAtEnd()) {
+ it.Set(this->GetBackgroundValue());
+ ++it;
+ }
+ }
+ }
+}
+//--------------------------------------------------------------------
+
#endif //#define CLITKBOOLEANOPERATORLABELIMAGEFILTER_TXX
git://git.creatis.insa-lyon.fr / clitk.git / blobdiff