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
template <class ImageType>
clitk::ExtractMediastinumFilter<ImageType>::
ExtractMediastinumFilter():
- clitk::FilterBase(),
- clitk::FilterWithAnatomicalFeatureDatabaseManagement(),
- itk::ImageToImageFilter<ImageType, MaskImageType>()
+ clitk::StructuresExtractionFilter<ImageType>()
{
this->SetNumberOfRequiredInputs(1);
-
SetBackgroundValuePatient(0);
SetBackgroundValueLung(0);
SetBackgroundValueBones(0);
SetForegroundValueLeftLung(1);
SetForegroundValueRightLung(2);
- SetBackgroundValue(0);
- SetForegroundValue(1);
-
- SetIntermediateSpacing(6);
- SetFuzzyThreshold1(0.5);
- SetFuzzyThreshold2(0.6);
- SetFuzzyThreshold3(0.05);
-
- SetDistanceMaxToAnteriorPartOfTheSpine(10);
- SetOutputMediastinumFilename("mediastinum.mhd");
-
+ SetDistanceMaxToAnteriorPartOfTheVertebralBody(10);
+ SetOutputMediastinumFilename("mediastinum.mha");
UseBonesOff();
}
//--------------------------------------------------------------------
}
//--------------------------------------------------------------------
-
//--------------------------------------------------------------------
template <class ImageType>
-template<class ArgsInfoType>
void
clitk::ExtractMediastinumFilter<ImageType>::
-SetArgsInfo(ArgsInfoType mArgsInfo)
+SetFuzzyThreshold(std::string tag, double value)
{
- SetVerboseOption_GGO(mArgsInfo);
- SetVerboseStep_GGO(mArgsInfo);
- SetWriteStep_GGO(mArgsInfo);
- SetVerboseWarningOff_GGO(mArgsInfo);
-
- SetIntermediateSpacing_GGO(mArgsInfo);
- SetFuzzyThreshold1_GGO(mArgsInfo);
- SetFuzzyThreshold2_GGO(mArgsInfo);
- SetFuzzyThreshold3_GGO(mArgsInfo);
-
- SetAFDBFilename_GGO(mArgsInfo);
- SetDistanceMaxToAnteriorPartOfTheSpine_GGO(mArgsInfo);
- SetUseBones_GGO(mArgsInfo);
+ m_FuzzyThreshold[tag] = value;
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+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;
+ }
- SetLowerThreshold_GGO(mArgsInfo);
- SetUpperThreshold_GGO(mArgsInfo);
+ return m_FuzzyThreshold[tag];
}
//--------------------------------------------------------------------
clitk::ExtractMediastinumFilter<ImageType>::
GenerateInputRequestedRegion()
{
- //DD("GenerateInputRequestedRegion");
+ // DD("GenerateInputRequestedRegion");
// Do not call default
- // Superclass::GenerateInputRequestedRegion();
- // DD("End GenerateInputRequestedRegion");
+ // Superclass::GenerateInputRequestedRegion();
+ // DD("End GenerateInputRequestedRegion");
}
-
//--------------------------------------------------------------------
//--------------------------------------------------------------------
+
//--------------------------------------------------------------------
template <class ImageType>
void
clitk::ExtractMediastinumFilter<ImageType>::
GenerateOutputInformation() {
- // DD("GenerateOutputInformation");
// Do not call default
- // Superclass::GenerateOutputInformation();
+ // Superclass::GenerateOutputInformation();
//--------------------------------------------------------------------
// Get input pointers
- LoadAFDB();
+ this->LoadAFDB();
ImageConstPointer input = dynamic_cast<const ImageType*>(itk::ProcessObject::GetInput(0));
- MaskImagePointer patient = GetAFDB()->template GetImage <MaskImageType>("patient");
- MaskImagePointer lung = GetAFDB()->template GetImage <MaskImageType>("lungs");
- MaskImagePointer bones = GetAFDB()->template GetImage <MaskImageType>("bones");
- MaskImagePointer trachea = GetAFDB()->template GetImage <MaskImageType>("trachea");
-
+ 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 1: Crop support (patient) to lung extend in RL
- StartNewStep("Crop support like lungs along LR");
+ // 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);
this->template StopCurrentStep<MaskImageType>(output);
//--------------------------------------------------------------------
- // Step 2: Crop support (previous) to bones extend in AP
+ // Step : Crop support (previous) to bones extend in AP
if (GetUseBones()) {
- StartNewStep("Crop support like bones along AP");
+ 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
}
//--------------------------------------------------------------------
- // Step 3: patient minus lungs, minus bones, minus trachea
- StartNewStep("Patient contours minus lungs, bones, trachea");
+ // 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();
output = boolFilter->GetOutput();
// Auto crop to gain support area
- output = clitk::AutoCrop<MaskImageType>(output, GetBackgroundValue());
+ output = clitk::AutoCrop<MaskImageType>(output, this->GetBackgroundValue());
this->template StopCurrentStep<MaskImageType>(output);
//--------------------------------------------------------------------
- // Step 4: LR limits from lung (need separate lung ?)
+ // 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 ? )
- StartNewStep("Left/Right limits with lungs");
- /*
- ImagePointer right_lung = clitk::SetBackground<MaskImageType, MaskImageType>(lung, lung, 2, 0);
- ImagePointer left_lung = clitk::SetBackground<MaskImageType, MaskImageType>(lung, lung, 1, 0);
- writeImage<MaskImageType>(right_lung, "right.mhd");
- writeImage<MaskImageType>(left_lung, "left.mhd");
- */
- typedef clitk::AddRelativePositionConstraintToLabelImageFilter<MaskImageType> RelPosFilterType;
- typename RelPosFilterType::Pointer relPosFilter = RelPosFilterType::New();
- relPosFilter->SetCurrentStepBaseId(this->GetCurrentStepId());
- relPosFilter->VerboseStepOff();
- relPosFilter->WriteStepOff();
- relPosFilter->SetInput(output);
- //relPosFilter->SetInputObject(left_lung);
- relPosFilter->SetInputObject(lung);
- relPosFilter->SetOrientationType(RelPosFilterType::LeftTo); // warning left lung is at right ;)
- relPosFilter->SetIntermediateSpacing(GetIntermediateSpacing());
- relPosFilter->SetFuzzyThreshold(GetFuzzyThreshold1());
- relPosFilter->Update();
- output = relPosFilter->GetOutput();
- //writeImage<MaskImageType>(right_lung, "step4-left.mhd");
-
- relPosFilter->SetInput(output);
- relPosFilter->SetCurrentStepBaseId(this->GetCurrentStepId());
- relPosFilter->VerboseStepOff();
- relPosFilter->WriteStepOff();
- //relPosFilter->SetInputObject(right_lung);
- relPosFilter->SetInputObject(lung);
- relPosFilter->SetOrientationType(RelPosFilterType::RightTo);
- relPosFilter->SetIntermediateSpacing(GetIntermediateSpacing());
- relPosFilter->SetFuzzyThreshold(GetFuzzyThreshold1());
- relPosFilter->Update();
- output = relPosFilter->GetOutput();
+ 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 5: AP limits from bones
+ // Step : AP limits from bones
// Separate the bones in the ant-post middle
- MaskImageConstPointer bones_ant;
- MaskImageConstPointer bones_post;
- MaskImagePointType middle_AntPost__position;
+ MaskImagePointer bones_ant;
+ MaskImagePointer bones_post;
+ MaskImagePointType middle_AntPost_position;
+ middle_AntPost_position.Fill(NumericTraits<MaskImagePointType::ValueType>::Zero);
if (GetUseBones()) {
- StartNewStep("Ant/Post limits with bones");
- 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;
+ 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;
- bones->TransformPhysicalPointToIndex(middle_AntPost__position, 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;
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;
+ // 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);
roiFilter->ReleaseDataFlagOff();
roiFilter->Update();
bones_ant = roiFilter->GetOutput();
- writeImage<MaskImageType>(bones_ant, "b_ant.mhd");
+ // writeImage<MaskImageType>(bones_ant, "b_ant.mhd");
// POST part
roiFilter = ROIFilterType::New();
index[1] = bones->GetLargestPossibleRegion().GetIndex()[1] + size[1]-1;
roiFilter->ReleaseDataFlagOff();
roiFilter->Update();
bones_post = roiFilter->GetOutput();
- writeImage<MaskImageType>(bones_post, "b_post.mhd");
-
- // Go !
- relPosFilter->SetCurrentStepNumber(0);
- relPosFilter->ResetPipeline();// = RelPosFilterType::New();
- relPosFilter->SetCurrentStepBaseId(this->GetCurrentStepId());
- relPosFilter->VerboseStepOff();
- relPosFilter->WriteStepOff();
- relPosFilter->SetInput(output);
- relPosFilter->SetInputObject(bones_post);
- relPosFilter->SetOrientationType(RelPosFilterType::AntTo);
- relPosFilter->SetIntermediateSpacing(GetIntermediateSpacing());
- relPosFilter->SetFuzzyThreshold(GetFuzzyThreshold2());
- relPosFilter->Update();
- output = relPosFilter->GetOutput();
- writeImage<MaskImageType>(output, "post.mhd");
-
- relPosFilter->SetInput(relPosFilter->GetOutput());
- relPosFilter->SetCurrentStepBaseId(this->GetCurrentStepId());
- relPosFilter->VerboseStepOff();
- relPosFilter->WriteStepOff();
- relPosFilter->SetInput(output);
- relPosFilter->SetInputObject(bones_ant);
- relPosFilter->SetOrientationType(RelPosFilterType::PostTo);
- relPosFilter->SetIntermediateSpacing(GetIntermediateSpacing());
- relPosFilter->SetFuzzyThreshold(GetFuzzyThreshold2());
- relPosFilter->Update();
- output = relPosFilter->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);
}
//--------------------------------------------------------------------
- // Step 6: Get CCL
- StartNewStep("Keep main connected component");
- output = clitk::Labelize<MaskImageType>(output, GetBackgroundValue(), false, 500);
- // output = RemoveLabels<MaskImageType>(output, BG, param->GetLabelsToRemove());
- output = clitk::KeepLabels<MaskImageType>(output, GetBackgroundValue(),
- GetForegroundValue(), 1, 1, 0);
+ // 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);
//--------------------------------------------------------------------
- // Step 7 : Slice by Slice to optimize posterior part
- // Warning slice does not necesseraly correspond between 'output' and 'bones'
- typedef clitk::ExtractSliceFilter<MaskImageType> ExtractSliceFilterType;
- typename ExtractSliceFilterType::Pointer extractSliceFilter = ExtractSliceFilterType::New();
- typedef typename ExtractSliceFilterType::SliceType SliceType;
- std::vector<typename SliceType::Pointer> mSlices;
- if (GetUseBones()) {
- StartNewStep("Rafine posterior part according to vertebral body");
- extractSliceFilter->SetInput(bones_post);
- extractSliceFilter->SetDirection(2);
- extractSliceFilter->Update();
- std::vector<double> mVertebralAntPositionBySlice;
- extractSliceFilter->GetOutputSlices(mSlices);
- for(unsigned int i=0; i<mSlices.size(); i++) {
- mSlices[i] = Labelize<SliceType>(mSlices[i], 0, true, 10);
- mSlices[i] = KeepLabels<SliceType>(mSlices[i],
- GetBackgroundValue(),
- GetForegroundValue(), 1, 2, true); // keep two first
- // Find most anterior point (start of the vertebral)
- typename itk::ImageRegionIteratorWithIndex<SliceType>
- iter(mSlices[i], mSlices[i]->GetLargestPossibleRegion());
- iter.GoToBegin();
- bool stop = false;
- while (!stop) {
- if (iter.Get() != GetBackgroundValue())
- stop = true; // not foreground because we keep two main label
- ++iter;
- if (iter.IsAtEnd()) stop = true;
- }
- if (!iter.IsAtEnd()) {
- typename SliceType::PointType p;
- mSlices[i]->TransformIndexToPhysicalPoint(iter.GetIndex(),p);
- mVertebralAntPositionBySlice.push_back(p[1]);
- }
- else {
- mVertebralAntPositionBySlice.push_back(bones_post->GetOrigin()[1]+(bones->GetLargestPossibleRegion().GetSize()[1]*bones->GetSpacing()[1]));
- DD(mVertebralAntPositionBySlice.back());
- DD("ERROR ?? NO FG in bones here ?");
- }
- }
+ // Generic RelativePosition processes
+ output = this->ApplyRelativePositionList("Mediastinum", output);
- // Cut Post position slice by slice
- {
- MaskImageRegionType region;
- MaskImageSizeType size;
- MaskImageIndexType start;
- size[2] = 1;
- start[0] = output->GetLargestPossibleRegion().GetIndex()[0];
- for(unsigned int i=0; i<mSlices.size(); i++) {
- // Compute index
- MaskImagePointType point;
- point[0] = 0;
-
- //TODO 10 mm OPTION
-
- point[1] = mVertebralAntPositionBySlice[i]+GetDistanceMaxToAnteriorPartOfTheSpine();// ADD ONE CM
- point[2] = bones_post->GetOrigin()[2]+(bones_post->GetLargestPossibleRegion().GetIndex()[2]+i)*bones_post->GetSpacing()[2];
- MaskImageIndexType index;
- output->TransformPhysicalPointToIndex(point, index);
- // Compute region
- start[2] = index[2];
- start[1] = output->GetLargestPossibleRegion().GetIndex()[1]+index[1];
- size[0] = output->GetLargestPossibleRegion().GetSize()[0];
- size[1] = output->GetLargestPossibleRegion().GetSize()[1]-start[1];
- region.SetSize(size);
- region.SetIndex(start);
- // Fill Region
- if (output->GetLargestPossibleRegion().IsInside(start)) {
- itk::ImageRegionIteratorWithIndex<MaskImageType> it(output, region);
- it.GoToBegin();
- while (!it.IsAtEnd()) {
- it.Set(GetBackgroundValue());
- ++it;
- }
- }
- }
- }
- this->template StopCurrentStep<MaskImageType>(output);
- }
//--------------------------------------------------------------------
- // Step 8: Trial segmentation KMeans
- if (0) {
- StartNewStep("K means");
- // Take input, crop like current mask
- typedef CropLikeImageFilter<ImageType> CropLikeFilterType;
- typename CropLikeFilterType::Pointer cropLikeFilter = CropLikeFilterType::New();
- cropLikeFilter->SetInput(input);
- cropLikeFilter->SetCropLikeImage(output);
- cropLikeFilter->Update();
- ImagePointer working_input = cropLikeFilter->GetOutput();
- writeImage<ImageType>(working_input, "crop-input.mhd");
- // Set bG at -1000
- working_input = clitk::SetBackground<ImageType, MaskImageType>(working_input, output, GetBackgroundValue(), -1000);
- writeImage<ImageType>(working_input, "crop-input2.mhd");
- // Kmeans
- typedef itk::ScalarImageKmeansImageFilter<ImageType> KMeansFilterType;
- typename KMeansFilterType::Pointer kmeansFilter = KMeansFilterType::New();
- kmeansFilter->SetInput(working_input);
- // const unsigned int numberOfInitialClasses = 3;
- // const unsigned int useNonContiguousLabels = 0;
- kmeansFilter->AddClassWithInitialMean(-1000);
- kmeansFilter->AddClassWithInitialMean(30);
- kmeansFilter->AddClassWithInitialMean(-40); // ==> I want this one
- DD("Go!");
- kmeansFilter->Update();
- DD("End");
- typename KMeansFilterType::ParametersType estimatedMeans = kmeansFilter->GetFinalMeans();
- const unsigned int numberOfClasses = estimatedMeans.Size();
- for ( unsigned int i = 0 ; i < numberOfClasses ; ++i ) {
- std::cout << "cluster[" << i << "] ";
- std::cout << " estimated mean : " << estimatedMeans[i] << std::endl;
- }
- MaskImageType::Pointer kmeans = kmeansFilter->GetOutput();
- kmeans = clitk::SetBackground<MaskImageType, MaskImageType>(kmeans, kmeans,
- 1, GetBackgroundValue());
- writeImage<MaskImageType>(kmeans, "kmeans.mhd");
- // Get final results, and remove from current mask
- boolFilter = BoolFilterType::New();
- boolFilter->InPlaceOn();
- boolFilter->SetInput1(output);
- boolFilter->SetInput2(kmeans);
- boolFilter->SetOperationType(BoolFilterType::And);
- boolFilter->Update();
- output = boolFilter->GetOutput();
- writeImage<MaskImageType>(output, "out-kmean.mhd");
- this->template StopCurrentStep<MaskImageType>(output);
-
- // TODO -> FillMASK ?
- // comment speed ? mask ? 2 class ?
-
-
- //TODO
- // Confidence connected ?
-
+ // 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 8: Lower limits from lung (need separate lung ?)
- if (1) {
- // StartNewStep("Trial : minus segmented struct");
- // MaskImagePointer heart = GetAFDB()->template GetImage <MaskImageType>("heart");
- // boolFilter = BoolFilterType::New();
- // boolFilter->InPlaceOn();
- // boolFilter->SetInput1(output);
- // boolFilter->SetInput2(heart);
- // boolFilter->SetOperationType(BoolFilterType::AndNot);
- // boolFilter->Update();
- // output = boolFilter->GetOutput(); // not needed because InPlace
-
- // Not below the heart
- // relPosFilter = RelPosFilterType::New();
- // relPosFilter->SetCurrentStepBaseId(this->GetCurrentStepId());
- // relPosFilter->VerboseStepOff();
- // relPosFilter->WriteStepOff();
- // relPosFilter->SetInput(output);
- // relPosFilter->SetInputObject(heart);
- // relPosFilter->SetOrientationType(RelPosFilterType::SupTo);
- // relPosFilter->SetIntermediateSpacing(GetIntermediateSpacing());
- // relPosFilter->SetFuzzyThreshold(GetFuzzyThreshold3());
- // relPosFilter->Update();
- // output = relPosFilter->GetOutput();
- }
+ // 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 8: Lower limits from lung (need separate lung ?)
- if (0) {
- StartNewStep("Lower limits with lungs");
- // TODO BOFFF ????
- relPosFilter = RelPosFilterType::New();
- relPosFilter->SetCurrentStepBaseId(this->GetCurrentStepId());
- relPosFilter->VerboseStepOff();
- relPosFilter->WriteStepOff();
- relPosFilter->SetInput(output);
- // relPosFilter->SetInputObject(left_lung);
- relPosFilter->SetInputObject(lung);
- relPosFilter->SetOrientationType(RelPosFilterType::SupTo);
- relPosFilter->SetIntermediateSpacing(GetIntermediateSpacing());
- relPosFilter->SetFuzzyThreshold(GetFuzzyThreshold3());
- relPosFilter->Update();
- output = relPosFilter->GetOutput();
- 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 10: Slice by Slice CCL
- StartNewStep("Slice by Slice keep only one component");
+ // Step: Slice by Slice CCL
+ this->StartNewStep("[Mediastinum] Slice by Slice keep only one component");
typedef clitk::ExtractSliceFilter<MaskImageType> ExtractSliceFilterType;
// typename ExtractSliceFilterType::Pointer
- extractSliceFilter = ExtractSliceFilterType::New();
+ ExtractSliceFilterType::Pointer extractSliceFilter = ExtractSliceFilterType::New();
extractSliceFilter->SetInput(output);
extractSliceFilter->SetDirection(2);
extractSliceFilter->Update();
typedef typename ExtractSliceFilterType::SliceType SliceType;
- // std::vector<typename SliceType::Pointer>
- mSlices.clear();
+ 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);
output = joinFilter->GetOutput();
this->template StopCurrentStep<MaskImageType>(output);
- //--------------------------------------------------------------------
- // Step 9: Binarize to remove too high HU
- // --> warning CCL slice by slice must be done before
- StartNewStep("Remove hypersignal (bones and injected part");
- // Crop initial ct like current support
- typedef CropLikeImageFilter<ImageType> CropLikeFilterType;
- typename CropLikeFilterType::Pointer cropLikeFilter = CropLikeFilterType::New();
- cropLikeFilter->SetInput(input);
- cropLikeFilter->SetCropLikeImage(output);
- cropLikeFilter->Update();
- ImagePointer working_input = cropLikeFilter->GetOutput();
- writeImage<ImageType>(working_input, "crop-ct.mhd");
- // Binarize
- typedef itk::BinaryThresholdImageFilter<ImageType, MaskImageType> InputBinarizeFilterType;
- typename InputBinarizeFilterType::Pointer binarizeFilter=InputBinarizeFilterType::New();
- binarizeFilter->SetInput(working_input);
- binarizeFilter->SetLowerThreshold(GetLowerThreshold());
- binarizeFilter->SetUpperThreshold(GetUpperThreshold());
- binarizeFilter->SetInsideValue(this->GetBackgroundValue()); // opposite
- binarizeFilter->SetOutsideValue(this->GetForegroundValue()); // opposite
- binarizeFilter->Update();
- MaskImagePointer working_bin = binarizeFilter->GetOutput();
- writeImage<MaskImageType>(working_bin, "bin.mhd");
- // Remove from support
- boolFilter = BoolFilterType::New();
- boolFilter->InPlaceOn();
- boolFilter->SetInput1(output);
- boolFilter->SetInput2(working_bin);
- boolFilter->SetOperationType(BoolFilterType::AndNot);
- boolFilter->Update();
- output = boolFilter->GetOutput();
- StopCurrentStep<MaskImageType>(output);
-
//--------------------------------------------------------------------
// Step 10 : AutoCrop
- StartNewStep("AutoCrop");
- output = clitk::AutoCrop<MaskImageType>(output, GetBackgroundValue());
+ this->StartNewStep("[Mediastinum] AutoCrop");
+ output = clitk::AutoCrop<MaskImageType>(output, this->GetBackgroundValue());
this->template StopCurrentStep<MaskImageType>(output);
- // Bones ? pb with RAM ? FillHoles ?
-
- // how to do with post part ? spine /lung ?
- // POST the spine (should be separated from the rest)
- /// DO THAT ---->>
- // histo Y on the whole bones_post (3D) -> result is the Y center on the spine (?)
- // by slice on bones_post
- // find the most ant point in the center
- // from this point go to post until out of bones.
- //
-
-
// End, set the real size
this->GetOutput(0)->SetRegions(output->GetLargestPossibleRegion());
this->GetOutput(0)->SetLargestPossibleRegion(output->GetLargestPossibleRegion());
clitk::ExtractMediastinumFilter<ImageType>::
GenerateData()
{
- DD("GenerateData");
this->GraftOutput(output);
// Store image filenames into AFDB
- GetAFDB()->SetImageFilename("mediastinum", this->GetOutputMediastinumFilename());
- WriteAFDB();
+ 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