From 907b0bad00cbf772fbf362879c74d673253f97bb Mon Sep 17 00:00:00 2001 From: dsarrut Date: Thu, 24 Mar 2011 07:33:42 +0000 Subject: [PATCH] Add S3A, correct minors for S7 --- segmentation/clitkExtractLung.ggo | 6 +- segmentation/clitkExtractLymphStation_3A.txx | 66 +++ segmentation/clitkExtractLymphStation_3P.txx | 63 ++- segmentation/clitkExtractLymphStation_7.txx | 466 +++++++++++------- segmentation/clitkExtractLymphStation_8.txx | 27 +- segmentation/clitkExtractLymphStations.ggo | 7 + .../clitkExtractLymphStationsFilter.h | 26 +- .../clitkExtractLymphStationsFilter.txx | 78 ++- ...clitkExtractLymphStationsGenericFilter.txx | 8 + 9 files changed, 511 insertions(+), 236 deletions(-) create mode 100644 segmentation/clitkExtractLymphStation_3A.txx diff --git a/segmentation/clitkExtractLung.ggo b/segmentation/clitkExtractLung.ggo index a5cf62a..5f53b4d 100644 --- a/segmentation/clitkExtractLung.ggo +++ b/segmentation/clitkExtractLung.ggo @@ -16,9 +16,9 @@ option "verboseMemory" - "Display memory usage" flag off section "I/O" option "input" i "Input CT image filename" string yes -option "afdb" a "Output Anatomical Feature DB (Carina position)" string no -option "output" o "Output lungs mask filename" string yes -option "outputTrachea" t "Output trachea mask filename" string no +option "afdb" a "Output Anatomical Feature DB (Carina position)" string no default="default.afdb" +option "output" o "Output lungs mask filename" string no default="lung.mhd" +option "outputTrachea" t "Output trachea mask filename" string no default="trachea.mhd" section "Step 1 : Air remove" diff --git a/segmentation/clitkExtractLymphStation_3A.txx b/segmentation/clitkExtractLymphStation_3A.txx new file mode 100644 index 0000000..3e70059 --- /dev/null +++ b/segmentation/clitkExtractLymphStation_3A.txx @@ -0,0 +1,66 @@ + +#include +#include + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +ExtractStation_3A_SetDefaultValues() +{ +} +//-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +ExtractStation_3A_SI_Limits() +{ + // Apex of the chest or Sternum & Carina. + StartNewStep("[Station 3A] Inf/Sup limits with Sternum and Carina"); + + // Get Carina position (has been determined in Station8) + m_CarinaZ = GetAFDB()->GetDouble("CarinaZ"); + + // Get Sternum and search for the upper position + MaskImagePointer Sternum = GetAFDB()->template GetImage("Sternum"); + + // Search most sup point + MaskImagePointType ps = Sternum->GetOrigin(); // initialise to avoid warning + clitk::FindExtremaPointInAGivenDirection(Sternum, GetBackgroundValue(), 2, false, ps); + double m_SternumZ = ps[2]+Sternum->GetSpacing()[2]; // One more slice, because it is below this point + + //* Crop support : + m_Working_Support = + clitk::CropImageAlongOneAxis(m_Working_Support, 2, + m_CarinaZ, m_SternumZ, true, + GetBackgroundValue()); + + StopCurrentStep(m_Working_Support); + m_ListOfStations["3A"] = m_Working_Support; +} +//-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +ExtractStation_3A_Ant_Limits() +{ + StartNewStep("[Station 3A] Ant limits with Sternum"); + + // Get Sternum, keep posterior part. + MaskImagePointer Sternum = GetAFDB()->template GetImage("Sternum"); + m_Working_Support = + clitk::SliceBySliceRelativePosition(m_Working_Support, Sternum, 2, + 0.5, "PostTo", + false, 3, true, false); + StopCurrentStep(m_Working_Support); + m_ListOfStations["3A"] = m_Working_Support; +} +//-------------------------------------------------------------------- + + diff --git a/segmentation/clitkExtractLymphStation_3P.txx b/segmentation/clitkExtractLymphStation_3P.txx index d1be710..eb50f76 100644 --- a/segmentation/clitkExtractLymphStation_3P.txx +++ b/segmentation/clitkExtractLymphStation_3P.txx @@ -65,23 +65,13 @@ ExtractStation_3P_Remove_Structures() StartNewStep("[Station 3P] Remove some structures."); - MaskImagePointer Aorta = GetAFDB()->template GetImage("Aorta"); - clitk::AndNot(m_Working_Support, Aorta, GetBackgroundValue()); - - MaskImagePointer VertebralBody = GetAFDB()->template GetImage("VertebralBody"); - clitk::AndNot(m_Working_Support, VertebralBody); - - MaskImagePointer SubclavianArtery = GetAFDB()->template GetImage("SubclavianArtery"); - clitk::AndNot(m_Working_Support, SubclavianArtery); - - MaskImagePointer Esophagus = GetAFDB()->template GetImage("Esophagus"); - clitk::AndNot(m_Working_Support, Esophagus); - - MaskImagePointer AzygousVein = GetAFDB()->template GetImage("AzygousVein"); - clitk::AndNot(m_Working_Support, AzygousVein); - - MaskImagePointer Thyroid = GetAFDB()->template GetImage("Thyroid"); - clitk::AndNot(m_Working_Support, Thyroid); + Remove_Structures("Aorta"); + Remove_Structures("VertebralBody"); + Remove_Structures("SubclavianArtery"); + Remove_Structures("Esophagus"); + Remove_Structures("Azygousvein"); + Remove_Structures("Thyroid"); + Remove_Structures("VertebralArtery"); StopCurrentStep(m_Working_Support); m_ListOfStations["3P"] = m_Working_Support; @@ -430,6 +420,29 @@ ExtractStation_3P_LR_sup_Limits() } //-------------------------------------------------------------------- +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +ExtractStation_3P_LR_sup_Limits_2() +{ + /* + Use VertebralArtery to limit. + + + StartNewStep("[Station 3P] Left/Right limits with VertebralArtery"); + + // Load structures + MaskImagePointer VertebralArtery = GetAFDB()->template GetImage("VertebralArtery"); + + clitk::AndNot(m_Working_Support, VertebralArtery); + + StopCurrentStep(m_Working_Support); + m_ListOfStations["3P"] = m_Working_Support; + */ +} +//-------------------------------------------------------------------- + //-------------------------------------------------------------------- template void @@ -444,10 +457,11 @@ ExtractStation_3P_LR_inf_Limits() inf : not Right to Azygousvein */ - StartNewStep("[Station 3P] Left/Right limits (inferior part) "); + StartNewStep("[Station 3P] Left/Right limits (inferior part) with Azygousvein and Aorta"); // Load structures MaskImagePointer AzygousVein = GetAFDB()->template GetImage("AzygousVein"); + MaskImagePointer Aorta = GetAFDB()->template GetImage("Aorta"); typedef clitk::AddRelativePositionConstraintToLabelImageFilter RelPosFilterType; typename RelPosFilterType::Pointer relPosFilter = RelPosFilterType::New(); @@ -465,6 +479,19 @@ ExtractStation_3P_LR_inf_Limits() relPosFilter->Update(); m_Working_Support = relPosFilter->GetOutput(); + writeImage(m_Working_Support, "before-L-aorta.mhd"); + relPosFilter->SetInput(m_Working_Support); + relPosFilter->SetInputObject(Aorta); + relPosFilter->AddOrientationTypeString("L"); + relPosFilter->SetInverseOrientationFlag(true); + // relPosFilter->SetIntermediateSpacing(3); + relPosFilter->SetIntermediateSpacingFlag(false); + relPosFilter->SetFuzzyThreshold(0.7); + relPosFilter->AutoCropFlagOn(); + relPosFilter->Update(); + m_Working_Support = relPosFilter->GetOutput(); + writeImage(m_Working_Support, "after-L-aorta.mhd"); + StopCurrentStep(m_Working_Support); m_ListOfStations["3P"] = m_Working_Support; } diff --git a/segmentation/clitkExtractLymphStation_7.txx b/segmentation/clitkExtractLymphStation_7.txx index 278efc8..725c852 100644 --- a/segmentation/clitkExtractLymphStation_7.txx +++ b/segmentation/clitkExtractLymphStation_7.txx @@ -1,32 +1,85 @@ +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +ExtractStation_7_SetDefaultValues() +{ + SetFuzzyThresholdForS7("Bronchi", 0.1); + SetFuzzyThresholdForS7("LeftSuperiorPulmonaryVein", 0.3); + SetFuzzyThresholdForS7("RightSuperiorPulmonaryVein", 0.2); + SetFuzzyThresholdForS7("RightPulmonaryArtery", 0.3); + SetFuzzyThresholdForS7("LeftPulmonaryArtery", 0.5); + SetFuzzyThresholdForS7("SVC", 0.2); +} +//-------------------------------------------------------------------- + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +SetFuzzyThresholdForS7(std::string tag, double value) +{ + m_FuzzyThresholdForS7[tag] = value; +} +//-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template +double +clitk::ExtractLymphStationsFilter:: +GetFuzzyThresholdForS7(std::string tag) +{ + if (m_FuzzyThresholdForS7.find(tag) != m_FuzzyThresholdForS7.end()) { + return m_FuzzyThresholdForS7[tag]; + } + else { + clitkExceptionMacro("Could not find options "+tag+" in the m_FuzzyThresholdForS7 list"); + } +} +//-------------------------------------------------------------------- + + //-------------------------------------------------------------------- template void clitk::ExtractLymphStationsFilter:: ExtractStation_7_SI_Limits() { + StartNewStep("[Station7] Inf/Sup mediastinum limits with carina/LLLBronchus"); // Get Inputs MaskImagePointer Trachea = GetAFDB()->template GetImage ("Trachea"); - double m_CarinaZ = GetAFDB()->GetPoint3D("Carina", 2); - DD(m_CarinaZ); - double m_OriginOfRightMiddleLobeBronchusZ = GetAFDB()->GetPoint3D("OriginOfRightMiddleLobeBronchus", 2); - DD(m_OriginOfRightMiddleLobeBronchusZ); - - /* Crop support : - Superior limit = carina - Inferior limit = origin right middle lobe bronchus */ - StartNewStep("[Station7] Inf/Sup mediastinum limits with carina/bronchus"); + + // We suppoe that CarinaZ was already computed (S8) + double m_CarinaZ = GetAFDB()->GetDouble("CarinaZ"); + + // double m_OriginOfRightMiddleLobeBronchusZ = GetAFDB()->GetPoint3D("OriginOfRightMiddleLobeBronchus", 2); + // DD(m_OriginOfRightMiddleLobeBronchusZ); + MaskImagePointer UpperBorderOfLLLBronchus = GetAFDB()->template GetImage("UpperBorderOfLLLBronchus"); + + // Search most inf point (WHY ? IS IT THE RIGHT STRUCTURE ??) + MaskImagePointType ps = UpperBorderOfLLLBronchus->GetOrigin(); // initialise to avoid warning + clitk::FindExtremaPointInAGivenDirection(UpperBorderOfLLLBronchus, GetBackgroundValue(), 2, true, ps); + double m_UpperBorderOfLLLBronchusZ = ps[2]; + + /* + std::vector centroids; + clitk::ComputeCentroids(UpperBorderOfLLLBronchus, GetBackgroundValue(), centroids); + double m_UpperBorderOfLLLBronchusZ = centroids[1][2]; + DD(m_UpperBorderOfLLLBronchusZ) + */ + + /* Crop support */ m_Working_Support = clitk::CropImageAlongOneAxis(m_Working_Support, 2, - m_OriginOfRightMiddleLobeBronchusZ, + m_UpperBorderOfLLLBronchusZ, m_CarinaZ, true, GetBackgroundValue()); - /* Crop trachea - Superior limit = carina - Inferior limit = origin right middle lobe bronchus*/ - m_working_trachea = + /* Crop trachea */ + m_Working_Trachea = clitk::CropImageAlongOneAxis(Trachea, 2, - m_OriginOfRightMiddleLobeBronchusZ, + m_UpperBorderOfLLLBronchusZ, m_CarinaZ, true, GetBackgroundValue()); @@ -43,98 +96,190 @@ clitk::ExtractLymphStationsFilter:: ExtractStation_7_RL_Limits() { // ---------------------------------------------------------------- - // Separate trachea in two CCL - StartNewStep("[Station7] Separate trachea under carina"); - - // Labelize and consider two main labels - m_working_trachea = Labelize(m_working_trachea, 0, true, 1); - - // Carina position must at the first slice that separate the two - // main bronchus (not superiorly) Check that upper slice is composed - // of at least two labels - typedef itk::ImageSliceIteratorWithIndex SliceIteratorType; - SliceIteratorType iter(m_working_trachea, m_working_trachea->GetLargestPossibleRegion()); - iter.SetFirstDirection(0); - iter.SetSecondDirection(1); - iter.GoToReverseBegin(); // Start from the end (because image is IS not SI) - int maxLabel=0; - while (!iter.IsAtReverseEndOfSlice()) { - while (!iter.IsAtReverseEndOfLine()) { - if (iter.Get() > maxLabel) maxLabel = iter.Get(); - --iter; - } - iter.PreviousLine(); + StartNewStep("[Station7] Limits with bronchus : RightTo the left bronchus"); + + // First consider bronchus and keep the main CCL by slice + m_RightBronchus = GetAFDB()->template GetImage ("RightBronchus"); + m_LeftBronchus = GetAFDB()->template GetImage ("LeftBronchus"); + + // Extract slices, Label, compute centroid, keep most central connected component + std::vector slices_leftbronchus; + std::vector slices_rightbronchus; + clitk::ExtractSlices(m_LeftBronchus, 2, slices_leftbronchus); + clitk::ExtractSlices(m_RightBronchus, 2, slices_rightbronchus); + + // Loop on slices + for(uint i=0; i(slices_leftbronchus[i], 0, false, 10); + std::vector c; + clitk::ComputeCentroids(slices_leftbronchus[i], GetBackgroundValue(), c); + if (c.size() > 1) { + double most_at_left = c[1][0]; + int most_at_left_index=1; + for(uint j=1; j(slices_leftbronchus[i], most_at_left_index, + most_at_left_index, GetBackgroundValue(), GetForegroundValue()); + } // end c.size } - if (maxLabel < 2) { - clitkExceptionMacro("First slice form Carina does not seems to seperate the two main bronchus. Abort"); + + for(uint i=0; i(slices_rightbronchus[i], 0, false, 10); + std::vector c; + clitk::ComputeCentroids(slices_rightbronchus[i], GetBackgroundValue(), c); + if (c.size() > 1) { + double most_at_right = c[1][0]; + int most_at_right_index=1; + for(uint j=1; j most_at_right) { + most_at_right = c[j][0]; + most_at_right_index = j; + } + } + // Put all other CCL to Background + slices_rightbronchus[i] = + clitk::Binarize(slices_rightbronchus[i], most_at_right_index, + most_at_right_index, GetBackgroundValue(), GetForegroundValue()); + } // end c.size } + + // Joint slices + m_LeftBronchus = clitk::JoinSlices(slices_leftbronchus, m_LeftBronchus, 2); + m_RightBronchus = clitk::JoinSlices(slices_rightbronchus, m_RightBronchus, 2); - // Compute 3D centroids of both parts to identify the left from the - // right bronchus - std::vector c; - clitk::ComputeCentroids(m_working_trachea, GetBackgroundValue(), c); - ImagePointType C1 = c[1]; - ImagePointType C2 = c[2]; - - ImagePixelType leftLabel; - ImagePixelType rightLabel; - if (C1[0] < C2[0]) { leftLabel = 1; rightLabel = 2; } - else { leftLabel = 2; rightLabel = 1; } - - StopCurrentStep(m_working_trachea); - - //----------------------------------------------------- - // Select LeftLabel (set one label to Backgroundvalue) - m_LeftBronchus = - SetBackground(m_working_trachea, m_working_trachea, - rightLabel, GetBackgroundValue(), false); - m_RightBronchus = - SetBackground(m_working_trachea, m_working_trachea, - leftLabel, GetBackgroundValue(), false); - - StartNewStep("[Station7] Limits with bronchus (slice by slice) : RightTo left bronchus"); - m_Working_Support = - clitk::SliceBySliceRelativePosition(m_Working_Support, - m_LeftBronchus, 2, - GetFuzzyThreshold(), "RightTo", - true, 4); + writeImage(m_LeftBronchus, "step-left.mhd"); + writeImage(m_RightBronchus, "step-right.mhd"); - StartNewStep("[Station7] Limits with bronchus (slice by slice) : LeftTo right bronchus"); m_Working_Support = - clitk::SliceBySliceRelativePosition(m_Working_Support, - m_RightBronchus, - 2, GetFuzzyThreshold(), "LeftTo", - true, 4); + clitk::SliceBySliceRelativePosition(m_Working_Support, m_LeftBronchus, 2, + GetFuzzyThresholdForS7("Bronchi"), "RightTo", + false, 3, false); + StopCurrentStep(m_Working_Support); - StartNewStep("[Station7] Limits with bronchus (slice by slice) : not AntTo left bronchus"); - m_Working_Support = - clitk::SliceBySliceRelativePosition(m_Working_Support, - m_LeftBronchus, - 2, GetFuzzyThreshold(), "AntTo", - true, 4, true); // NOT - StartNewStep("[Station7] Limits with bronchus (slice by slice) : not AntTo right bronchus"); + // ---------------------------------------------------------------- + StartNewStep("[Station7] Limits with bronchus : LeftTo the right bronchus"); m_Working_Support = - clitk::SliceBySliceRelativePosition(m_Working_Support, - m_RightBronchus, - 2, GetFuzzyThreshold(), "AntTo", - true, 4, true); + clitk::SliceBySliceRelativePosition(m_Working_Support, m_RightBronchus, 2, + GetFuzzyThresholdForS7("Bronchi"), "LeftTo", + false, 3, false); + StopCurrentStep(m_Working_Support); - StartNewStep("[Station7] Limits with bronchus (slice by slice) : not PostTo left bronchus"); - m_Working_Support = - clitk::SliceBySliceRelativePosition(m_Working_Support, - m_LeftBronchus, - 2, GetFuzzyThreshold(), "PostTo", - true, 4, true); - StartNewStep("[Station7] Limits with bronchus (slice by slice) : not PostTo right bronchus"); - m_Working_Support = - clitk::SliceBySliceRelativePosition(m_Working_Support, - m_RightBronchus, - 2, GetFuzzyThreshold(), "PostTo", - true, 4, true); - m_Station7 = m_Working_Support; - StopCurrentStep(m_Station7); + // ---------------------------------------------------------------- + StartNewStep("[Station7] Limits with LeftSuperiorPulmonaryVein"); + try { + MaskImagePointer LeftSuperiorPulmonaryVein = GetAFDB()->template GetImage("LeftSuperiorPulmonaryVein"); + typedef SliceBySliceRelativePositionFilter SliceRelPosFilterType; + typename SliceRelPosFilterType::Pointer sliceRelPosFilter = SliceRelPosFilterType::New(); + sliceRelPosFilter->SetInput(m_Working_Support); + sliceRelPosFilter->SetInputObject(LeftSuperiorPulmonaryVein); + sliceRelPosFilter->SetDirection(2); + sliceRelPosFilter->SetFuzzyThreshold(GetFuzzyThresholdForS7("LeftSuperiorPulmonaryVein")); + sliceRelPosFilter->AddOrientationTypeString("NotLeftTo"); + sliceRelPosFilter->AddOrientationTypeString("NotAntTo"); + sliceRelPosFilter->SetIntermediateSpacingFlag(true); + sliceRelPosFilter->SetIntermediateSpacing(3); + sliceRelPosFilter->SetUniqueConnectedComponentBySlice(false); + sliceRelPosFilter->SetAutoCropFlag(false); + sliceRelPosFilter->IgnoreEmptySliceObjectFlagOn(); + sliceRelPosFilter->Update(); + m_Working_Support = sliceRelPosFilter->GetOutput(); + StopCurrentStep(m_Working_Support); + } + catch (clitk::ExceptionObject e) { + std::cout << "Not LeftSuperiorPulmonaryVein, skip" << std::endl; + } + + // ---------------------------------------------------------------- + StartNewStep("[Station7] Limits with RightSuperiorPulmonaryVein"); + try { + MaskImagePointer RightSuperiorPulmonaryVein = GetAFDB()->template GetImage("RightSuperiorPulmonaryVein"); + typedef SliceBySliceRelativePositionFilter SliceRelPosFilterType; + typename SliceRelPosFilterType::Pointer sliceRelPosFilter = SliceRelPosFilterType::New(); + sliceRelPosFilter->SetInput(m_Working_Support); + sliceRelPosFilter->SetInputObject(RightSuperiorPulmonaryVein); + sliceRelPosFilter->SetDirection(2); + sliceRelPosFilter->SetFuzzyThreshold(GetFuzzyThresholdForS7("RightSuperiorPulmonaryVein")); + sliceRelPosFilter->AddOrientationTypeString("NotRightTo"); + sliceRelPosFilter->AddOrientationTypeString("NotAntTo"); + sliceRelPosFilter->AddOrientationTypeString("NotPostTo"); + sliceRelPosFilter->SetIntermediateSpacingFlag(true); + sliceRelPosFilter->SetIntermediateSpacing(3); + sliceRelPosFilter->SetUniqueConnectedComponentBySlice(false); + sliceRelPosFilter->SetAutoCropFlag(false); + sliceRelPosFilter->IgnoreEmptySliceObjectFlagOn(); + sliceRelPosFilter->Update(); + m_Working_Support = sliceRelPosFilter->GetOutput(); + StopCurrentStep(m_Working_Support); + } + catch (clitk::ExceptionObject e) { + std::cout << "Not RightSuperiorPulmonaryVein, skip" << std::endl; + } + + // ---------------------------------------------------------------- + StartNewStep("[Station7] Limits with RightPulmonaryArtery"); + MaskImagePointer RightPulmonaryArtery = GetAFDB()->template GetImage("RightPulmonaryArtery"); + typedef SliceBySliceRelativePositionFilter SliceRelPosFilterType; + typename SliceRelPosFilterType::Pointer sliceRelPosFilter = SliceRelPosFilterType::New(); + sliceRelPosFilter->SetInput(m_Working_Support); + sliceRelPosFilter->SetInputObject(RightPulmonaryArtery); + sliceRelPosFilter->SetDirection(2); + sliceRelPosFilter->SetFuzzyThreshold(GetFuzzyThresholdForS7("RightPulmonaryArtery")); + sliceRelPosFilter->AddOrientationTypeString("NotAntTo"); + sliceRelPosFilter->SetIntermediateSpacingFlag(true); + sliceRelPosFilter->SetIntermediateSpacing(3); + sliceRelPosFilter->SetUniqueConnectedComponentBySlice(false); + sliceRelPosFilter->SetAutoCropFlag(false); + sliceRelPosFilter->IgnoreEmptySliceObjectFlagOn(); + sliceRelPosFilter->Update(); + m_Working_Support = sliceRelPosFilter->GetOutput(); + StopCurrentStep(m_Working_Support); + + // ---------------------------------------------------------------- + StartNewStep("[Station7] Limits with LeftPulmonaryArtery"); + MaskImagePointer LeftPulmonaryArtery = GetAFDB()->template GetImage("LeftPulmonaryArtery"); + sliceRelPosFilter = SliceRelPosFilterType::New(); + sliceRelPosFilter->SetInput(m_Working_Support); + sliceRelPosFilter->SetInputObject(LeftPulmonaryArtery); + sliceRelPosFilter->SetDirection(2); + sliceRelPosFilter->SetFuzzyThreshold(GetFuzzyThresholdForS7("LeftPulmonaryArtery")); + sliceRelPosFilter->AddOrientationTypeString("NotAntTo"); + sliceRelPosFilter->SetIntermediateSpacingFlag(true); + sliceRelPosFilter->SetIntermediateSpacing(3); + sliceRelPosFilter->SetUniqueConnectedComponentBySlice(false); + sliceRelPosFilter->SetAutoCropFlag(false); + sliceRelPosFilter->IgnoreEmptySliceObjectFlagOn(); + sliceRelPosFilter->Update(); + m_Working_Support = sliceRelPosFilter->GetOutput(); + StopCurrentStep(m_Working_Support); + + StartNewStep("[Station7] Limits with SVC"); + MaskImagePointer SVC = GetAFDB()->template GetImage("SVC"); + sliceRelPosFilter = SliceRelPosFilterType::New(); + sliceRelPosFilter->SetInput(m_Working_Support); + sliceRelPosFilter->SetInputObject(SVC); + sliceRelPosFilter->SetDirection(2); + sliceRelPosFilter->SetFuzzyThreshold(GetFuzzyThresholdForS7("SVC")); + sliceRelPosFilter->AddOrientationTypeString("NotRightTo"); + sliceRelPosFilter->AddOrientationTypeString("NotAntTo"); + sliceRelPosFilter->SetIntermediateSpacingFlag(true); + sliceRelPosFilter->SetIntermediateSpacing(3); + sliceRelPosFilter->SetUniqueConnectedComponentBySlice(false); + sliceRelPosFilter->SetAutoCropFlag(true); + sliceRelPosFilter->IgnoreEmptySliceObjectFlagOn(); + sliceRelPosFilter->Update(); + m_Working_Support = sliceRelPosFilter->GetOutput(); + StopCurrentStep(m_Working_Support); + + // End + m_ListOfStations["7"] = m_Working_Support; } //-------------------------------------------------------------------- @@ -145,20 +290,31 @@ void clitk::ExtractLymphStationsFilter:: ExtractStation_7_Posterior_Limits() { - StartNewStep("[Station7] Posterior limits -> must be AntTo post wall of the bronchi"); + StartNewStep("[Station7] Posterior limits -> must be AntTo post wall of the bronchi (OLD CLASSIF)"); // Search for points that are the most left/post/ant and most // right/post/ant of the left and right bronchus // extract, loop slices, label/keep, find extrema x 3 - FindExtremaPointsInBronchus(m_LeftBronchus, 0, 15, - m_RightMostInLeftBronchus, - m_AntMostInLeftBronchus, - m_PostMostInLeftBronchus); - FindExtremaPointsInBronchus(m_RightBronchus, 1, 15, - m_LeftMostInRightBronchus, - m_AntMostInRightBronchus, - m_PostMostInRightBronchus); + /* FindExtremaPointsInBronchus(m_LeftBronchus, 0, 15, m_RightMostInLeftBronchus, + m_AntMostInLeftBronchus, m_PostMostInLeftBronchus); + FindExtremaPointsInBronchus(m_RightBronchus, 1, 15, m_LeftMostInRightBronchus, + m_AntMostInRightBronchus, m_PostMostInRightBronchus); + */ + + // First cut bronchus to the correct sup/inf support + MaskImagePointer RightBronchus = clitk::ResizeImageLike(m_RightBronchus, m_Working_Support, GetBackgroundValue()); + MaskImagePointer LeftBronchus = clitk::ResizeImageLike(m_LeftBronchus, m_Working_Support, GetBackgroundValue()); + + // Find extrema points + FindExtremaPointsInBronchus(RightBronchus, 0, 10, m_LeftMostInRightBronchus, + m_AntMostInRightBronchus, m_PostMostInRightBronchus); + + FindExtremaPointsInBronchus(LeftBronchus, 1, 10, m_RightMostInLeftBronchus, + m_AntMostInLeftBronchus, m_PostMostInLeftBronchus); + + + // DEBUG std::ofstream osrl; openFileForWriting(osrl, "osrl.txt"); osrl << "LANDMARKS1" << std::endl; std::ofstream osal; openFileForWriting(osal, "osal.txt"); osal << "LANDMARKS1" << std::endl; @@ -174,7 +330,9 @@ ExtractStation_7_Posterior_Limits() << " " << m_AntMostInLeftBronchus[i][2] << " 0 0 " << std::endl; ospl << i << " " << m_PostMostInLeftBronchus[i][0] << " " << m_PostMostInLeftBronchus[i][1] << " " << m_PostMostInLeftBronchus[i][2] << " 0 0 " << std::endl; + } + for(uint i=0; i to accelerate, start with formula, when change sign -> stop and fill - */ - typedef itk::ImageSliceIteratorWithIndex SliceIteratorType; - SliceIteratorType iter = SliceIteratorType(m_Working_Support, - m_Working_Support->GetLargestPossibleRegion()); - iter.SetFirstDirection(0); - iter.SetSecondDirection(1); - iter.GoToBegin(); - int i=0; - MaskImageType::PointType A; - MaskImageType::PointType B; - MaskImageType::PointType C; - while (!iter.IsAtEnd()) { - A = m_PostMostInLeftBronchus[i]; - B = m_PostMostInRightBronchus[i]; - C = A; - C[1] -= 10; // I know I must keep this point - double s = (B[0] - A[0]) * (C[1] - A[1]) - (B[1] - A[1]) * (C[0] - A[0]); - bool isPositive = s<0; - while (!iter.IsAtEndOfSlice()) { - while (!iter.IsAtEndOfLine()) { - // Very slow, I know ... but image should be very small - m_Working_Support->TransformIndexToPhysicalPoint(iter.GetIndex(), C); - double s = (B[0] - A[0]) * (C[1] - A[1]) - (B[1] - A[1]) * (C[0] - A[0]); - if (s == 0) iter.Set(2); - if (isPositive) { - if (s > 0) iter.Set(GetBackgroundValue()); - } - else { - if (s < 0) iter.Set(GetBackgroundValue()); - } - ++iter; - } - iter.NextLine(); - } - iter.NextSlice(); - ++i; - } - - //----------------------------------------------------- - // StartNewStep("[Station7] Anterior limits"); - + clitk::SliceBySliceSetBackgroundFromLineSeparation(m_Working_Support, + m_PostMostInRightBronchus, + m_PostMostInLeftBronchus, + GetBackgroundValue(), 1, -10); + // If needed -> can do the same with AntMost. - // MISSING FROM NOW - - // Station 4R, Station 4L, the right pulmonary artery, and/or the - // left superior pulmonary vein + // End + StopCurrentStep(m_Working_Support); + m_ListOfStations["7"] = m_Working_Support; +} +//-------------------------------------------------------------------- +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +ExtractStation_7_Remove_Structures() +{ - //----------------------------------------------------- - //----------------------------------------------------- - // ALSO SUBSTRACT ARTERY/VEIN (initially in the support) + //-------------------------------------------------------------------- + StartNewStep("[Station7] remove some structures"); + Remove_Structures("AzygousVein"); + Remove_Structures("Aorta"); + Remove_Structures("Esophagus"); + Remove_Structures("RightPulmonaryArtery"); + Remove_Structures("LeftPulmonaryArtery"); + Remove_Structures("LeftSuperiorPulmonaryVein"); + Remove_Structures("PulmonaryTrunk"); + Remove_Structures("VertebralBody"); - // Set output - m_Station7 = m_Working_Support; + // END + StopCurrentStep(m_Working_Support); + m_ListOfStations["7"] = m_Working_Support; + return; } //-------------------------------------------------------------------- + diff --git a/segmentation/clitkExtractLymphStation_8.txx b/segmentation/clitkExtractLymphStation_8.txx index dd3fdca..545780e 100644 --- a/segmentation/clitkExtractLymphStation_8.txx +++ b/segmentation/clitkExtractLymphStation_8.txx @@ -462,8 +462,8 @@ ExtractStation_8_Ant_Inf_Limits() relPosFilter->WriteStepFlagOff(); relPosFilter->SetInput(m_Working_Support); relPosFilter->SetInputObject(m_Esophagus); - relPosFilter->AddOrientationTypeString("P"); - relPosFilter->InverseOrientationFlagOff(); + relPosFilter->AddOrientationTypeString("PostTo"); + // relPosFilter->InverseOrientationFlagOff(); relPosFilter->SetDirection(2); // Z axis relPosFilter->UniqueConnectedComponentBySliceOff(); relPosFilter->SetIntermediateSpacing(3); @@ -635,7 +635,7 @@ clitk::ExtractLymphStationsFilter:: ExtractStation_8_Single_CCL_Limits() { //-------------------------------------------------------------------- - StartNewStep("[Station8] Single CCL"); + StartNewStep("[Station8 or 3P] Slice by slice, keep a single CCL (the closest to VertebralBody)"); // Consider slices std::vector slices; @@ -654,15 +654,11 @@ ExtractStation_8_Single_CCL_Limits() std::vector c; clitk::ComputeCentroids(slices[i], GetBackgroundValue(), centroids); clitk::ComputeCentroids(slices_vb[i], GetBackgroundValue(), c); - // DD(c.size()); - // DD(centroids.size()); if ((c.size() > 1) && (centroids.size() > 1)) { // keep the one which is the closer to vertebralbody centroid double distance=1000000; int index=0; - // DD(c[1]); for(uint j=1; j c; clitk::ComputeCentroids(eso_slices[i], GetBackgroundValue(), c); if (c.size() >1) { - // DD(c[1]); eso_slices[i] = clitk::CropImageAbove(eso_slices[i], 1, c[1][1], false, GetBackgroundValue()); eso_slices[i] = @@ -1041,7 +1036,6 @@ ExtractStation_8_LR_Limits() j++; } } - // DD(j); sp_MostAntVertebralBody[1] += GetDistanceMaxToAnteriorPartOfTheSpine(); // Consider offset // Crop the vertebralbody below this most post line @@ -1065,8 +1059,6 @@ ExtractStation_8_LR_Limits() clitk::PointsUtils::Convert2DTo3D(sp_MostLeftVertebralBody, VertebralBody, i, p); p_MostLeftVertebralBody.push_back(p); - // DD(p_MostLeftVertebralBody.back()); - /* ------------------------------------------------------------------------------------- Find most Left point in Esophagus */ @@ -1082,14 +1074,12 @@ ExtractStation_8_LR_Limits() // p_MostLeftEsophagus.push_back(p); // sp_MostLeftEsophagus = sp_temp; // Retrieve previous 2D position // found = true; - // DD(p); // } } else { clitk::PointsUtils::Convert2DTo3D(sp_MostLeftEsophagus, EsophagusForSlice, i, p); p_MostLeftEsophagus.push_back(p); // sp_temp = sp_MostLeftEsophagus; // Store previous 2D position - // // DD(p_MostLeftEsophagus.back()); } /* ------------------------------------------------------------------------------------- @@ -1141,11 +1131,8 @@ ExtractStation_8_Remove_Structures() //-------------------------------------------------------------------- StartNewStep("[Station8] remove some structures"); - MaskImagePointer Aorta = GetAFDB()->template GetImage("Aorta"); - clitk::AndNot(m_Working_Support, Aorta, GetBackgroundValue()); - - MaskImagePointer Esophagus = GetAFDB()->template GetImage("Esophagus"); - clitk::AndNot(m_Working_Support, Esophagus, GetBackgroundValue()); + Remove_Structures("Aorta"); + Remove_Structures("Esophagus"); // END StopCurrentStep(m_Working_Support); @@ -1243,8 +1230,8 @@ ExtractStation_8_LR_Limits_old() relPosFilter->WriteStepFlagOff(); relPosFilter->SetInput(m_Working_Support); relPosFilter->SetInputObject(Esophagus); - relPosFilter->AddOrientationTypeString("L"); - relPosFilter->InverseOrientationFlagOn(); // Not Left to + relPosFilter->AddOrientationTypeString("NotLeftTo"); + // relPosFilter->InverseOrientationFlagOn(); // Not Left to relPosFilter->SetDirection(2); // Z axis relPosFilter->UniqueConnectedComponentBySliceOff(); // important relPosFilter->SetIntermediateSpacing(4); diff --git a/segmentation/clitkExtractLymphStations.ggo b/segmentation/clitkExtractLymphStations.ggo index c19c607..57cf24a 100644 --- a/segmentation/clitkExtractLymphStations.ggo +++ b/segmentation/clitkExtractLymphStations.ggo @@ -26,3 +26,10 @@ option "esophagusDilatationForRight" - "Dilatation of esophagus, in mm, for 'rig option "fuzzyThresholdForS8" - "Threshold for 'Post' to dilated Esophagus" double default="0.5" no option "injectedThresholdForS8" - "Threshold injected areas in the ct" double default="150" no +section "Options for Station 7" +option "tS7_Bronchi" - "Threshold for Left/Right bronchi" double default="0.1" no +option "tS7_LeftSuperiorPulmonaryVein" - "Threshold for LeftSuperiorPulmonaryVein" double default="0.3" no +option "tS7_RightSuperiorPulmonaryVein" - "Threshold for RightSuperiorPulmonaryVein" double default="0.2" no +option "tS7_RightPulmonaryArtery" - "Threshold for RightPulmonaryArtery" double default="0.3" no +option "tS7_LeftPulmonaryArtery" - "Threshold for LeftPulmonaryArtery (NOT USEFUL YET)" double default="0.5" no +option "tS7_SVC" - "Threshold for SVC" double default="0.2" no \ No newline at end of file diff --git a/segmentation/clitkExtractLymphStationsFilter.h b/segmentation/clitkExtractLymphStationsFilter.h index 773e3bd..d67dfa2 100644 --- a/segmentation/clitkExtractLymphStationsFilter.h +++ b/segmentation/clitkExtractLymphStationsFilter.h @@ -96,11 +96,10 @@ namespace clitk { itkGetConstMacro(InjectedThresholdForS8, double); // Station 7 - itkSetMacro(FuzzyThreshold, double); - itkGetConstMacro(FuzzyThreshold, double); - itkSetMacro(Station7Filename, std::string); - itkGetConstMacro(Station7Filename, std::string); - + void SetFuzzyThresholdForS7(std::string tag, double value); + double GetFuzzyThresholdForS7(std::string tag); + + // All stations bool GetComputeStation(std::string s); void AddComputeStation(std::string station) ; @@ -121,6 +120,7 @@ namespace clitk { std::map m_ComputeStationMap; bool CheckForStation(std::string station); + void Remove_Structures(std::string s); // Station 8 double m_DistanceMaxToAnteriorPartOfTheSpine; @@ -156,19 +156,26 @@ namespace clitk { void ExtractStation_3P_Ant_Limits(); void ExtractStation_3P_Post_Limits(); void ExtractStation_3P_LR_sup_Limits(); + void ExtractStation_3P_LR_sup_Limits_2(); void ExtractStation_3P_LR_inf_Limits(); + // Station 3A + void ExtractStation_3A(); + void ExtractStation_3A_SetDefaultValues(); + void ExtractStation_3A_SI_Limits(); + void ExtractStation_3A_Ant_Limits(); + // Station 7 void ExtractStation_7(); + void ExtractStation_7_SetDefaultValues(); 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; + void ExtractStation_7_Remove_Structures(); + MaskImagePointer m_Working_Trachea; + std::map m_FuzzyThresholdForS7; MaskImagePointer m_LeftBronchus; MaskImagePointer m_RightBronchus; - MaskImagePointer m_Station7; typedef std::vector ListOfPointsType; ListOfPointsType m_RightMostInLeftBronchus; ListOfPointsType m_AntMostInLeftBronchus; @@ -205,6 +212,7 @@ namespace clitk { #include "clitkExtractLymphStationsFilter.txx" #include "clitkExtractLymphStation_8.txx" #include "clitkExtractLymphStation_3P.txx" +#include "clitkExtractLymphStation_3A.txx" #include "clitkExtractLymphStation_7.txx" #include "clitkExtractLymphStation_4RL.txx" #endif diff --git a/segmentation/clitkExtractLymphStationsFilter.txx b/segmentation/clitkExtractLymphStationsFilter.txx index ea9e2b6..aec05d7 100644 --- a/segmentation/clitkExtractLymphStationsFilter.txx +++ b/segmentation/clitkExtractLymphStationsFilter.txx @@ -55,10 +55,8 @@ ExtractLymphStationsFilter(): // Default values ExtractStation_8_SetDefaultValues(); ExtractStation_3P_SetDefaultValues(); - - // Station 7 - SetFuzzyThreshold(0.5); - SetStation7Filename("station7.mhd"); + ExtractStation_3A_SetDefaultValues(); + ExtractStation_7_SetDefaultValues(); } //-------------------------------------------------------------------- @@ -85,13 +83,19 @@ GenerateOutputInformation() { ExtractStation_3P(); StopSubStep(); - if (0) { // temporary suppress - // Extract Station7 - StartNewStep("Station 7"); - StartSubStep(); - ExtractStation_7(); - StopSubStep(); + // Extract Station3A + StartNewStep("Station 3A"); + StartSubStep(); + ExtractStation_3A(); + StopSubStep(); + + // Extract Station7 + StartNewStep("Station 7"); + StartSubStep(); + ExtractStation_7(); + StopSubStep(); + if (0) { // temporary suppress // Extract Station4RL StartNewStep("Station 4RL"); StartSubStep(); @@ -131,8 +135,6 @@ template void clitk::ExtractLymphStationsFilter:: GenerateData() { - DD("GenerateData, graft output"); - // Final Step -> graft output (if SetNthOutput => redo) this->GraftOutput(m_ListOfStations["8"]); } @@ -228,11 +230,13 @@ ExtractStation_3P() { if (CheckForStation("3P")) { ExtractStation_3P_SI_Limits(); - ExtractStation_3P_Remove_Structures(); ExtractStation_3P_Ant_Limits(); ExtractStation_3P_Post_Limits(); ExtractStation_3P_LR_sup_Limits(); + // ExtractStation_3P_LR_sup_Limits_2(); ExtractStation_3P_LR_inf_Limits(); + ExtractStation_8_Single_CCL_Limits(); // YES 8 ! + ExtractStation_3P_Remove_Structures(); // after CCL // Store image filenames into AFDB writeImage(m_ListOfStations["3P"], "seg/Station3P.mhd"); GetAFDB()->SetImageFilename("Station3P", "seg/Station3P.mhd"); @@ -247,14 +251,35 @@ template void clitk::ExtractLymphStationsFilter:: ExtractStation_7() { - if (m_ListOfStations["7"]) { - DD("Station 7 support already exist -> use it"); - m_Working_Support = m_ListOfStations["7"]; - } - else m_Working_Support = m_Mediastinum; + if (CheckForStation("7")) { ExtractStation_7_SI_Limits(); ExtractStation_7_RL_Limits(); ExtractStation_7_Posterior_Limits(); + // ExtractStation_8_Single_CCL_Limits(); // Yes the same than 8 + ExtractStation_7_Remove_Structures(); + // Store image filenames into AFDB + writeImage(m_ListOfStations["7"], "seg/Station7.mhd"); + GetAFDB()->SetImageFilename("Station7", "seg/Station7.mhd"); + WriteAFDB(); + } +} +//-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +ExtractStation_3A() +{ + if (CheckForStation("3A")) { + ExtractStation_3A_SI_Limits(); + ExtractStation_3A_Ant_Limits(); + // Store image filenames into AFDB + writeImage(m_ListOfStations["3A"], "seg/Station3A.mhd"); + GetAFDB()->SetImageFilename("Station3A", "seg/Station3A.mhd"); + WriteAFDB(); + } } //-------------------------------------------------------------------- @@ -264,6 +289,8 @@ template void clitk::ExtractLymphStationsFilter:: ExtractStation_4RL() { + DD("TODO"); + exit(0); /* WARNING ONLY 4R FIRST !!! (not same inf limits) */ @@ -275,6 +302,21 @@ ExtractStation_4RL() { //-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +Remove_Structures(std::string s) +{ + try { + StartNewStep("[Station7] Remove "+s); + MaskImagePointer Structure = GetAFDB()->template GetImage(s); + clitk::AndNot(m_Working_Support, Structure, GetBackgroundValue()); + } + catch(clitk::ExceptionObject e) { + std::cout << s << " not found, skip." << std::endl; + } +} +//-------------------------------------------------------------------- diff --git a/segmentation/clitkExtractLymphStationsGenericFilter.txx b/segmentation/clitkExtractLymphStationsGenericFilter.txx index c5c6378..0dafe7a 100644 --- a/segmentation/clitkExtractLymphStationsGenericFilter.txx +++ b/segmentation/clitkExtractLymphStationsGenericFilter.txx @@ -108,6 +108,14 @@ SetOptionsFromArgsInfoToFilter(FilterType * f) for(uint i=0; iAddComputeStation(mArgsInfo.station_arg[i]); + + // Station 7 + f->SetFuzzyThresholdForS7("Bronchi", mArgsInfo.tS7_Bronchi_arg); + f->SetFuzzyThresholdForS7("LeftSuperiorPulmonaryVein", mArgsInfo.tS7_LeftSuperiorPulmonaryVein_arg); + f->SetFuzzyThresholdForS7("RightSuperiorPulmonaryVein", mArgsInfo.tS7_RightSuperiorPulmonaryVein_arg); + f->SetFuzzyThresholdForS7("RightPulmonaryArtery", mArgsInfo.tS7_RightPulmonaryArtery_arg); + f->SetFuzzyThresholdForS7("LeftPulmonaryArtery", mArgsInfo.tS7_LeftPulmonaryArtery_arg); + f->SetFuzzyThresholdForS7("SVC", mArgsInfo.tS7_SVC_arg); } //-------------------------------------------------------------------- -- 2.45.2