From 5e2af376544fce0c6dc46bb3c3227d35b501c1f1 Mon Sep 17 00:00:00 2001 From: dsarrut Date: Wed, 15 Dec 2010 08:54:58 +0000 Subject: [PATCH] various update --- .../clitkAnatomicalFeatureDatabase.cxx | 13 +- segmentation/clitkAnatomicalFeatureDatabase.h | 1 + .../clitkExtractAirwaysTreeInfoFilter.h | 24 +- .../clitkExtractAirwaysTreeInfoFilter.txx | 280 +++++++- segmentation/clitkExtractBones.ggo | 4 + segmentation/clitkExtractBonesFilter.h | 12 +- segmentation/clitkExtractBonesFilter.txx | 20 +- segmentation/clitkExtractLungFilter.txx | 18 +- segmentation/clitkExtractLymphStations.ggo | 11 +- .../clitkExtractLymphStationsFilter.h | 121 ++-- .../clitkExtractLymphStationsFilter.txx | 266 ++----- .../clitkExtractLymphStationsGenericFilter.h | 12 +- ...clitkExtractLymphStationsGenericFilter.txx | 53 +- segmentation/clitkExtractMediastinum.ggo | 23 +- segmentation/clitkExtractMediastinumFilter.h | 140 ++-- .../clitkExtractMediastinumFilter.txx | 647 +++++++++++++----- .../clitkExtractMediastinumGenericFilter.txx | 33 +- 17 files changed, 1068 insertions(+), 610 deletions(-) diff --git a/segmentation/clitkAnatomicalFeatureDatabase.cxx b/segmentation/clitkAnatomicalFeatureDatabase.cxx index 007aee5..539e170 100644 --- a/segmentation/clitkAnatomicalFeatureDatabase.cxx +++ b/segmentation/clitkAnatomicalFeatureDatabase.cxx @@ -22,8 +22,8 @@ // std #include #include -#include -#include +#include +#include //-------------------------------------------------------------------- clitk::AnatomicalFeatureDatabase::AnatomicalFeatureDatabase() @@ -97,6 +97,15 @@ void clitk::AnatomicalFeatureDatabase::SetPoint3D(std::string tag, PointType3D & //-------------------------------------------------------------------- +//-------------------------------------------------------------------- +double clitk::AnatomicalFeatureDatabase::GetPoint3D(std::string tag, int dim) +{ + PointType3D p; + GetPoint3D(tag, p); + return p[dim]; +} +//-------------------------------------------------------------------- + //-------------------------------------------------------------------- void clitk::AnatomicalFeatureDatabase::GetPoint3D(std::string tag, PointType3D & p) { diff --git a/segmentation/clitkAnatomicalFeatureDatabase.h b/segmentation/clitkAnatomicalFeatureDatabase.h index dcf7efd..f1324cb 100644 --- a/segmentation/clitkAnatomicalFeatureDatabase.h +++ b/segmentation/clitkAnatomicalFeatureDatabase.h @@ -50,6 +50,7 @@ namespace clitk { typedef itk::Point PointType3D; void SetPoint3D(TagType tag, PointType3D & p); void GetPoint3D(TagType tag, PointType3D & p); + double GetPoint3D(std::string tag, int dim); // Set Get image void SetImageFilename(TagType tag, std::string f); diff --git a/segmentation/clitkExtractAirwaysTreeInfoFilter.h b/segmentation/clitkExtractAirwaysTreeInfoFilter.h index 220d8ce..2b4d8fc 100644 --- a/segmentation/clitkExtractAirwaysTreeInfoFilter.h +++ b/segmentation/clitkExtractAirwaysTreeInfoFilter.h @@ -53,8 +53,10 @@ namespace clitk { _l1 = l1; _l2 = l2; } + IndexType index; PointType point; + double weight; PixelType l; PixelType l1; PixelType l2; @@ -109,7 +111,23 @@ namespace clitk { typedef typename InternalImageType::Pointer InternalImagePointer; typedef typename InternalImageType::IndexType InternalIndexType; typedef LabelizeParameters LabelParamType; - + + // Data Structures for trees + struct FullTreeNodeType { + ImageIndexType index; + ImagePointType point; + double weight; + }; + struct StructuralTreeNodeType { + ImageIndexType index; + ImagePointType point; + double weight; + }; + typedef tree FullTreeType; + typedef tree StructuralTreeType; + FullTreeType mFullSkeletonTree; + StructuralTreeType mStructuralSkeletonTree; + /** Connect inputs */ void SetInput(const ImageType * image); @@ -148,6 +166,10 @@ namespace clitk { virtual void GenerateData(); TreeType m_SkeletonTree; + void TrackFromThisIndex2(ImageIndexType index, ImagePointer skeleton, + ImagePixelType label, + typename FullTreeType::iterator currentNode, + typename StructuralTreeType::iterator currentSNode); void TrackFromThisIndex(std::vector & listOfBifurcations, ImagePointer skeleton, ImageIndexType index, diff --git a/segmentation/clitkExtractAirwaysTreeInfoFilter.txx b/segmentation/clitkExtractAirwaysTreeInfoFilter.txx index 6145ad7..da32e04 100644 --- a/segmentation/clitkExtractAirwaysTreeInfoFilter.txx +++ b/segmentation/clitkExtractAirwaysTreeInfoFilter.txx @@ -108,14 +108,30 @@ void clitk::ExtractAirwaysTreeInfoFilter:: GenerateData() { - StartNewStep("Thinning filter (skeleton)"); + + // Crop just abov lungs ? + + // read skeleton if already exist + bool foundSkeleton = true; + try { + skeleton = GetAFDB()->template GetImage ("skeleton"); + } + catch (clitk::ExceptionObject e) { + DD("did not found skeleton"); + foundSkeleton = false; + } + // Extract skeleton - typedef itk::BinaryThinningImageFilter3D ThinningFilterType; - typename ThinningFilterType::Pointer thinningFilter = ThinningFilterType::New(); - thinningFilter->SetInput(input); // input = trachea - thinningFilter->Update(); - skeleton = thinningFilter->GetOutput(); - StopCurrentStep(skeleton); + if (!foundSkeleton) { + StartNewStep("Thinning filter (skeleton)"); + typedef itk::BinaryThinningImageFilter3D ThinningFilterType; + typename ThinningFilterType::Pointer thinningFilter = ThinningFilterType::New(); + thinningFilter->SetInput(input); // input = trachea + thinningFilter->Update(); + skeleton = thinningFilter->GetOutput(); + StopCurrentStep(skeleton); + writeImage(skeleton, "skeleton.mhd"); + } // Find first point for tracking StartNewStep("Find first point for tracking"); @@ -131,16 +147,18 @@ GenerateData() typename ImageType::IndexType index = it.GetIndex(); DD(index); - // Initialize neighborhooditerator - typedef itk::NeighborhoodIterator NeighborhoodIteratorType; - typename NeighborhoodIteratorType::SizeType radius; - radius.Fill(1); - NeighborhoodIteratorType nit(radius, skeleton, skeleton->GetLargestPossibleRegion()); + if (0) { + // Initialize neighborhooditerator + typedef itk::NeighborhoodIterator NeighborhoodIteratorType; + typename NeighborhoodIteratorType::SizeType radius; + radius.Fill(1); + NeighborhoodIteratorType nit(radius, skeleton, skeleton->GetLargestPossibleRegion()); + } - // Find first label number (must be different from BG and FG) + // Find a label number different from BG and FG typename ImageType::PixelType label = GetForegroundValue()+1; while ((label == GetBackgroundValue()) || (label == GetForegroundValue())) { label++; } - DD(label); + DD((int)label); /* Tracking ? @@ -151,20 +169,138 @@ GenerateData() -> follow at Left */ - - // Track from the first point + // NEW Track from the first point StartNewStep("Start tracking"); - std::vector listOfBifurcations; - m_SkeletonTree.set_head(index); - TrackFromThisIndex(listOfBifurcations, skeleton, index, label, m_SkeletonTree.begin()); - DD("end track"); - - // Convert index into physical point coordinates - for(unsigned int i=0; iTransformIndexToPhysicalPoint(listOfBifurcations[i].index, - listOfBifurcations[i].point); + FullTreeNodeType n; + n.index = index; + skeleton->TransformIndexToPhysicalPoint(n.index, n.point); // à mettre dans TrackFromThisIndex2 + mFullSkeletonTree.set_head(n); + StructuralTreeNodeType sn; + sn.index = index; + skeleton->TransformIndexToPhysicalPoint(sn.index, sn.point); + mStructuralSkeletonTree.set_head(sn); + TrackFromThisIndex2(index, skeleton, label, mFullSkeletonTree.begin(), mStructuralSkeletonTree.begin()); + StopCurrentStep(); + + // Reput FG instead of label in the skeleton image + skeleton = clitk::SetBackground(skeleton, skeleton, label, GetForegroundValue()); + + // Debug + typename StructuralTreeType::iterator sit = mStructuralSkeletonTree.begin(); + while (sit != mStructuralSkeletonTree.end()) { + DD(sit->point); + ++sit; + } + + // compute weight : n longueurs à chaque bifurfaction. + // parcours de FullTreeNodeType, à partir de leaf, remonter de proche en proche la distance eucl. + // par post order + + // Init + typename FullTreeType::iterator fit = mFullSkeletonTree.begin(); + while (fit != mFullSkeletonTree.end()) { + fit->weight = 0.0; + ++fit; + } + + DD("compute weight"); + typename FullTreeType::post_order_iterator pit = mFullSkeletonTree.begin_post(); + while (pit != mFullSkeletonTree.end_post()) { + //DD(pit->point); + /* + if (pit != mFullSkeletonTree.begin()) { + typename FullTreeType::iterator parent = mFullSkeletonTree.parent(pit); + double d = pit->point.EuclideanDistanceTo(parent->point); + // DD(parent->weight); + //DD(d); + parent->weight += d; + if (pit.number_of_children() > 1) { + DD(pit.number_of_children()); + DD(pit->point); + DD(pit->weight); + DD(parent->weight); + DD(mFullSkeletonTree.child(pit, 0)->weight); + DD(mFullSkeletonTree.child(pit, 1)->weight); + } + } + */ + double previous = pit->weight; + for(uint i=0; ipoint.EuclideanDistanceTo(mFullSkeletonTree.child(pit, i)->point); + // pit->weight = pit->weight + mFullSkeletonTree.child(pit, i)->weight + d; + if (i==0) + pit->weight = pit->weight + mFullSkeletonTree.child(pit, i)->weight + d; + if (i>0) { + DD(pit.number_of_children()); + DD(pit->point); + DD(pit->weight); + DD(mFullSkeletonTree.child(pit, 0)->weight); + DD(mFullSkeletonTree.child(pit, 1)->weight); + pit->weight = std::max(pit->weight, previous+mFullSkeletonTree.child(pit, i)->weight + d); + } + } + ++pit; + } + + DD("end"); + fit = mFullSkeletonTree.begin(); + while (fit != mFullSkeletonTree.end()) { + std::cout << "p = " << fit->point + << " " << fit->weight + << " child= " << fit.number_of_children() + << " -> "; + for(uint i=0; iweight; + } + std::cout << std::endl; + ++fit; + } + + /* Selection criteria ? + - at least 2 child + - from top to bottom + - equilibr ? + */ + DD("========================================="); + fit = mFullSkeletonTree.begin(); + while (fit != mFullSkeletonTree.end()) { + if (fit.number_of_children() > 1) { + std::cout << "ppp = " << fit->point + << " " << fit->weight << std::endl; + for(uint i=1; iweight; + double w2 = mFullSkeletonTree.child(fit, i)->weight; + // if (w1 <= 0.1) break; + // if (w2 <= 0.1) break; + DD(w1);DD(w2); + if (w1>w2) { + double sw=w1; + w1=w2; w2=sw; + } + DD(w1/w2); + if (w1/w2 > 0.7) { + DD("KEEP IT"); + } + } + } + ++fit; } + + if (0) { + // Track from the first point + StartNewStep("Start tracking OLD"); + std::vector listOfBifurcations; + m_SkeletonTree.set_head(index); + TrackFromThisIndex(listOfBifurcations, skeleton, index, label, m_SkeletonTree.begin()); + DD("end track"); + + // Convert index into physical point coordinates + for(unsigned int i=0; iTransformIndexToPhysicalPoint(listOfBifurcations[i].index, + listOfBifurcations[i].point); + } + // Search for the first slice that separate the bronchus // (carina). Labelize slice by slice, stop when the two points of // the skeleton ar not in the same connected component @@ -183,8 +319,35 @@ GenerateData() DD("REDO !!!!!!!!!!!!"); /** => chercher la bif qui a les plus important sous-arbres + - iterate sur m_SkeletonTree depth first, remonter au parent, add value = distance + + NONNNNN m_SkeletonTree n'est pas la structure mais l'ensemble avec ts le skeleton + + TreeIterator itt = m_SkeletonTree.begin(); + while (itt != m_SkeletonTree.end()) { + itt->weight = 0.0; + ++itt; + } + + typename tree::post_order_iterator pit = m_SkeletonTree.begin_post(); + while (pit != m_SkeletonTree.end_post()) { + typename tree::iterator parent = m_SkeletonTree.parent(pit); + double d = pit->point.EuclideanDistanceTo(parent); + DD(parent.weight); + DD(d); + parent.weight += d; + ++it; + } + + itt = m_SkeletonTree.begin(); + while (itt != m_SkeletonTree.end()) { + DD(itt->weight); + ++itt; + } **/ + + // search for carina bool stop = false; int slice_index = listOfBifurcations[0].index[2]; // first slice from carina in skeleton int i=0; @@ -240,12 +403,13 @@ GenerateData() if (GetVerboseStep()) { std::cout << "\t Found carina at " << carina_position << " mm" << std::endl; } - GetAFDB()->SetPoint3D("Carina", carina_position); + GetAFDB()->SetPoint3D("carina", carina_position); // Write bifurcation (debug) for(uint i=0; iSetPoint3D("Bif"+toString(i), listOfBifurcations[i].point); } + } // Set the output (skeleton); this->GraftOutput(skeleton); // not SetNthOutput @@ -253,6 +417,70 @@ GenerateData() //-------------------------------------------------------------------- +//-------------------------------------------------------------------- +/** + Start from the pixel "index" in the image "skeleton" and track the + path from neighbor pixels. Put the tracked path in the tree at the + currentNode position. Label is used to mark the FG pixel as already + visited. Progress recursively when several neighbors are found. + **/ +template +void +clitk::ExtractAirwaysTreeInfoFilter:: +TrackFromThisIndex2(ImageIndexType index, ImagePointer skeleton, + ImagePixelType label, + typename FullTreeType::iterator currentNode, + typename StructuralTreeType::iterator currentSNode) +{ + // Create NeighborhoodIterator + typedef itk::NeighborhoodIterator NeighborhoodIteratorType; + typename NeighborhoodIteratorType::SizeType radius; + radius.Fill(1); + NeighborhoodIteratorType nit(radius, skeleton, skeleton->GetLargestPossibleRegion()); + + // Track + std::vector listOfTrackedPoint; + bool stop = false; + while (!stop) { + nit.SetLocation(index); + nit.SetCenterPixel(label); + listOfTrackedPoint.clear(); + for(unsigned int i=0; iTransformIndexToPhysicalPoint(n.index, n.point); + currentNode = mFullSkeletonTree.append_child(currentNode, n); + skeleton->SetPixel(n.index, label); // change label in skeleton image + } + else { + if (listOfTrackedPoint.size() >= 2) { + for(uint i=0; iTransformIndexToPhysicalPoint(n.index, n.point); + typename FullTreeType::iterator node = mFullSkeletonTree.append_child(currentNode, n); + StructuralTreeNodeType sn; + sn.index = listOfTrackedPoint[i]; + skeleton->TransformIndexToPhysicalPoint(sn.index, sn.point); + typename StructuralTreeType::iterator snode = mStructuralSkeletonTree.append_child(currentSNode, sn); + TrackFromThisIndex2(listOfTrackedPoint[i], skeleton, label, node, snode); + } + } + stop = true; // this it the end of the tracking + } // end else + } // end while stop +} +//-------------------------------------------------------------------- + //-------------------------------------------------------------------- template void diff --git a/segmentation/clitkExtractBones.ggo b/segmentation/clitkExtractBones.ggo index cd737e7..fa265ef 100644 --- a/segmentation/clitkExtractBones.ggo +++ b/segmentation/clitkExtractBones.ggo @@ -40,4 +40,8 @@ option "upper2" - "Upper threshold for RG" double no default="1500" option "radius2" - "Neighborhood radius" int no multiple default="1" option "sampleRate2" - "Sample rate of label image for RG: number of voxels to skip between seeds" int no default="0" +section "Fill holes" +option "doNotFillHoles" - "Do not fill holes if set" flag on +option "dir" d "Directions (axes) to perform filling (defaults to 2,1,0)" int multiple no + option "noAutoCrop" - "If set : do no crop final mask to BoundingBox" flag on diff --git a/segmentation/clitkExtractBonesFilter.h b/segmentation/clitkExtractBonesFilter.h index e803e5a..2205ea3 100644 --- a/segmentation/clitkExtractBonesFilter.h +++ b/segmentation/clitkExtractBonesFilter.h @@ -157,7 +157,13 @@ namespace clitk { itkGetConstMacro(SampleRate2, int); GGO_DefineOption(sampleRate2, SetSampleRate2, int); - // Final Step + // Step fill holes + itkSetMacro(FillHoles, bool); + itkGetConstMacro(FillHoles, bool); + itkBooleanMacro(FillHoles); + GGO_DefineOption_Flag(doNotFillHoles, SetFillHoles); + + // Step Auto Crop itkSetMacro(AutoCrop, bool); itkGetConstMacro(AutoCrop, bool); itkBooleanMacro(AutoCrop); @@ -193,6 +199,10 @@ namespace clitk { InputImageSizeType m_Radius2; int m_SampleRate2; + // Step + bool m_FillHoles; + InputImageSizeType m_FillHolesDirections; + virtual void GenerateOutputInformation(); virtual void GenerateData(); diff --git a/segmentation/clitkExtractBonesFilter.txx b/segmentation/clitkExtractBonesFilter.txx index 3e8cc49..7f59e92 100644 --- a/segmentation/clitkExtractBonesFilter.txx +++ b/segmentation/clitkExtractBonesFilter.txx @@ -24,6 +24,7 @@ #include "clitkSetBackgroundImageFilter.h" #include "clitkSegmentationUtils.h" #include "clitkAutoCropFilter.h" +#include "clitkFillMaskFilter.h" // itk #include "itkBinaryThresholdImageFilter.h" @@ -63,6 +64,7 @@ ExtractBonesFilter(): SetRadius2(s); SetSampleRate2(0); AutoCropOn(); + FillHolesOn(); } //-------------------------------------------------------------------- @@ -90,6 +92,7 @@ SetArgsInfo(ArgsInfoType mArgsInfo) SetWriteStep_GGO(mArgsInfo); SetVerboseWarningOff_GGO(mArgsInfo); + SetAFDBFilename_GGO(mArgsInfo); SetOutputBonesFilename_GGO(mArgsInfo); SetInitialSmoothing_GGO(mArgsInfo); @@ -108,8 +111,7 @@ SetArgsInfo(ArgsInfoType mArgsInfo) SetRadius2_GGO(mArgsInfo); SetSampleRate2_GGO(mArgsInfo); SetAutoCrop_GGO(mArgsInfo); - - SetAFDBFilename_GGO(mArgsInfo); + SetFillHoles_GGO(mArgsInfo); } //-------------------------------------------------------------------- @@ -262,6 +264,20 @@ GenerateOutputInformation() { output = setBackgroundFilter->GetOutput(); + //-------------------------------------------------------------------- + //-------------------------------------------------------------------- + // Fill Bones + if (GetFillHoles()) { + StartNewStep("Fill Holes"); + typedef clitk::FillMaskFilter FillMaskFilterType; + typename FillMaskFilterType::Pointer fillMaskFilter = FillMaskFilterType::New(); + fillMaskFilter->SetInput(output); + fillMaskFilter->AddDirection(2); + fillMaskFilter->Update(); + output = fillMaskFilter->GetOutput(); + StopCurrentStep(output); + } + //-------------------------------------------------------------------- //-------------------------------------------------------------------- // [Optional] diff --git a/segmentation/clitkExtractLungFilter.txx b/segmentation/clitkExtractLungFilter.txx index 2147b88..52276f4 100644 --- a/segmentation/clitkExtractLungFilter.txx +++ b/segmentation/clitkExtractLungFilter.txx @@ -190,6 +190,7 @@ GenerateOutputInformation() cropFilter->SetCropLikeImage(patient); cropFilter->Update(); working_input = cropFilter->GetOutput(); + DD(working_input->GetLargestPossibleRegion()); StopCurrentStep(working_input); //-------------------------------------------------------------------- @@ -344,6 +345,7 @@ GenerateOutputInformation() autocropFilter2->SetInput(working_image); autocropFilter2->Update(); working_image = autocropFilter2->GetOutput(); + DD(working_image->GetLargestPossibleRegion()); StopCurrentStep(working_image); //-------------------------------------------------------------------- @@ -382,13 +384,14 @@ GenerateOutputInformation() // Fill Lungs if (GetFillHoles()) { StartNewStep("Fill Holes"); - /* + typedef clitk::FillMaskFilter FillMaskFilterType; typename FillMaskFilterType::Pointer fillMaskFilter = FillMaskFilterType::New(); - fillMaskFilter(working_image); + fillMaskFilter->SetInput(working_image); + fillMaskFilter->AddDirection(2); + //fillMaskFilter->AddDirection(1); fillMaskFilter->Update(); working_image = fillMaskFilter->GetOutput(); StopCurrentStep(working_image); - */ } //-------------------------------------------------------------------- @@ -411,13 +414,14 @@ GenerateOutputInformation() // Decompose the first label static const unsigned int Dim = ImageType::ImageDimension; if (initialNumberOfLabels<2) { + DD(initialNumberOfLabels); // Structuring element radius typename ImageType::SizeType radius; for (unsigned int i=0;i DecomposeAndReconstructFilterType; typename DecomposeAndReconstructFilterType::Pointer decomposeAndReconstructFilter=DecomposeAndReconstructFilterType::New(); decomposeAndReconstructFilter->SetInput(working_image); - decomposeAndReconstructFilter->SetVerbose(false); + decomposeAndReconstructFilter->SetVerbose(true); decomposeAndReconstructFilter->SetRadius(radius); decomposeAndReconstructFilter->SetMaximumNumberOfLabels(2); decomposeAndReconstructFilter->SetMinimumObjectSize(this->GetMinimalComponentSize()); @@ -534,7 +538,8 @@ TracheaRegionGrowing() f->SetLower(-2000); f->SetUpper(GetUpperThresholdForTrachea()); f->SetMinimumLowerThreshold(-2000); - f->SetMaximumUpperThreshold(0); + // f->SetMaximumUpperThreshold(0); // MAYBE TO CHANGE ??? + f->SetMaximumUpperThreshold(-800); // MAYBE TO CHANGE ??? f->SetAdaptLowerBorder(false); f->SetAdaptUpperBorder(true); f->SetMinimumSize(5000); @@ -542,12 +547,15 @@ TracheaRegionGrowing() f->SetMultiplier(GetMultiplierForTrachea()); f->SetThresholdStepSize(GetThresholdStepSizeForTrachea()); f->SetMinimumThresholdStepSize(1); + f->VerboseOn(); for(unsigned int i=0; iAddSeed(m_Seeds[i]); // DD(m_Seeds[i]); } f->Update(); + writeImage(f->GetOutput(), "trg.mhd"); + // take first (main) connected component trachea_tmp = Labelize(f->GetOutput(), GetBackgroundValue(), diff --git a/segmentation/clitkExtractLymphStations.ggo b/segmentation/clitkExtractLymphStations.ggo index 3e493a1..1f3f546 100644 --- a/segmentation/clitkExtractLymphStations.ggo +++ b/segmentation/clitkExtractLymphStations.ggo @@ -14,14 +14,9 @@ option "verboseWarningOff" - "Do not display warning" flag off section "I/O" option "afdb" a "Input Anatomical Feature DB" string no -option "mediastinum" m "Input mediastinum mask filename" string yes -option "trachea" t "Input trachea mask filename" string yes -option "output" o "Output lungs mask filename" string yes - -option "carenaZposition" c "Sup-Inf position of the carena (in mm, with origin)" double no -option "middleLobeBronchusZposition" b "Sup-Inf position of the middle lobe bronchus (in mm, with origin)" double no -option "spacing" - "Intermediate resampling spacing" double no default="6" -option "fuzzy1" - "Fuzzy relative position threshold" double no default="0.6" +option "input" i "Input filename" string no +option "output" o "Output lungs mask filename" string no + diff --git a/segmentation/clitkExtractLymphStationsFilter.h b/segmentation/clitkExtractLymphStationsFilter.h index 6db26a2..d130157 100644 --- a/segmentation/clitkExtractLymphStationsFilter.h +++ b/segmentation/clitkExtractLymphStationsFilter.h @@ -28,10 +28,7 @@ namespace clitk { //-------------------------------------------------------------------- /* Try to extract the LymphStations part of a thorax CT. - Inputs : - - Patient label image - - Lungs label image - - Bones label image + Need a set of Anatomical Features (AFDB) */ //-------------------------------------------------------------------- @@ -39,12 +36,12 @@ namespace clitk { class ITK_EXPORT ExtractLymphStationsFilter: public virtual clitk::FilterBase, public clitk::FilterWithAnatomicalFeatureDatabaseManagement, - public itk::ImageToImageFilter + public itk::ImageToImageFilter > { public: /** Standard class typedefs. */ - typedef itk::ImageToImageFilter Superclass; + typedef itk::ImageToImageFilter > Superclass; typedef ExtractLymphStationsFilter Self; typedef itk::SmartPointer Pointer; typedef itk::SmartPointer ConstPointer; @@ -53,8 +50,7 @@ namespace clitk { itkNewMacro(Self); /** Run-time type information (and related methods). */ - itkTypeMacro(ExtractLymphStationsFilter, InPlaceImageFilter); - FILTERBASE_INIT; + itkTypeMacro(ExtractLymphStationsFilter, ImageToImageFilter); /** Some convenient typedefs. */ typedef TImageType ImageType; @@ -66,46 +62,33 @@ namespace clitk { typedef typename ImageType::IndexType ImageIndexType; typedef typename ImageType::PointType ImagePointType; - /** Connect inputs */ - void SetInputMediastinumLabelImage(const TImageType * image, ImagePixelType bg=0); - void SetInputTracheaLabelImage(const TImageType * image, ImagePixelType bg=0); - - /** ImageDimension constants */ - itkStaticConstMacro(ImageDimension, unsigned int, TImageType::ImageDimension); + typedef uchar MaskImagePixelType; + typedef itk::Image MaskImageType; + typedef typename MaskImageType::Pointer MaskImagePointer; + typedef typename MaskImageType::RegionType MaskImageRegionType; + typedef typename MaskImageType::SizeType MaskImageSizeType; + typedef typename MaskImageType::IndexType MaskImageIndexType; + typedef typename MaskImageType::PointType MaskImagePointType; - // Set all options at a time - template - void SetArgsInfo(ArgsInfoType arg); + /** ImageDimension constants */ + itkStaticConstMacro(ImageDimension, unsigned int, ImageType::ImageDimension); + FILTERBASE_INIT; - // Background / Foreground - itkSetMacro(BackgroundValueMediastinum, ImagePixelType); - itkGetConstMacro(BackgroundValueMediastinum, ImagePixelType); - //GGO_DefineOption(MediastinumBG, SetBackgroundValueMediastinum, ImagePixelType); + /** Main options (from ggo) */ + template + void SetArgsInfo(ArgsInfoType & argsinfo); + + itkGetConstMacro(BackgroundValue, MaskImagePixelType); + itkGetConstMacro(ForegroundValue, MaskImagePixelType); + itkSetMacro(BackgroundValue, MaskImagePixelType); + itkSetMacro(ForegroundValue, MaskImagePixelType); - itkSetMacro(BackgroundValueTrachea, ImagePixelType); - itkGetConstMacro(BackgroundValueTrachea, ImagePixelType); - //GGO_DefineOption(TracheaBG, SetBackgroundValueTrachea, ImagePixelType); - - itkGetConstMacro(BackgroundValue, ImagePixelType); - itkGetConstMacro(ForegroundValue, ImagePixelType); - - //itkSetMacro(CarinaZPositionInMM, double); - void SetCarinaZPositionInMM(double d) { m_CarinaZPositionInMM = d; Modified(); m_CarinaZPositionInMMIsSet = true; } - itkGetConstMacro(CarinaZPositionInMM, double); - GGO_DefineOption(carenaZposition, SetCarinaZPositionInMM, double); + // Station 7 + itkSetMacro(FuzzyThreshold, double); + itkGetConstMacro(FuzzyThreshold, double); + itkSetMacro(Station7Filename, std::string); + itkGetConstMacro(Station7Filename, std::string); - itkSetMacro(MiddleLobeBronchusZPositionInMM, double); - itkGetConstMacro(MiddleLobeBronchusZPositionInMM, double); - GGO_DefineOption(middleLobeBronchusZposition, SetMiddleLobeBronchusZPositionInMM, double); - - itkSetMacro(IntermediateSpacing, double); - itkGetConstMacro(IntermediateSpacing, double); - GGO_DefineOption(spacing, SetIntermediateSpacing, double); - - itkSetMacro(FuzzyThreshold1, double); - itkGetConstMacro(FuzzyThreshold1, double); - GGO_DefineOption(fuzzy1, SetFuzzyThreshold1, double); - protected: ExtractLymphStationsFilter(); virtual ~ExtractLymphStationsFilter() {} @@ -113,27 +96,35 @@ namespace clitk { virtual void GenerateOutputInformation(); virtual void GenerateInputRequestedRegion(); virtual void GenerateData(); - - itkSetMacro(BackgroundValue, ImagePixelType); - itkSetMacro(ForegroundValue, ImagePixelType); - - ImageConstPointer m_mediastinum; - ImageConstPointer m_trachea; - ImagePointer m_working_image; - ImagePointer m_working_trachea; - ImagePointer m_output; - - ImagePixelType m_BackgroundValueMediastinum; - ImagePixelType m_BackgroundValueTrachea; - ImagePixelType m_BackgroundValue; - ImagePixelType m_ForegroundValue; - - double m_CarinaZPositionInMM; - bool m_CarinaZPositionInMMIsSet; - double m_MiddleLobeBronchusZPositionInMM; - double m_IntermediateSpacing; - double m_FuzzyThreshold1; + ImageConstPointer m_Input; + MaskImagePointer m_Support; + MaskImagePointer m_Working_Support; + MaskImagePointer m_Output; + MaskImagePixelType m_BackgroundValue; + MaskImagePixelType m_ForegroundValue; + + // Common + MaskImagePointer m_Trachea; + + // Station 7 + void ExtractStation_7(); + void ExtractStation_7_SI_Limits(); + void ExtractStation_7_RL_Limits(); + void ExtractStation_7_Posterior_Limits(); + std::string m_Station7Filename; + MaskImagePointer m_working_trachea; + double m_FuzzyThreshold; + MaskImagePointer m_LeftBronchus; + MaskImagePointer m_RightBronchus; + MaskImagePointer m_Station7; + + // Station 4RL + void ExtractStation_4RL(); + void ExtractStation_4RL_SI_Limits(); + void ExtractStation_4RL_LR_Limits(); + MaskImagePointer m_Station4RL; + private: ExtractLymphStationsFilter(const Self&); //purposely not implemented void operator=(const Self&); //purposely not implemented @@ -146,6 +137,8 @@ namespace clitk { #ifndef ITK_MANUAL_INSTANTIATION #include "clitkExtractLymphStationsFilter.txx" +#include "clitkExtractLymphStation_7.txx" +#include "clitkExtractLymphStation_4RL.txx" #endif #endif diff --git a/segmentation/clitkExtractLymphStationsFilter.txx b/segmentation/clitkExtractLymphStationsFilter.txx index b81ae07..7b2fd2f 100644 --- a/segmentation/clitkExtractLymphStationsFilter.txx +++ b/segmentation/clitkExtractLymphStationsFilter.txx @@ -34,6 +34,7 @@ #include #include #include +#include #include // itk ENST @@ -45,27 +46,26 @@ clitk::ExtractLymphStationsFilter:: ExtractLymphStationsFilter(): clitk::FilterBase(), clitk::FilterWithAnatomicalFeatureDatabaseManagement(), - itk::ImageToImageFilter() + itk::ImageToImageFilter() { this->SetNumberOfRequiredInputs(1); - SetBackgroundValueMediastinum(0); - SetBackgroundValueTrachea(0); SetBackgroundValue(0); SetForegroundValue(1); - SetIntermediateSpacing(6); - SetFuzzyThreshold1(0.6); - m_CarinaZPositionInMMIsSet = false; + + // Station 7 + SetFuzzyThreshold(0.5); + SetStation7Filename("station7.mhd"); } //-------------------------------------------------------------------- //-------------------------------------------------------------------- template +template void clitk::ExtractLymphStationsFilter:: -SetInputMediastinumLabelImage(const TImageType * image, ImagePixelType bg) { - this->SetNthInput(0, const_cast(image)); - SetBackgroundValueMediastinum(bg); +SetArgsInfo(ArgsInfoType & argsinfo) { + DD("SetArgsInfo"); } //-------------------------------------------------------------------- @@ -74,30 +74,43 @@ SetInputMediastinumLabelImage(const TImageType * image, ImagePixelType bg) { template void clitk::ExtractLymphStationsFilter:: -SetInputTracheaLabelImage(const TImageType * image, ImagePixelType bg) { - this->SetNthInput(1, const_cast(image)); - SetBackgroundValueTrachea(bg); +GenerateOutputInformation() { + // Get inputs + LoadAFDB(); + m_Input = dynamic_cast(itk::ProcessObject::GetInput(0)); + m_Support = GetAFDB()->template GetImage ("mediastinum"); + + // + typedef clitk::BooleanOperatorLabelImageFilter BFilter; + BFilter::Pointer merge = BFilter::New(); + + // Extract Station7 + ExtractStation_7(); + m_Output = m_Station7; + + // Extract Station4RL + ExtractStation_4RL(); + + writeImage(m_Station4RL, "s4rl.mhd"); + // writeImage(m_Output, "ouput.mhd"); + //writeImage(m_Working_Support, "ws.mhd"); + /*merge->SetInput1(m_Station7); + merge->SetInput2(m_Station4RL); // support + merge->SetOperationType(BFilter::AndNot); CHANGE OPERATOR + merge->SetForegroundValue(4); + merge->Update(); + m_Output = merge->GetOutput(); + */ } //-------------------------------------------------------------------- //-------------------------------------------------------------------- template -template void clitk::ExtractLymphStationsFilter:: -SetArgsInfo(ArgsInfoType mArgsInfo) -{ - SetVerboseOption_GGO(mArgsInfo); - SetVerboseStep_GGO(mArgsInfo); - SetWriteStep_GGO(mArgsInfo); - SetVerboseWarningOff_GGO(mArgsInfo); - SetAFDBFilename_GGO(mArgsInfo); - SetCarinaZPositionInMM_GGO(mArgsInfo); - m_CarinaZPositionInMMIsSet = false; - SetMiddleLobeBronchusZPositionInMM_GGO(mArgsInfo); - SetIntermediateSpacing_GGO(mArgsInfo); - SetFuzzyThreshold1_GGO(mArgsInfo); +GenerateInputRequestedRegion() { + DD("GenerateInputRequestedRegion (nothing?)"); } //-------------------------------------------------------------------- @@ -106,189 +119,24 @@ SetArgsInfo(ArgsInfoType mArgsInfo) template void clitk::ExtractLymphStationsFilter:: -GenerateOutputInformation() { - // Superclass::GenerateOutputInformation(); - - // Get input - m_mediastinum = dynamic_cast(itk::ProcessObject::GetInput(0)); - m_trachea = dynamic_cast(itk::ProcessObject::GetInput(1)); - - // Get landmarks information - if (!m_CarinaZPositionInMMIsSet) { - ImagePointType carina; - LoadAFDB(); - GetAFDB()->GetPoint3D("Carina", carina); - DD(carina); - m_CarinaZPositionInMM = carina[2]; - } - DD(m_CarinaZPositionInMM); - - // ---------------------------------------------------------------- - // ---------------------------------------------------------------- - // Superior limit = carina - // Inferior limit = origin right middle lobe bronchus - StartNewStep("Inf/Sup mediastinum limits with carina/bronchus"); - ImageRegionType region = m_mediastinum->GetLargestPossibleRegion(); DD(region); - ImageSizeType size = region.GetSize(); - ImageIndexType index = region.GetIndex(); - DD(m_CarinaZPositionInMM); - DD(m_MiddleLobeBronchusZPositionInMM); - index[2] = floor((m_MiddleLobeBronchusZPositionInMM - m_mediastinum->GetOrigin()[2]) / m_mediastinum->GetSpacing()[2]); - size[2] = 1+ceil((m_CarinaZPositionInMM-m_MiddleLobeBronchusZPositionInMM) / m_mediastinum->GetSpacing()[2]); // +1 to - region.SetSize(size); - region.SetIndex(index); - DD(region); - typedef itk::RegionOfInterestImageFilter CropFilterType; - typename CropFilterType::Pointer cropFilter = CropFilterType::New(); - cropFilter->SetInput(m_mediastinum); - cropFilter->SetRegionOfInterest(region); - cropFilter->Update(); - m_working_image = cropFilter->GetOutput(); - // Auto Crop (because following rel pos is faster) - m_working_image = clitk::AutoCrop(m_working_image, 0); - StopCurrentStep(m_working_image); - m_output = m_working_image; - - // ---------------------------------------------------------------- - // ---------------------------------------------------------------- - // Separate trachea in two CCL - StartNewStep("Separate trachea under carina"); - // DD(region); - ImageRegionType trachea_region = m_trachea->GetLargestPossibleRegion(); - for(int i=0; i<3; i++) { - index[i] = floor(((index[i]*m_mediastinum->GetSpacing()[i])+m_mediastinum->GetOrigin()[i] - -m_trachea->GetOrigin()[i])/m_trachea->GetSpacing()[i]); - // DD(index[i]); - size[i] = ceil((size[i]*m_mediastinum->GetSpacing()[i])/m_trachea->GetSpacing()[i]); - // DD(size[i]); - if (index[i] < 0) { - size[i] += index[i]; - index[i] = 0; - } - if (size[i]+index[i] > (trachea_region.GetSize()[i] + trachea_region.GetIndex()[i])) { - size[i] = trachea_region.GetSize()[i] + trachea_region.GetIndex()[i] - index[i]; - } - } - // DD(index); - // DD(size); - region.SetIndex(index); - region.SetSize(size); - // typedef itk::RegionOfInterestImageFilter CropFilterType; - // typename CropFilterType::Pointer - cropFilter = CropFilterType::New(); - // m_trachea.Print(std::cout); - cropFilter->SetInput(m_trachea); - cropFilter->SetRegionOfInterest(region); - cropFilter->Update(); - m_working_trachea = cropFilter->GetOutput(); - - // Labelize and consider two main labels - m_working_trachea = Labelize(m_working_trachea, 0, true, 1); - - // Detect wich label is at Left - typedef itk::ImageSliceConstIteratorWithIndex SliceIteratorType; - SliceIteratorType iter(m_working_trachea, m_working_trachea->GetLargestPossibleRegion()); - iter.SetFirstDirection(0); - iter.SetSecondDirection(1); - iter.GoToBegin(); - bool stop = false; - ImagePixelType leftLabel; - ImagePixelType rightLabel; - while (!stop) { - if (iter.Get() != m_BackgroundValueTrachea) { - // DD(iter.GetIndex()); - // DD((int)iter.Get()); - leftLabel = iter.Get(); - stop = true; - } - ++iter; - } - if (leftLabel == 1) rightLabel = 2; - else rightLabel = 1; - DD((int)leftLabel); - DD((int)rightLabel); - - // End step - StopCurrentStep(m_working_trachea); - - //----------------------------------------------------- - /* DD("TEST Skeleton"); - typedef itk::BinaryThinningImageFilter SkeletonFilterType; - typename SkeletonFilterType::Pointer skeletonFilter = SkeletonFilterType::New(); - skeletonFilter->SetInput(m_working_trachea); - skeletonFilter->Update(); - writeImage(skeletonFilter->GetOutput(), "skel.mhd"); - writeImage(skeletonFilter->GetThinning(), "skel2.mhd"); - */ - - //----------------------------------------------------- - StartNewStep("Left limits with bronchus (slice by slice)"); - // Select LeftLabel (set right label to 0) - ImagePointer temp = SetBackground(m_working_trachea, m_working_trachea, rightLabel, 0); - writeImage(temp, "temp1.mhd"); - - typedef clitk::SliceBySliceRelativePositionFilter SliceRelPosFilterType; - typename SliceRelPosFilterType::Pointer sliceRelPosFilter = SliceRelPosFilterType::New(); - sliceRelPosFilter->SetCurrentStepBaseId(this->GetCurrentStepId()); - sliceRelPosFilter->VerboseStepOn(); - sliceRelPosFilter->WriteStepOn(); - sliceRelPosFilter->SetInput(m_working_image); - sliceRelPosFilter->SetInputObject(temp); - sliceRelPosFilter->SetDirection(2); - sliceRelPosFilter->SetFuzzyThreshold(0.6); - sliceRelPosFilter->SetOrientationType(SliceRelPosFilterType::RelPosFilterType::RightTo); - sliceRelPosFilter->Update(); - m_working_image = sliceRelPosFilter->GetOutput(); - writeImage(m_working_image, "afterslicebyslice.mhd"); - - - //----------------------------------------------------- - StartNewStep("Right limits with bronchus (slice by slice)"); - // Select LeftLabel (set right label to 0) - temp = SetBackground(m_working_trachea, m_working_trachea, leftLabel, 0); - writeImage(temp, "temp2.mhd"); - - sliceRelPosFilter = SliceRelPosFilterType::New(); - sliceRelPosFilter->SetCurrentStepBaseId(this->GetCurrentStepId()); - sliceRelPosFilter->VerboseStepOn(); - sliceRelPosFilter->WriteStepOn(); - sliceRelPosFilter->SetInput(m_working_image); - sliceRelPosFilter->SetInputObject(temp); - sliceRelPosFilter->SetDirection(2); - sliceRelPosFilter->SetFuzzyThreshold(0.6); - sliceRelPosFilter->SetOrientationType(SliceRelPosFilterType::RelPosFilterType::LeftTo); - sliceRelPosFilter->Update(); - m_working_image = sliceRelPosFilter->GetOutput(); - writeImage(m_working_image, "afterslicebyslice.mhd"); - - - DD("end"); - m_output = m_working_image; - StopCurrentStep(m_output); +GenerateData() { + DD("GenerateData, graft output"); - // Set output image information (required) - ImagePointer outputImage = this->GetOutput(0); - outputImage->SetRegions(m_working_image->GetLargestPossibleRegion()); - outputImage->SetOrigin(m_working_image->GetOrigin()); - outputImage->SetRequestedRegion(m_working_image->GetLargestPossibleRegion()); - DD("end2"); + // Final Step -> graft output (if SetNthOutput => redo) + this->GraftOutput(m_Output); } //-------------------------------------------------------------------- - + //-------------------------------------------------------------------- template void clitk::ExtractLymphStationsFilter:: -GenerateInputRequestedRegion() { - DD("GenerateInputRequestedRegion"); - // Call default - Superclass::GenerateInputRequestedRegion(); - // Following needed because output region can be greater than input (trachea) - ImagePointer mediastinum = dynamic_cast(itk::ProcessObject::GetInput(0)); - ImagePointer trachea = dynamic_cast(itk::ProcessObject::GetInput(1)); - mediastinum->SetRequestedRegion(mediastinum->GetLargestPossibleRegion()); - trachea->SetRequestedRegion(trachea->GetLargestPossibleRegion()); +ExtractStation_7() { + DD("ExtractStation_7"); + ExtractStation_7_SI_Limits(); + ExtractStation_7_RL_Limits(); + ExtractStation_7_Posterior_Limits(); } //-------------------------------------------------------------------- @@ -297,12 +145,20 @@ GenerateInputRequestedRegion() { template void clitk::ExtractLymphStationsFilter:: -GenerateData() { - DD("GenerateData"); - // Final Step -> graft output (if SetNthOutput => redo) - this->GraftOutput(m_output); +ExtractStation_4RL() { + DD("ExtractStation_4RL"); + writeImage(m_Support, "essai.mhd"); // OK + + /* + WARNING ONLY 4R FIRST !!! (not same inf limits) + */ + + ExtractStation_4RL_SI_Limits(); + ExtractStation_4RL_LR_Limits(); + } //-------------------------------------------------------------------- - + + #endif //#define CLITKBOOLEANOPERATORLABELIMAGEFILTER_TXX diff --git a/segmentation/clitkExtractLymphStationsGenericFilter.h b/segmentation/clitkExtractLymphStationsGenericFilter.h index fa2d18c..c5a9949 100644 --- a/segmentation/clitkExtractLymphStationsGenericFilter.h +++ b/segmentation/clitkExtractLymphStationsGenericFilter.h @@ -37,10 +37,10 @@ namespace clitk ExtractLymphStationsGenericFilter(); //-------------------------------------------------------------------- - typedef ExtractLymphStationsGenericFilter Self; typedef ImageToImageGenericFilter > Superclass; - typedef itk::SmartPointer Pointer; - typedef itk::SmartPointer ConstPointer; + typedef ExtractLymphStationsGenericFilter Self; + typedef itk::SmartPointer Pointer; + typedef itk::SmartPointer ConstPointer; //-------------------------------------------------------------------- itkNewMacro(Self); @@ -48,6 +48,8 @@ namespace clitk //-------------------------------------------------------------------- void SetArgsInfo(const ArgsInfoType & a); + template + void SetOptionsFromArgsInfoToFilter(FilterType * f) ; //-------------------------------------------------------------------- // Main function called each time the filter is updated @@ -58,6 +60,10 @@ namespace clitk template void InitializeImageType(); ArgsInfoType mArgsInfo; + private: + ExtractLymphStationsGenericFilter(const Self&); //purposely not implemented + void operator=(const Self&); //purposely not implemented + }; // end class //-------------------------------------------------------------------- diff --git a/segmentation/clitkExtractLymphStationsGenericFilter.txx b/segmentation/clitkExtractLymphStationsGenericFilter.txx index fb6d4da..9dfa850 100644 --- a/segmentation/clitkExtractLymphStationsGenericFilter.txx +++ b/segmentation/clitkExtractLymphStationsGenericFilter.txx @@ -28,7 +28,7 @@ clitk::ExtractLymphStationsGenericFilter::ExtractLymphStationsGene { // Default values cmdline_parser_clitkExtractLymphStations_init(&mArgsInfo); - InitializeImageType<3>(); + InitializeImageType<3>(); // Only for 3D images } //-------------------------------------------------------------------- @@ -38,10 +38,7 @@ template template void clitk::ExtractLymphStationsGenericFilter::InitializeImageType() { - ADD_IMAGE_TYPE(Dim, uchar); - // ADD_IMAGE_TYPE(Dim, short); - // ADD_IMAGE_TYPE(Dim, int); - // ADD_IMAGE_TYPE(Dim, float); + ADD_IMAGE_TYPE(Dim, short); // Can add float later } //-------------------------------------------------------------------- @@ -53,9 +50,23 @@ void clitk::ExtractLymphStationsGenericFilter::SetArgsInfo(const A mArgsInfo=a; SetIOVerbose(mArgsInfo.verbose_flag); if (mArgsInfo.imagetypes_flag) this->PrintAvailableImageTypes(); - if (mArgsInfo.mediastinum_given) AddInputFilename(mArgsInfo.mediastinum_arg); - if (mArgsInfo.trachea_given) AddInputFilename(mArgsInfo.trachea_arg); - if (mArgsInfo.output_given) AddOutputFilename(mArgsInfo.output_arg); + if (mArgsInfo.input_given) AddInputFilename(mArgsInfo.input_arg); + if (mArgsInfo.output_given) AddOutputFilename(mArgsInfo.output_arg); +} +//-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template +template +void +clitk::ExtractLymphStationsGenericFilter:: +SetOptionsFromArgsInfoToFilter(FilterType * f) +{ + f->SetVerboseOption(mArgsInfo.verbose_flag); + f->SetVerboseStep(mArgsInfo.verboseStep_flag); + f->SetWriteStep(mArgsInfo.writeStep_flag); + f->SetAFDBFilename(mArgsInfo.afdb_arg); } //-------------------------------------------------------------------- @@ -68,34 +79,24 @@ template void clitk::ExtractLymphStationsGenericFilter::UpdateWithInputImageType() { // Reading input - typename ImageType::Pointer mediastinum = this->template GetInput(0); - typename ImageType::Pointer trachea = this->template GetInput(1); + typename ImageType::Pointer input = this->template GetInput(0); // Create filter typedef clitk::ExtractLymphStationsFilter FilterType; typename FilterType::Pointer filter = FilterType::New(); // Set global Options - filter->SetInputMediastinumLabelImage(mediastinum, 0); // change 0 with BG - filter->SetInputTracheaLabelImage(trachea, 0); // change 0 with BG - filter->SetArgsInfo(mArgsInfo); + filter->SetInput(input); + SetOptionsFromArgsInfoToFilter(filter); // Go ! - // try { - filter->Update(); - // } - // catch(std::runtime_error e) { - - // // Check if error - // if (filter->HasError()) { - // SetLastError(filter->GetLastError()); - // // No output - // return; - // } + filter->Update(); // Write/Save results - typename ImageType::Pointer output = filter->GetOutput(); - this->template SetNextOutput(output); + typedef uchar MaskImagePixelType; + typedef itk::Image OutputImageType; + typename OutputImageType::Pointer output = filter->GetOutput(); + this->template SetNextOutput(output); } //-------------------------------------------------------------------- diff --git a/segmentation/clitkExtractMediastinum.ggo b/segmentation/clitkExtractMediastinum.ggo index 4a42029..48139f5 100644 --- a/segmentation/clitkExtractMediastinum.ggo +++ b/segmentation/clitkExtractMediastinum.ggo @@ -13,26 +13,23 @@ option "verboseWarningOff" - "Do not display warning" flag off section "I/O" -option "patient" p "Input patient mask filename" string yes -option "patientBG" - "Patient Background" int default="0" no -option "bones" b "Input bones mask filename" string yes -option "bonesBG" - "Bones Background" int default="0" no -option "lung" l "Input lung mask filename" string yes -option "lungBG" - "Lung Background" int default="0" no -option "lungLeft" - "Lung Left value" int default="1" no -option "lungRight" - "Lung Right value" int default="2" no -option "trachea" t "Input trachea mask filename" string yes -option "tracheaBG" - "Trachea Background" int default="0" no -option "afdb" a "Input Anatomical Feature DB" string no +option "input" i "Input CT filename" string yes +option "afdb" a "Input Anatomical Feature DB" string no +option "useBones" - "If set : do use bones mask (when image is not injected)" flag off option "output" o "Output lungs mask filename" string yes section "Step 1 : Left/Right limits with lungs" option "spacing" - "Intermediate resampling spacing" double no default="6" -option "fuzzy1" - "Fuzzy relative position threshold" double no default="0.4" +option "fuzzy1" - "Fuzzy relative position threshold" double no default="0.5" section "Step 2 : Ant/Post limits with bones" option "fuzzy2" - "Fuzzy relative position threshold" double no default="0.6" +option "antSpine" - "Distance max to anterior part of the spine in mm" double no default="10" section "Step 3 : Inf limits with Lung" -option "fuzzy3" - "Fuzzy relative position threshold" double no default="0.2" +option "fuzzy3" - "Fuzzy relative position threshold" double no default="0.05" + +section "Step x : threshold for removing bones and injected parts" +option "upper" - "Upper threshold" double no default="150" +option "lower" - "Lower threshold" double no default="-200" \ No newline at end of file diff --git a/segmentation/clitkExtractMediastinumFilter.h b/segmentation/clitkExtractMediastinumFilter.h index 4b64d13..f0d24b4 100644 --- a/segmentation/clitkExtractMediastinumFilter.h +++ b/segmentation/clitkExtractMediastinumFilter.h @@ -39,12 +39,32 @@ namespace clitk { class ITK_EXPORT ExtractMediastinumFilter: public virtual clitk::FilterBase, public clitk::FilterWithAnatomicalFeatureDatabaseManagement, - public itk::ImageToImageFilter + public itk::ImageToImageFilter > { public: + /** Some convenient typedefs. */ + typedef TImageType ImageType; + typedef typename ImageType::ConstPointer ImageConstPointer; + typedef typename ImageType::Pointer ImagePointer; + typedef typename ImageType::RegionType ImageRegionType; + typedef typename ImageType::PixelType ImagePixelType; + typedef typename ImageType::SizeType ImageSizeType; + typedef typename ImageType::IndexType ImageIndexType; + typedef typename ImageType::PointType ImagePointType; + + /** Some convenient typedefs. */ + typedef uchar MaskImagePixelType; + typedef itk::Image MaskImageType; + typedef typename MaskImageType::ConstPointer MaskImageConstPointer; + typedef typename MaskImageType::Pointer MaskImagePointer; + typedef typename MaskImageType::RegionType MaskImageRegionType; + typedef typename MaskImageType::SizeType MaskImageSizeType; + typedef typename MaskImageType::IndexType MaskImageIndexType; + typedef typename MaskImageType::PointType MaskImagePointType; + /** Standard class typedefs. */ - typedef itk::ImageToImageFilter Superclass; + typedef itk::ImageToImageFilter Superclass; typedef ExtractMediastinumFilter Self; typedef itk::SmartPointer Pointer; typedef itk::SmartPointer ConstPointer; @@ -56,57 +76,52 @@ namespace clitk { itkTypeMacro(ExtractMediastinumFilter, InPlaceImageFilter); FILTERBASE_INIT; - /** Some convenient typedefs. */ - typedef TImageType ImageType; - typedef typename ImageType::ConstPointer ImageConstPointer; - typedef typename ImageType::Pointer ImagePointer; - typedef typename ImageType::RegionType ImageRegionType; - typedef typename ImageType::PixelType ImagePixelType; - typedef typename ImageType::SizeType ImageSizeType; - typedef typename ImageType::IndexType ImageIndexType; - typedef typename ImageType::PointType ImagePointType; - /** Connect inputs */ - void SetInputPatientLabelImage(const TImageType * image, ImagePixelType bg=0); - void SetInputLungLabelImage(const TImageType * image, ImagePixelType bg=0, - ImagePixelType fgLeftLung=1, ImagePixelType fgRightLung=2); - void SetInputBonesLabelImage(const TImageType * image, ImagePixelType bg=0); - void SetInputTracheaLabelImage(const TImageType * image, ImagePixelType bg=0); + void SetInput(const ImageType * image); + void SetInputPatientLabelImage(const MaskImageType * image, MaskImagePixelType bg=0); + void SetInputLungLabelImage(const MaskImageType * image, MaskImagePixelType bg=0, + MaskImagePixelType fgLeftLung=1, MaskImagePixelType fgRightLung=2); + void SetInputBonesLabelImage(const MaskImageType * image, MaskImagePixelType bg=0); + void SetInputTracheaLabelImage(const MaskImageType * image, MaskImagePixelType bg=0); + // Output filename (for AFBD) + itkSetMacro(OutputMediastinumFilename, std::string); + itkGetConstMacro(OutputMediastinumFilename, std::string); + /** ImageDimension constants */ - itkStaticConstMacro(ImageDimension, unsigned int, TImageType::ImageDimension); + itkStaticConstMacro(ImageDimension, unsigned int, ImageType::ImageDimension); // Set all options at a time template void SetArgsInfo(ArgsInfoType arg); // Background / Foreground - itkSetMacro(BackgroundValuePatient, ImagePixelType); - itkGetConstMacro(BackgroundValuePatient, ImagePixelType); - GGO_DefineOption(patientBG, SetBackgroundValuePatient, ImagePixelType); + itkSetMacro(BackgroundValuePatient, MaskImagePixelType); + itkGetConstMacro(BackgroundValuePatient, MaskImagePixelType); + // GGO_DefineOption(patientBG, SetBackgroundValuePatient, MaskImagePixelType); - itkSetMacro(BackgroundValueLung, ImagePixelType); - itkGetConstMacro(BackgroundValueLung, ImagePixelType); - GGO_DefineOption(lungBG, SetBackgroundValueLung, ImagePixelType); + itkSetMacro(BackgroundValueLung, MaskImagePixelType); + itkGetConstMacro(BackgroundValueLung, MaskImagePixelType); + // GGO_DefineOption(lungBG, SetBackgroundValueLung, MaskImagePixelType); - itkSetMacro(BackgroundValueBones, ImagePixelType); - itkGetConstMacro(BackgroundValueBones, ImagePixelType); - GGO_DefineOption(bonesBG, SetBackgroundValueBones, ImagePixelType); + itkSetMacro(BackgroundValueBones, MaskImagePixelType); + itkGetConstMacro(BackgroundValueBones, MaskImagePixelType); + // GGO_DefineOption(bonesBG, SetBackgroundValueBones, MaskImagePixelType); - itkGetConstMacro(BackgroundValue, ImagePixelType); - itkGetConstMacro(ForegroundValue, ImagePixelType); + itkGetConstMacro(BackgroundValue, MaskImagePixelType); + itkGetConstMacro(ForegroundValue, MaskImagePixelType); - itkSetMacro(ForegroundValueLeftLung, ImagePixelType); - itkGetConstMacro(ForegroundValueLeftLung, ImagePixelType); - GGO_DefineOption(lungLeft, SetForegroundValueLeftLung, ImagePixelType); + itkSetMacro(ForegroundValueLeftLung, MaskImagePixelType); + itkGetConstMacro(ForegroundValueLeftLung, MaskImagePixelType); + // GGO_DefineOption(lungLeft, SetForegroundValueLeftLung, MaskImagePixelType); - itkSetMacro(ForegroundValueRightLung, ImagePixelType); - itkGetConstMacro(ForegroundValueRightLung, ImagePixelType); - GGO_DefineOption(lungRight, SetForegroundValueRightLung, ImagePixelType); + itkSetMacro(ForegroundValueRightLung, MaskImagePixelType); + itkGetConstMacro(ForegroundValueRightLung, MaskImagePixelType); + // GGO_DefineOption(lungRight, SetForegroundValueRightLung, MaskImagePixelType); - itkSetMacro(BackgroundValueTrachea, ImagePixelType); - itkGetConstMacro(BackgroundValueTrachea, ImagePixelType); - GGO_DefineOption(lungBG, SetBackgroundValueTrachea, ImagePixelType); + itkSetMacro(BackgroundValueTrachea, MaskImagePixelType); + itkGetConstMacro(BackgroundValueTrachea, MaskImagePixelType); + // GGO_DefineOption(lungBG, SetBackgroundValueTrachea, MaskImagePixelType); itkSetMacro(IntermediateSpacing, double); itkGetConstMacro(IntermediateSpacing, double); @@ -124,6 +139,23 @@ namespace clitk { itkGetConstMacro(FuzzyThreshold3, double); GGO_DefineOption(fuzzy3, SetFuzzyThreshold3, double); + itkSetMacro(DistanceMaxToAnteriorPartOfTheSpine, double); + itkGetConstMacro(DistanceMaxToAnteriorPartOfTheSpine, double); + GGO_DefineOption(antSpine, SetDistanceMaxToAnteriorPartOfTheSpine, double); + + itkBooleanMacro(UseBones); + itkSetMacro(UseBones, bool); + itkGetConstMacro(UseBones, bool); + GGO_DefineOption_Flag(useBones, SetUseBones); + + itkSetMacro(UpperThreshold, double); + itkGetConstMacro(UpperThreshold, double); + GGO_DefineOption(upper, SetUpperThreshold, double); + + itkSetMacro(LowerThreshold, double); + itkGetConstMacro(LowerThreshold, double); + GGO_DefineOption(lower, SetLowerThreshold, double); + protected: ExtractMediastinumFilter(); virtual ~ExtractMediastinumFilter() {} @@ -132,23 +164,31 @@ namespace clitk { virtual void GenerateInputRequestedRegion(); virtual void GenerateData(); - itkSetMacro(BackgroundValue, ImagePixelType); - itkSetMacro(ForegroundValue, ImagePixelType); - - ImagePixelType m_BackgroundValuePatient; - ImagePixelType m_BackgroundValueLung; - ImagePixelType m_BackgroundValueBones; - ImagePixelType m_BackgroundValueTrachea; - ImagePixelType m_ForegroundValueLeftLung; - ImagePixelType m_ForegroundValueRightLung; - - ImagePixelType m_BackgroundValue; - ImagePixelType m_ForegroundValue; + itkSetMacro(BackgroundValue, MaskImagePixelType); + itkSetMacro(ForegroundValue, MaskImagePixelType); + MaskImagePixelType m_BackgroundValuePatient; + MaskImagePixelType m_BackgroundValueLung; + MaskImagePixelType m_BackgroundValueBones; + MaskImagePixelType m_BackgroundValueTrachea; + MaskImagePixelType m_ForegroundValueLeftLung; + MaskImagePixelType m_ForegroundValueRightLung; + + MaskImagePixelType m_BackgroundValue; + MaskImagePixelType m_ForegroundValue; + + typename MaskImageType::Pointer output; + double m_IntermediateSpacing; double m_FuzzyThreshold1; double m_FuzzyThreshold2; double m_FuzzyThreshold3; + double m_DistanceMaxToAnteriorPartOfTheSpine; + bool m_UseBones; + double m_UpperThreshold; + double m_LowerThreshold; + + std::string m_OutputMediastinumFilename; private: ExtractMediastinumFilter(const Self&); //purposely not implemented diff --git a/segmentation/clitkExtractMediastinumFilter.txx b/segmentation/clitkExtractMediastinumFilter.txx index b76e617..7c77546 100644 --- a/segmentation/clitkExtractMediastinumFilter.txx +++ b/segmentation/clitkExtractMediastinumFilter.txx @@ -23,6 +23,7 @@ #include "clitkCommon.h" #include "clitkExtractMediastinumFilter.h" #include "clitkAddRelativePositionConstraintToLabelImageFilter.h" +#include "clitkSliceBySliceRelativePositionFilter.h" #include "clitkSegmentationUtils.h" #include "clitkExtractAirwaysTreeInfoFilter.h" #include "clitkCropLikeImageFilter.h" @@ -35,6 +36,7 @@ #include "itkLabelImageToStatisticsLabelMapFilter.h" #include "itkRegionOfInterestImageFilter.h" #include "itkBinaryThresholdImageFilter.h" +#include "itkScalarImageKmeansImageFilter.h" // itk ENST #include "RelativePositionPropImageFilter.h" @@ -45,9 +47,9 @@ clitk::ExtractMediastinumFilter:: ExtractMediastinumFilter(): clitk::FilterBase(), clitk::FilterWithAnatomicalFeatureDatabaseManagement(), - itk::ImageToImageFilter() + itk::ImageToImageFilter() { - this->SetNumberOfRequiredInputs(4); + this->SetNumberOfRequiredInputs(1); SetBackgroundValuePatient(0); SetBackgroundValueLung(0); @@ -58,9 +60,14 @@ ExtractMediastinumFilter(): SetForegroundValue(1); SetIntermediateSpacing(6); - SetFuzzyThreshold1(0.4); + SetFuzzyThreshold1(0.5); SetFuzzyThreshold2(0.6); - SetFuzzyThreshold3(0.2); + SetFuzzyThreshold3(0.05); + + SetDistanceMaxToAnteriorPartOfTheSpine(10); + SetOutputMediastinumFilename("mediastinum.mhd"); + + UseBonesOff(); } //-------------------------------------------------------------------- @@ -69,9 +76,9 @@ ExtractMediastinumFilter(): template void clitk::ExtractMediastinumFilter:: -SetInputPatientLabelImage(const ImageType * image, ImagePixelType bg) +SetInputPatientLabelImage(const MaskImageType * image, MaskImagePixelType bg) { - this->SetNthInput(0, const_cast(image)); + this->SetNthInput(0, const_cast(image)); m_BackgroundValuePatient = bg; } //-------------------------------------------------------------------- @@ -81,10 +88,10 @@ SetInputPatientLabelImage(const ImageType * image, ImagePixelType bg) template void clitk::ExtractMediastinumFilter:: -SetInputLungLabelImage(const ImageType * image, ImagePixelType bg, - ImagePixelType fgLeft, ImagePixelType fgRight) +SetInputLungLabelImage(const MaskImageType * image, MaskImagePixelType bg, + MaskImagePixelType fgLeft, MaskImagePixelType fgRight) { - this->SetNthInput(1, const_cast(image)); + this->SetNthInput(1, const_cast(image)); m_BackgroundValueLung = bg; m_ForegroundValueLeftLung = fgLeft; m_ForegroundValueRightLung = fgRight; @@ -96,9 +103,9 @@ SetInputLungLabelImage(const ImageType * image, ImagePixelType bg, template void clitk::ExtractMediastinumFilter:: -SetInputBonesLabelImage(const ImageType * image, ImagePixelType bg) +SetInputBonesLabelImage(const MaskImageType * image, MaskImagePixelType bg) { - this->SetNthInput(2, const_cast(image)); + this->SetNthInput(2, const_cast(image)); m_BackgroundValueBones = bg; } //-------------------------------------------------------------------- @@ -108,9 +115,9 @@ SetInputBonesLabelImage(const ImageType * image, ImagePixelType bg) template void clitk::ExtractMediastinumFilter:: -SetInputTracheaLabelImage(const ImageType * image, ImagePixelType bg) +SetInputTracheaLabelImage(const MaskImageType * image, MaskImagePixelType bg) { - this->SetNthInput(3, const_cast(image)); + this->SetNthInput(3, const_cast(image)); m_BackgroundValueTrachea = bg; } //-------------------------------------------------------------------- @@ -128,19 +135,17 @@ SetArgsInfo(ArgsInfoType mArgsInfo) SetWriteStep_GGO(mArgsInfo); SetVerboseWarningOff_GGO(mArgsInfo); - SetBackgroundValuePatient_GGO(mArgsInfo); - SetBackgroundValueLung_GGO(mArgsInfo); - SetBackgroundValueTrachea_GGO(mArgsInfo); - - SetForegroundValueLeftLung_GGO(mArgsInfo); - SetForegroundValueRightLung_GGO(mArgsInfo); - SetIntermediateSpacing_GGO(mArgsInfo); SetFuzzyThreshold1_GGO(mArgsInfo); SetFuzzyThreshold2_GGO(mArgsInfo); SetFuzzyThreshold3_GGO(mArgsInfo); - SetAFDBFilename_GGO(mArgsInfo); + SetAFDBFilename_GGO(mArgsInfo); + SetDistanceMaxToAnteriorPartOfTheSpine_GGO(mArgsInfo); + SetUseBones_GGO(mArgsInfo); + + SetLowerThreshold_GGO(mArgsInfo); + SetUpperThreshold_GGO(mArgsInfo); } //-------------------------------------------------------------------- @@ -149,11 +154,14 @@ SetArgsInfo(ArgsInfoType mArgsInfo) template void clitk::ExtractMediastinumFilter:: -GenerateOutputInformation() { - ImagePointer input = dynamic_cast(itk::ProcessObject::GetInput(0)); - ImagePointer outputImage = this->GetOutput(0); - outputImage->SetRegions(outputImage->GetLargestPossibleRegion()); +GenerateInputRequestedRegion() +{ + //DD("GenerateInputRequestedRegion"); + // Do not call default + // Superclass::GenerateInputRequestedRegion(); + // DD("End GenerateInputRequestedRegion"); } + //-------------------------------------------------------------------- @@ -161,22 +169,9 @@ GenerateOutputInformation() { template void clitk::ExtractMediastinumFilter:: -GenerateInputRequestedRegion() +SetInput(const ImageType * image) { - // Call default - Superclass::GenerateInputRequestedRegion(); - - // Get input pointers - LoadAFDB(); - ImagePointer patient = GetAFDB()->template GetImage ("patient"); - ImagePointer lung = GetAFDB()->template GetImage ("lungs"); - ImagePointer bones = GetAFDB()->template GetImage ("bones"); - ImagePointer trachea = GetAFDB()->template GetImage ("trachea"); - - patient->SetRequestedRegion(patient->GetLargestPossibleRegion()); - lung->SetRequestedRegion(lung->GetLargestPossibleRegion()); - bones->SetRequestedRegion(bones->GetLargestPossibleRegion()); - trachea->SetRequestedRegion(trachea->GetLargestPossibleRegion()); + this->SetNthInput(0, const_cast(image)); } //-------------------------------------------------------------------- @@ -185,69 +180,83 @@ GenerateInputRequestedRegion() template void clitk::ExtractMediastinumFilter:: -GenerateData() -{ +GenerateOutputInformation() { + // DD("GenerateOutputInformation"); + // Do not call default + // Superclass::GenerateOutputInformation(); + + //-------------------------------------------------------------------- // Get input pointers - ImagePointer patient = GetAFDB()->template GetImage ("patient"); - ImagePointer lung = GetAFDB()->template GetImage ("lungs"); - ImagePointer bones = GetAFDB()->template GetImage ("bones"); - ImagePointer trachea = GetAFDB()->template GetImage ("trachea"); + LoadAFDB(); + ImageConstPointer input = dynamic_cast(itk::ProcessObject::GetInput(0)); + MaskImagePointer patient = GetAFDB()->template GetImage ("patient"); + MaskImagePointer lung = GetAFDB()->template GetImage ("lungs"); + MaskImagePointer bones = GetAFDB()->template GetImage ("bones"); + MaskImagePointer trachea = GetAFDB()->template GetImage ("trachea"); - // Get output pointer - ImagePointer output; - - // Step 0: Crop support (patient) to lung extend in RL + //-------------------------------------------------------------------- + // Step 1: Crop support (patient) to lung extend in RL StartNewStep("Crop support like lungs along LR"); - typedef clitk::CropLikeImageFilter CropFilterType; + typedef clitk::CropLikeImageFilter 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(output); - - // Step 0: Crop support (previous) to bones extend in AP - StartNewStep("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(output); - - // Step 1: patient minus lungs, minus bones - StartNewStep("Patient contours minus lungs and minus bones"); - typedef clitk::BooleanOperatorLabelImageFilter BoolFilterType; + this->template StopCurrentStep(output); + + //-------------------------------------------------------------------- + // Step 2: Crop support (previous) to bones extend in AP + if (GetUseBones()) { + StartNewStep("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(output); + } + + //-------------------------------------------------------------------- + // Step 3: patient minus lungs, minus bones, minus trachea + StartNewStep("Patient contours minus lungs, bones, trachea"); + typedef clitk::BooleanOperatorLabelImageFilter BoolFilterType; typename BoolFilterType::Pointer boolFilter = BoolFilterType::New(); boolFilter->InPlaceOn(); 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(); // Auto crop to gain support area - output = clitk::AutoCrop(output, GetBackgroundValue()); - this->template StopCurrentStep(output); - - // Step 2: LR limits from lung (need separate lung ?) + output = clitk::AutoCrop(output, GetBackgroundValue()); + this->template StopCurrentStep(output); + + //-------------------------------------------------------------------- + // Step 4: 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"); - /* - // WE DO NOT NEED THE FOLLOWING ? - // 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). - ImagePointer right_lung = clitk::SetBackground(lung, lung, 2, 0); - ImagePointer left_lung = clitk::SetBackground(lung, lung, 1, 0); - writeImage(right_lung, "right.mhd"); - writeImage(left_lung, "left.mhd"); + ImagePointer right_lung = clitk::SetBackground(lung, lung, 2, 0); + ImagePointer left_lung = clitk::SetBackground(lung, lung, 1, 0); + writeImage(right_lung, "right.mhd"); + writeImage(left_lung, "left.mhd"); */ - - typedef clitk::AddRelativePositionConstraintToLabelImageFilter RelPosFilterType; + typedef clitk::AddRelativePositionConstraintToLabelImageFilter RelPosFilterType; typename RelPosFilterType::Pointer relPosFilter = RelPosFilterType::New(); relPosFilter->SetCurrentStepBaseId(this->GetCurrentStepId()); relPosFilter->VerboseStepOff(); @@ -260,7 +269,7 @@ GenerateData() relPosFilter->SetFuzzyThreshold(GetFuzzyThreshold1()); relPosFilter->Update(); output = relPosFilter->GetOutput(); - // writeImage(right_lung, "step4-left.mhd"); + //writeImage(right_lung, "step4-left.mhd"); relPosFilter->SetInput(output); relPosFilter->SetCurrentStepBaseId(this->GetCurrentStepId()); @@ -273,118 +282,378 @@ GenerateData() relPosFilter->SetFuzzyThreshold(GetFuzzyThreshold1()); relPosFilter->Update(); output = relPosFilter->GetOutput(); - this->template StopCurrentStep(output); - - // Step 3: AP limits from bones - StartNewStep("Ant/Post limits with bones"); - ImageConstPointer bones_ant; - ImageConstPointer bones_post; - - // Find ant-post coordinate of trachea (assume the carena position is a - // good estimation of the ant-post position of the trachea) - ImagePointType carina; - LoadAFDB(); - GetAFDB()->GetPoint3D("Carina", carina); - DD(carina); - ImageIndexType index_trachea; - bones->TransformPhysicalPointToIndex(carina, index_trachea); - DD(index_trachea); + this->template StopCurrentStep(output); + + //-------------------------------------------------------------------- + // Step 5: AP limits from bones + // Separate the bones in the ant-post middle + MaskImageConstPointer bones_ant; + MaskImageConstPointer bones_post; + MaskImagePointType middle_AntPost__position; + 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; + MaskImageIndexType index_bones_middle; + bones->TransformPhysicalPointToIndex(middle_AntPost__position, index_bones_middle); - // Split bone image first into two parts (ant and post) - typedef itk::RegionOfInterestImageFilter ROIFilterType; - // typedef itk::ExtractImageFilter ROIFilterType; - typename ROIFilterType::Pointer roiFilter = ROIFilterType::New(); - ImageRegionType region = bones->GetLargestPossibleRegion(); - ImageSizeType size = region.GetSize(); - DD(size); - size[1] = index_trachea[1]; //size[1]/2.0; - DD(size); - region.SetSize(size); - roiFilter->SetInput(bones); - // roiFilter->SetExtractionRegion(region); - roiFilter->SetRegionOfInterest(region); - roiFilter->ReleaseDataFlagOff(); - roiFilter->Update(); - bones_ant = roiFilter->GetOutput(); - writeImage(bones_ant, "b_ant.mhd"); - - // roiFilter->ResetPipeline();// = ROIFilterType::New(); - roiFilter = ROIFilterType::New(); - ImageIndexType index = region.GetIndex(); - index[1] = bones->GetLargestPossibleRegion().GetIndex()[1] + size[1]-1; - size[1] = bones->GetLargestPossibleRegion().GetSize()[1] - size[1]; - DD(size); - region.SetIndex(index); - region.SetSize(size); - roiFilter->SetInput(bones); - // roiFilter->SetExtractionRegion(region); - roiFilter->SetRegionOfInterest(region); - roiFilter->ReleaseDataFlagOff(); - roiFilter->Update(); - bones_post = roiFilter->GetOutput(); - writeImage(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(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(); - this->template StopCurrentStep(output); - - // Get CCL + // Split bone image first into two parts (ant and post), and crop + // lateraly to get vertebral + typedef itk::RegionOfInterestImageFilter ROIFilterType; + // typedef itk::ExtractImageFilter 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(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(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(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(); + this->template StopCurrentStep(output); + } + + //-------------------------------------------------------------------- + // Step 6: Get CCL StartNewStep("Keep main connected component"); - output = clitk::Labelize(output, GetBackgroundValue(), true, 500); - // output = RemoveLabels(output, BG, param->GetLabelsToRemove()); - output = clitk::KeepLabels(output, GetBackgroundValue(), - GetForegroundValue(), 1, 1, 0); - this->template StopCurrentStep(output); - - // Step : Lower limits from lung (need separate lung ?) - StartNewStep("Lower limits with lungs"); - relPosFilter = RelPosFilterType::New(); - relPosFilter->SetCurrentStepBaseId(this->GetCurrentStepId()); - relPosFilter->VerboseStepOff(); - relPosFilter->WriteStepOff(); - relPosFilter->SetInput(output); - DD(output->GetLargestPossibleRegion().GetIndex()); - // relPosFilter->SetInputObject(left_lung); - relPosFilter->SetInputObject(lung); - relPosFilter->SetOrientationType(RelPosFilterType::SupTo); - relPosFilter->SetIntermediateSpacing(GetIntermediateSpacing()); - relPosFilter->SetFuzzyThreshold(GetFuzzyThreshold3()); - relPosFilter->Update(); - output = relPosFilter->GetOutput(); - DD(output->GetLargestPossibleRegion()); + output = clitk::Labelize(output, GetBackgroundValue(), false, 500); + // output = RemoveLabels(output, BG, param->GetLabelsToRemove()); + output = clitk::KeepLabels(output, GetBackgroundValue(), + GetForegroundValue(), 1, 1, 0); + this->template StopCurrentStep(output); + + //-------------------------------------------------------------------- + // Step 7 : Slice by Slice to optimize posterior part + // Warning slice does not necesseraly correspond between 'output' and 'bones' + typedef clitk::ExtractSliceFilter ExtractSliceFilterType; + typename ExtractSliceFilterType::Pointer extractSliceFilter = ExtractSliceFilterType::New(); + typedef typename ExtractSliceFilterType::SliceType SliceType; + std::vector mSlices; + if (GetUseBones()) { + StartNewStep("Rafine posterior part according to vertebral body"); + extractSliceFilter->SetInput(bones_post); + extractSliceFilter->SetDirection(2); + extractSliceFilter->Update(); + std::vector mVertebralAntPositionBySlice; + extractSliceFilter->GetOutputSlices(mSlices); + for(unsigned int i=0; i(mSlices[i], 0, true, 10); + mSlices[i] = KeepLabels(mSlices[i], + GetBackgroundValue(), + GetForegroundValue(), 1, 2, true); // keep two first + // Find most anterior point (start of the vertebral) + typename itk::ImageRegionIteratorWithIndex + 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 ?"); + } + } + + // 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; iGetOrigin()[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 it(output, region); + it.GoToBegin(); + while (!it.IsAtEnd()) { + it.Set(GetBackgroundValue()); + ++it; + } + } + } + } + this->template StopCurrentStep(output); + } + + //-------------------------------------------------------------------- + // Step 8: Trial segmentation KMeans + if (0) { + StartNewStep("K means"); + // Take input, crop like current mask + typedef CropLikeImageFilter CropLikeFilterType; + typename CropLikeFilterType::Pointer cropLikeFilter = CropLikeFilterType::New(); + cropLikeFilter->SetInput(input); + cropLikeFilter->SetCropLikeImage(output); + cropLikeFilter->Update(); + ImagePointer working_input = cropLikeFilter->GetOutput(); + writeImage(working_input, "crop-input.mhd"); + // Set bG at -1000 + working_input = clitk::SetBackground(working_input, output, GetBackgroundValue(), -1000); + writeImage(working_input, "crop-input2.mhd"); + // Kmeans + typedef itk::ScalarImageKmeansImageFilter 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(kmeans, kmeans, + 1, GetBackgroundValue()); + writeImage(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(output, "out-kmean.mhd"); + this->template StopCurrentStep(output); + + // TODO -> FillMASK ? + // comment speed ? mask ? 2 class ? + + + //TODO + // Confidence connected ? + + } + + //-------------------------------------------------------------------- + // Step 8: Lower limits from lung (need separate lung ?) + if (1) { + // StartNewStep("Trial : minus segmented struct"); + // MaskImagePointer heart = GetAFDB()->template GetImage ("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 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(output); + } + + //-------------------------------------------------------------------- + // Step 10: Slice by Slice CCL + StartNewStep("Slice by Slice keep only one component"); + typedef clitk::ExtractSliceFilter ExtractSliceFilterType; + // typename ExtractSliceFilterType::Pointer + extractSliceFilter = ExtractSliceFilterType::New(); + extractSliceFilter->SetInput(output); + extractSliceFilter->SetDirection(2); + extractSliceFilter->Update(); + typedef typename ExtractSliceFilterType::SliceType SliceType; + // std::vector + mSlices.clear(); + extractSliceFilter->GetOutputSlices(mSlices); + for(unsigned int i=0; i(mSlices[i], 0, true, 100); + mSlices[i] = KeepLabels(mSlices[i], 0, 1, 1, 1, true); + } + typedef itk::JoinSeriesImageFilter JoinSeriesFilterType; + typename JoinSeriesFilterType::Pointer joinFilter = JoinSeriesFilterType::New(); + joinFilter->SetOrigin(output->GetOrigin()[2]); + joinFilter->SetSpacing(output->GetSpacing()[2]); + for(unsigned int i=0; iPushBackInput(mSlices[i]); + } + joinFilter->Update(); + output = joinFilter->GetOutput(); + this->template StopCurrentStep(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 CropLikeFilterType; + typename CropLikeFilterType::Pointer cropLikeFilter = CropLikeFilterType::New(); + cropLikeFilter->SetInput(input); + cropLikeFilter->SetCropLikeImage(output); + cropLikeFilter->Update(); + ImagePointer working_input = cropLikeFilter->GetOutput(); + writeImage(working_input, "crop-ct.mhd"); + // Binarize + typedef itk::BinaryThresholdImageFilter 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(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(output); + + //-------------------------------------------------------------------- + // Step 10 : AutoCrop + StartNewStep("AutoCrop"); + output = clitk::AutoCrop(output, GetBackgroundValue()); + this->template StopCurrentStep(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()); + this->GetOutput(0)->SetRequestedRegion(output->GetLargestPossibleRegion()); + this->GetOutput(0)->SetBufferedRegion(output->GetLargestPossibleRegion()); +} +//-------------------------------------------------------------------- - output = clitk::AutoCrop(output, GetBackgroundValue()); - // roiFilter = ROIFilterType::New(); - //roiFilter->SetInput(output); - //roiFilter->Update(); - //output = roiFilter->GetOutput(); - // Final Step -> set output - this->SetNthOutput(0, output); +//-------------------------------------------------------------------- +template +void +clitk::ExtractMediastinumFilter:: +GenerateData() +{ + DD("GenerateData"); + this->GraftOutput(output); + // Store image filenames into AFDB + GetAFDB()->SetImageFilename("mediastinum", this->GetOutputMediastinumFilename()); + WriteAFDB(); } //-------------------------------------------------------------------- diff --git a/segmentation/clitkExtractMediastinumGenericFilter.txx b/segmentation/clitkExtractMediastinumGenericFilter.txx index b5cfaba..881ebba 100644 --- a/segmentation/clitkExtractMediastinumGenericFilter.txx +++ b/segmentation/clitkExtractMediastinumGenericFilter.txx @@ -38,7 +38,7 @@ template template void clitk::ExtractMediastinumGenericFilter::InitializeImageType() { - ADD_IMAGE_TYPE(Dim, uchar); + ADD_IMAGE_TYPE(Dim, short); // ADD_IMAGE_TYPE(Dim, short); // ADD_IMAGE_TYPE(Dim, int); // ADD_IMAGE_TYPE(Dim, float); @@ -53,10 +53,11 @@ void clitk::ExtractMediastinumGenericFilter::SetArgsInfo(const Arg mArgsInfo=a; SetIOVerbose(mArgsInfo.verbose_flag); if (mArgsInfo.imagetypes_flag) this->PrintAvailableImageTypes(); - if (mArgsInfo.patient_given) AddInputFilename(mArgsInfo.patient_arg); - if (mArgsInfo.lung_given) AddInputFilename(mArgsInfo.lung_arg); - if (mArgsInfo.bones_given) AddInputFilename(mArgsInfo.bones_arg); - if (mArgsInfo.trachea_given) AddInputFilename(mArgsInfo.trachea_arg); + if (mArgsInfo.input_given) AddInputFilename(mArgsInfo.input_arg); + //if (mArgsInfo.patient_given) AddInputFilename(mArgsInfo.patient_arg); + //if (mArgsInfo.lung_given) AddInputFilename(mArgsInfo.lung_arg); + //if (mArgsInfo.bones_given) AddInputFilename(mArgsInfo.bones_arg); + //if (mArgsInfo.trachea_given) AddInputFilename(mArgsInfo.trachea_arg); if (mArgsInfo.output_given) AddOutputFilename(mArgsInfo.output_arg); } //-------------------------------------------------------------------- @@ -70,28 +71,30 @@ template void clitk::ExtractMediastinumGenericFilter::UpdateWithInputImageType() { // Reading input - typename ImageType::Pointer patient = this->template GetInput(0); - typename ImageType::Pointer lung = this->template GetInput(1); - typename ImageType::Pointer bones = this->template GetInput(2); - typename ImageType::Pointer trachea = this->template GetInput(3); + typename ImageType::Pointer input = this->template GetInput(0); + // typename ImageType::Pointer patient = this->template GetInput(0); + // typename ImageType::Pointer lung = this->template GetInput(1); + // typename ImageType::Pointer bones = this->template GetInput(2); + // typename ImageType::Pointer trachea = this->template GetInput(3); // Create filter typedef clitk::ExtractMediastinumFilter FilterType; typename FilterType::Pointer filter = FilterType::New(); // Set global Options - filter->SetInputPatientLabelImage(patient, mArgsInfo.patientBG_arg); - filter->SetInputLungLabelImage(lung, mArgsInfo.lungBG_arg, mArgsInfo.lungRight_arg, mArgsInfo.lungLeft_arg); - filter->SetInputBonesLabelImage(bones, mArgsInfo.bonesBG_arg); - filter->SetInputTracheaLabelImage(trachea, mArgsInfo.tracheaBG_arg); + filter->SetInput(input); + // filter->SetInputPatientLabelImage(patient, mArgsInfo.patientBG_arg); + // filter->SetInputLungLabelImage(lung, mArgsInfo.lungBG_arg, mArgsInfo.lungRight_arg, mArgsInfo.lungLeft_arg); + // filter->SetInputBonesLabelImage(bones, mArgsInfo.bonesBG_arg); + filter->SetOutputMediastinumFilename(mArgsInfo.output_arg); filter->SetArgsInfo(mArgsInfo); // Go ! filter->Update(); // Write/Save results - typename ImageType::Pointer output = filter->GetOutput(); - this->template SetNextOutput(output); + typename FilterType::MaskImageType::Pointer output = filter->GetOutput(); + this->template SetNextOutput(output); } //-------------------------------------------------------------------- -- 2.47.1