From: dsarrut Date: Thu, 3 Mar 2011 10:15:11 +0000 (+0000) Subject: small improvement for S8 X-Git-Tag: v1.2.0~215 X-Git-Url: https://git.creatis.insa-lyon.fr/pubgit/?a=commitdiff_plain;h=a339bdc482ea9752ec53195bc9a47e8b05dba582;p=clitk.git small improvement for S8 --- diff --git a/segmentation/clitkAnatomicalFeatureDatabase.txx b/segmentation/clitkAnatomicalFeatureDatabase.txx index 7ec7f97..0a15c6b 100644 --- a/segmentation/clitkAnatomicalFeatureDatabase.txx +++ b/segmentation/clitkAnatomicalFeatureDatabase.txx @@ -71,7 +71,6 @@ ReleaseImage(std::string tag) } else { typename ImageType::Pointer image = GetImage(tag); - DD(image->GetReferenceCount()); image->SetReferenceCount(image->GetReferenceCount()-1); m_MapOfImage.erase(tag); } diff --git a/segmentation/clitkExtractLungFilter.txx b/segmentation/clitkExtractLungFilter.txx index fda0203..32d480c 100644 --- a/segmentation/clitkExtractLungFilter.txx +++ b/segmentation/clitkExtractLungFilter.txx @@ -392,7 +392,6 @@ GenerateOutputInformation() unsigned int initialNumberOfLabels = statisticsImageFilter->GetMaximum(); working_mask = statisticsImageFilter->GetOutput(); - DD(initialNumberOfLabels); PrintMemory(GetVerboseMemoryFlag(), "After count label"); // Decompose the first label @@ -404,7 +403,7 @@ GenerateOutputInformation() typedef clitk::DecomposeAndReconstructImageFilter DecomposeAndReconstructFilterType; typename DecomposeAndReconstructFilterType::Pointer decomposeAndReconstructFilter=DecomposeAndReconstructFilterType::New(); decomposeAndReconstructFilter->SetInput(working_mask); - decomposeAndReconstructFilter->SetVerbose(true); + decomposeAndReconstructFilter->SetVerbose(false); decomposeAndReconstructFilter->SetRadius(radius); decomposeAndReconstructFilter->SetMaximumNumberOfLabels(2); decomposeAndReconstructFilter->SetMinimumObjectSize(this->GetMinimalComponentSize()); diff --git a/segmentation/clitkExtractLymphStation_3P.txx b/segmentation/clitkExtractLymphStation_3P.txx index effcd4f..d1be710 100644 --- a/segmentation/clitkExtractLymphStation_3P.txx +++ b/segmentation/clitkExtractLymphStation_3P.txx @@ -11,6 +11,7 @@ ExtractStation_3P_SetDefaultValues() } //-------------------------------------------------------------------- + //-------------------------------------------------------------------- template void @@ -19,37 +20,452 @@ ExtractStation_3P_SI_Limits() { /* Apex of the chest & Carina. - */ + */ StartNewStep("[Station 3P] Inf/Sup limits with apex of the chest and carina"); - writeImage(m_Working_Support, "support.mhd"); - // Get Carina position (has been determined in Station8) m_CarinaZ = GetAFDB()->GetDouble("CarinaZ"); - DD(m_CarinaZ); - // Get Apex of the Chest + // Get Apex of the Chest. The "lungs" structure is autocroped so we + // consider the most superior point. MaskImagePointer Lungs = GetAFDB()->template GetImage("Lungs"); - DD("lung ok"); + MaskImageIndexType index = Lungs->GetLargestPossibleRegion().GetIndex(); + index += Lungs->GetLargestPossibleRegion().GetSize(); MaskImagePointType p; - bool found = clitk::FindExtremaPointInAGivenDirection(Lungs, - GetBackgroundValue(), - 1, true, p); - DD(found); - DD(p); + Lungs->TransformIndexToPhysicalPoint(index, p); + p[2] -= 5; // 5 mm below because the autocrop is slightly above the apex double m_ApexOfTheChest = p[2]; - DD(m_ApexOfTheChest); /* Crop support : Superior limit = carina Inferior limit = Apex of the chest */ m_Working_Support = clitk::CropImageAlongOneAxis(m_Working_Support, 2, - m_ApexOfTheChest, - m_CarinaZ, true, + m_CarinaZ, + m_ApexOfTheChest, true, GetBackgroundValue()); StopCurrentStep(m_Working_Support); m_ListOfStations["3P"] = m_Working_Support; } //-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +ExtractStation_3P_Remove_Structures() +{ + /* + 3 - (sup) remove Aorta, VB, subcl + not LR subcl ! -> a séparer LR ? + (inf) remove Eso, Aorta, Azygousvein + */ + + 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); + + StopCurrentStep(m_Working_Support); + m_ListOfStations["3P"] = m_Working_Support; +} +//-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +ExtractStation_3P_Ant_Limits() +{ + /* + Ant Post limit : + + Anteriorly, it is in contact with the posterior aspect of Stations + 1–2 superiorly (Fig. 2A–C) and with Stations 4R and 4L inferiorly + (Fig. 2D–I and 3A–C). The anterior limit of Station 3P is kept + posterior to the trachea, which is defined by an imaginary + horizontal line running along the posterior wall of the trachea + (Fig. 2B,E, red line). Posteriorly, it is delineated along the + anterior and lateral borders of the vertebral body until an + imaginary horizontal line running 1 cm posteriorly to the + anterior border of the vertebral body (Fig. 2D). + + 1 - post to the trachea : define an imaginary line just post the + Trachea and remove what is anterior. + + */ + StartNewStep("[Station 3P] Ant limits with Trachea "); + + // Load Trachea + MaskImagePointer Trachea = GetAFDB()->template GetImage("Trachea"); + + // Crop like current support (need by SliceBySliceSetBackgroundFromLineSeparation after) + Trachea = + clitk::ResizeImageLike(Trachea, m_Working_Support, GetBackgroundValue()); + + // Slice by slice, determine the most post point of the trachea (A) + // and consider a second point on the same horizontal line (B) + std::vector p_A; + std::vector p_B; + std::vector slices; + clitk::ExtractSlices(Trachea, 2, slices); + MaskImagePointType p; + typename MaskSliceType::PointType sp; + for(uint i=0; i(slices[i], 0, true, 100); + slices[i] = KeepLabels(slices[i], GetBackgroundValue(), + GetForegroundValue(), 1, 1, true); + // Find most posterior point + clitk::FindExtremaPointInAGivenDirection(slices[i], GetBackgroundValue(), + 1, false, sp); + clitk::PointsUtils::Convert2DTo3D(sp, Trachea, i, p); + p_A.push_back(p); + p[0] -= 20; + p_B.push_back(p); + } + clitk::WriteListOfLandmarks(p_A, "trachea-post.txt"); + + // Remove Ant part above those lines + clitk::SliceBySliceSetBackgroundFromLineSeparation(m_Working_Support, + p_A, p_B, + GetBackgroundValue(), + 1, 10); + // End, release memory + GetAFDB()->template ReleaseImage("Trachea"); + StopCurrentStep(m_Working_Support); + m_ListOfStations["3P"] = m_Working_Support; +} +//-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +ExtractStation_3P_Post_Limits() +{ + /* + Ant Post limit : + + Anteriorly, it is in contact with the posterior aspect of Stations + 1–2 superiorly (Fig. 2A–C) and with Stations 4R and 4L inferiorly + (Fig. 2D–I and 3A–C). The anterior limit of Station 3P is kept + posterior to the trachea, which is defined by an imaginary + horizontal line running along the posterior wall of the trachea + (Fig. 2B,E, red line). Posteriorly, it is delineated along the + anterior and lateral borders of the vertebral body until an + imaginary horizontal line running 1 cm posteriorly to the + anterior border of the vertebral body (Fig. 2D). + + 2 - post to the trachea : define an imaginary line just post the + Trachea and remove what is anterior. + + */ + StartNewStep("[Station 3P] Post limits with VertebralBody "); + + // Load VertebralBody + MaskImagePointer VertebralBody = GetAFDB()->template GetImage("VertebralBody"); + + // Crop like current support (need by SliceBySliceSetBackgroundFromLineSeparation after) + VertebralBody = clitk::ResizeImageLike(VertebralBody, m_Working_Support, GetBackgroundValue()); + + writeImage(VertebralBody, "vb.mhd"); + + // Compute VertebralBody most Ant position (again because slices + // changes). Slice by slice, determine the most post point of the + // trachea (A) and consider a second point on the same horizontal + // line (B) + std::vector p_A; + std::vector p_B; + std::vector slices; + clitk::ExtractSlices(VertebralBody, 2, slices); + MaskImagePointType p; + typename MaskSliceType::PointType sp; + for(uint i=0; i(slices[i], GetBackgroundValue(), + 1, true, sp); + + // If the VertebralBody stop superiorly before the end of + // m_Working_Support, we consider the same previous point. + if (!found) { + p = p_A.back(); + p[2] += VertebralBody->GetSpacing()[2]; + p_A.push_back(p); + p = p_B.back(); + p[2] += VertebralBody->GetSpacing()[2]; + p_B.push_back(p); + } + else { + clitk::PointsUtils::Convert2DTo3D(sp, VertebralBody, i, p); + p[1] += 10; // Consider 10 mm more post + p_A.push_back(p); + p[0] -= 20; + p_B.push_back(p); + } + } + clitk::WriteListOfLandmarks(p_A, "vb-ant.txt"); + + + // Remove Ant part above those lines + clitk::SliceBySliceSetBackgroundFromLineSeparation(m_Working_Support, + p_A, p_B, + GetBackgroundValue(), + 1, -10); + writeImage(m_Working_Support, "a.mhd"); + + StopCurrentStep(m_Working_Support); + m_ListOfStations["3P"] = m_Working_Support; +} +//-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +ExtractStation_3P_LR_sup_Limits() +{ + /* + "On the right side, the limit is defined by the air–soft-tissue + interface. On the left side, it is defined by the air–tissue + interface superiorly (Fig. 2A–C) and the aorta inferiorly + (Figs. 2D–I and 3A–C)." + + sup : + Resizelike support : Trachea, SubclavianArtery + Trachea, slice by slice, get centroid trachea + SubclavianArtery, slice by slice, CCL + prendre la 1ère à L et R, not at Left + + */ + StartNewStep("[Station 3P] Left/Right limits (superior part) "); + + // Load structures + MaskImagePointer Trachea = GetAFDB()->template GetImage("Trachea"); + MaskImagePointer SubclavianArtery = GetAFDB()->template GetImage("SubclavianArtery"); + + // Crop like current support + Trachea = clitk::ResizeImageLike(Trachea, m_Working_Support, GetBackgroundValue()); + SubclavianArtery = clitk::ResizeImageLike(SubclavianArtery, m_Working_Support, GetBackgroundValue()); + + writeImage(Trachea, "tr.mhd"); + writeImage(SubclavianArtery, "sca.mhd"); + + // Get list of slices + std::vector slices_support; + std::vector slices_trachea; + std::vector slices_subclavianartery; + clitk::ExtractSlices(m_Working_Support, 2, slices_support); + clitk::ExtractSlices(Trachea, 2, slices_trachea); + clitk::ExtractSlices(SubclavianArtery, 2, slices_subclavianartery); + + // Loop on slices + std::vector points; + MaskImagePointType p; + for(uint i=0; i centroids; + typename MaskSliceType::PointType c; + ComputeCentroids(slices_trachea[i], GetBackgroundValue(), centroids); + c = centroids[1]; + + // [debug] Store point + clitk::PointsUtils::Convert2DTo3D(centroids[1], Trachea, i, p); + points.push_back(p); + + // Get Right and Left CCL in SubclavianArtery + slices_subclavianartery[i] = Labelize(slices_subclavianartery[i], 0, true, 10); + ComputeCentroids(slices_subclavianartery[i], GetBackgroundValue(), centroids); + + if (centroids.size() > 1) { + // Determine the one at Right/Left -> first after Trachea + // centroid + typename MaskSliceType::PointType right; + typename MaskSliceType::PointType left; + int label_right=-1; + int label_left=-1; + right[0] = c[0]-100; + left[0] = c[0]+100; + for(uint j=1; j= right[0]) { + right = centroids[j]; + label_right = j; + } + } + if (centroids[j][0] > c[0]) { // At Left of Trachea centroid + if (centroids[j][0] <= left[0]) { + left = centroids[j]; + label_left = j; + } + } + } + + if (label_right != -1) { + + // Debug points + clitk::PointsUtils::Convert2DTo3D(centroids[label_right], SubclavianArtery, i, p); + points.push_back(p); + + // Set Background and ForegroundValue according to label_right + MaskSlicePointer object = + clitk::Binarize(slices_subclavianartery[i], label_right, label_right, + GetBackgroundValue(), GetForegroundValue()); + + // Relative Position : not at Right + typedef clitk::AddRelativePositionConstraintToLabelImageFilter RelPosFilterType; + typename RelPosFilterType::Pointer relPosFilter = RelPosFilterType::New(); + relPosFilter->VerboseStepFlagOff(); + relPosFilter->WriteStepFlagOff(); + relPosFilter->SetBackgroundValue(GetBackgroundValue()); + relPosFilter->SetInput(slices_support[i]); + relPosFilter->SetInputObject(object); + relPosFilter->AddOrientationTypeString("R"); + relPosFilter->SetInverseOrientationFlag(true); + // relPosFilter->SetIntermediateSpacing(3); + relPosFilter->SetIntermediateSpacingFlag(false); + relPosFilter->SetFuzzyThreshold(0.7); + relPosFilter->AutoCropFlagOff(); // important ! because we join the slices after this loop + relPosFilter->Update(); + slices_support[i] = relPosFilter->GetOutput(); + + // Relative Position : not Anterior + relPosFilter = RelPosFilterType::New(); + relPosFilter->VerboseStepFlagOff(); + relPosFilter->WriteStepFlagOff(); + relPosFilter->SetBackgroundValue(GetBackgroundValue()); + relPosFilter->SetInput(slices_support[i]); + relPosFilter->SetInputObject(object); + relPosFilter->AddOrientationTypeString("A"); + relPosFilter->SetInverseOrientationFlag(true); + // relPosFilter->SetIntermediateSpacing(3); + relPosFilter->SetIntermediateSpacingFlag(false); + relPosFilter->SetFuzzyThreshold(0.7); + relPosFilter->AutoCropFlagOff(); // important ! because we join the slices after this loop + relPosFilter->Update(); + slices_support[i] = relPosFilter->GetOutput(); + + } // End RelativePosition for Right + + + if (label_left != -1) { + + // Debug points + clitk::PointsUtils::Convert2DTo3D(centroids[label_left], SubclavianArtery, i, p); + points.push_back(p); + + // Set Background and ForegroundValue according to label_right + MaskSlicePointer object = + clitk::Binarize(slices_subclavianartery[i], label_left, label_left, + GetBackgroundValue(), GetForegroundValue()); + + // Relative Position : not at Right + typedef clitk::AddRelativePositionConstraintToLabelImageFilter RelPosFilterType; + typename RelPosFilterType::Pointer relPosFilter = RelPosFilterType::New(); + relPosFilter->VerboseStepFlagOff(); + relPosFilter->WriteStepFlagOff(); + relPosFilter->SetBackgroundValue(GetBackgroundValue()); + relPosFilter->SetInput(slices_support[i]); + relPosFilter->SetInputObject(object); + relPosFilter->AddOrientationTypeString("L"); + relPosFilter->SetInverseOrientationFlag(true); + // relPosFilter->SetIntermediateSpacing(3); + relPosFilter->SetIntermediateSpacingFlag(false); + relPosFilter->SetFuzzyThreshold(0.7); + relPosFilter->AutoCropFlagOff(); // important ! because we join the slices after this loop + relPosFilter->Update(); + slices_support[i] = relPosFilter->GetOutput(); + + // Relative Position : not Anterior + relPosFilter = RelPosFilterType::New(); + relPosFilter->VerboseStepFlagOff(); + relPosFilter->WriteStepFlagOff(); + relPosFilter->SetBackgroundValue(GetBackgroundValue()); + relPosFilter->SetInput(slices_support[i]); + relPosFilter->SetInputObject(object); + relPosFilter->AddOrientationTypeString("A"); + relPosFilter->SetInverseOrientationFlag(true); + // relPosFilter->SetIntermediateSpacing(3); + relPosFilter->SetIntermediateSpacingFlag(false); + relPosFilter->SetFuzzyThreshold(0.7); + relPosFilter->AutoCropFlagOff(); // important ! because we join the slices after this loop + relPosFilter->Update(); + slices_support[i] = relPosFilter->GetOutput(); + + } + + + } // if centroids.size > 1 + } // End loop slices + + // Joint slices + m_Working_Support = clitk::JoinSlices(slices_support, m_Working_Support, 2); + + // Save list points + clitk::WriteListOfLandmarks(points, "subcl-lr.txt"); + + + StopCurrentStep(m_Working_Support); + m_ListOfStations["3P"] = m_Working_Support; +} +//-------------------------------------------------------------------- + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +ExtractStation_3P_LR_inf_Limits() +{ + /* + "On the right side, the limit is defined by the air–soft-tissue + interface. On the left side, it is defined by the air–tissue + interface superiorly (Fig. 2A–C) and the aorta inferiorly + (Figs. 2D–I and 3A–C)." + + inf : not Right to Azygousvein + */ + StartNewStep("[Station 3P] Left/Right limits (inferior part) "); + + // Load structures + MaskImagePointer AzygousVein = GetAFDB()->template GetImage("AzygousVein"); + + typedef clitk::AddRelativePositionConstraintToLabelImageFilter RelPosFilterType; + typename RelPosFilterType::Pointer relPosFilter = RelPosFilterType::New(); + relPosFilter->VerboseStepFlagOff(); + relPosFilter->WriteStepFlagOff(); + relPosFilter->SetBackgroundValue(GetBackgroundValue()); + relPosFilter->SetInput(m_Working_Support); + relPosFilter->SetInputObject(AzygousVein); + relPosFilter->AddOrientationTypeString("R"); + relPosFilter->SetInverseOrientationFlag(true); + // relPosFilter->SetIntermediateSpacing(3); + relPosFilter->SetIntermediateSpacingFlag(false); + relPosFilter->SetFuzzyThreshold(0.7); + relPosFilter->AutoCropFlagOn(); + relPosFilter->Update(); + m_Working_Support = relPosFilter->GetOutput(); + + StopCurrentStep(m_Working_Support); + m_ListOfStations["3P"] = m_Working_Support; +} +//-------------------------------------------------------------------- diff --git a/segmentation/clitkExtractLymphStation_8.txx b/segmentation/clitkExtractLymphStation_8.txx index 8349261..dd3fdca 100644 --- a/segmentation/clitkExtractLymphStation_8.txx +++ b/segmentation/clitkExtractLymphStation_8.txx @@ -15,6 +15,7 @@ ExtractStation_8_SetDefaultValues() p[0] = 5; p[1] = 10; p[2] = 1; SetEsophagusDiltationForRight(p); SetFuzzyThresholdForS8(0.5); + SetInjectedThresholdForS8(150); } //-------------------------------------------------------------------- @@ -46,26 +47,16 @@ ExtractStation_8_SI_Limits() // Get Carina Z position MaskImagePointer Carina = GetAFDB()->template GetImage("Carina"); - // CHANGE TO -> OriginOfRightMiddleLobeBronchus ??? - std::vector centroids; clitk::ComputeCentroids(Carina, GetBackgroundValue(), centroids); m_CarinaZ = centroids[1][2]; + // DD(m_CarinaZ); // add one slice to include carina ? m_CarinaZ += m_Mediastinum->GetSpacing()[2]; // We dont need Carina structure from now Carina->Delete(); GetAFDB()->SetDouble("CarinaZ", m_CarinaZ); - // Get left lower lobe bronchus (L), take the upper border - // Get right bronchus intermedius (RML), take the lower border - - /* - double m_CarinaZ = GetAFDB()->GetPoint3D("Carina", 2) + - m_Mediastinum->GetSpacing()[2]; // add one slice to include carina - // double GOjunctionZ = GetAFDB()->GetPoint3D("GOjunction", 2); - */ - // Found most inferior part of the lung MaskImagePointer Lungs = GetAFDB()->template GetImage("Lungs"); // It should be already croped, so I took the origin and add 10mm above @@ -96,11 +87,16 @@ ExtractStation_8_SI_Limits() boolFilter->SetInput1(m_Working_Support); boolFilter->SetInput2(VertebralBody); boolFilter->SetOperationType(BoolFilterType::AndNot); - boolFilter->Update(); - boolFilter->SetInput1(boolFilter->GetOutput()); - boolFilter->SetInput2(Aorta); - boolFilter->SetOperationType(BoolFilterType::AndNot); - boolFilter->Update(); + boolFilter->Update(); + /* + + DO NOT REMOVE AORTA YET (LATER) + + boolFilter->SetInput1(boolFilter->GetOutput()); + boolFilter->SetInput2(Aorta); + boolFilter->SetOperationType(BoolFilterType::AndNot); + boolFilter->Update(); + */ m_Working_Support = boolFilter->GetOutput(); // Done. @@ -173,7 +169,7 @@ ExtractStation_8_Post_Limits() // DEBUG : write list of points clitk::WriteListOfLandmarks(vertebralAntPositions, - "vertebralMostAntPositions.txt"); + "S8-vertebralMostAntPositions-points.txt"); // Cut support posteriorly 1cm the most anterior point of the // VertebralBody. Do nothing below and above the VertebralBody. To @@ -221,7 +217,7 @@ ExtractStation_8_Post_Limits() template void clitk::ExtractLymphStationsFilter:: -ExtractStation_8_Ant_Limits() +ExtractStation_8_Ant_Sup_Limits() { //-------------------------------------------------------------------- StartNewStep("[Station8] Ant limits with S7 above Carina"); @@ -249,8 +245,8 @@ ExtractStation_8_Ant_Limits() GetBackgroundValue()); // Seprate into two main bronchi - MaskImagePointer LeftBronchus; MaskImagePointer RightBronchus; + MaskImagePointer LeftBronchus; // Labelize and consider the two first (main) labels m_Working_Trachea = Labelize(m_Working_Trachea, 0, true, 1); @@ -282,27 +278,35 @@ ExtractStation_8_Ant_Limits() 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; } + ImagePixelType rightLabel; + ImagePixelType leftLabel; + if (C1[0] < C2[0]) { rightLabel = 1; leftLabel = 2; } + else { rightLabel = 2; leftLabel = 1; } // Select LeftLabel (set one label to Backgroundvalue) - LeftBronchus = + RightBronchus = + clitk::Binarize(m_Working_Trachea, rightLabel, rightLabel, + GetBackgroundValue(), GetForegroundValue()); + /* SetBackground(m_Working_Trachea, m_Working_Trachea, - rightLabel, GetBackgroundValue(), false); - RightBronchus = + leftLabel, GetBackgroundValue(), false); + */ + LeftBronchus = clitk::Binarize(m_Working_Trachea, leftLabel, leftLabel, + GetBackgroundValue(), GetForegroundValue()); + /* SetBackground(m_Working_Trachea, m_Working_Trachea, - leftLabel, GetBackgroundValue(), false); + rightLabel, GetBackgroundValue(), false); + */ + // Crop images - LeftBronchus = clitk::AutoCrop(LeftBronchus, GetBackgroundValue()); RightBronchus = clitk::AutoCrop(RightBronchus, GetBackgroundValue()); + LeftBronchus = clitk::AutoCrop(LeftBronchus, GetBackgroundValue()); // Insert int AFDB if need after - GetAFDB()->template SetImage ("LeftBronchus", "seg/leftBronchus.mhd", - LeftBronchus, true); GetAFDB()->template SetImage ("RightBronchus", "seg/rightBronchus.mhd", RightBronchus, true); + GetAFDB()->template SetImage ("LeftBronchus", "seg/leftBronchus.mhd", + LeftBronchus, true); // Now crop below OriginOfRightMiddleLobeBronchusZ // It is not done before to keep entire bronchi. @@ -313,17 +317,17 @@ ExtractStation_8_Ant_Limits() clitk::ComputeCentroids(OriginOfRightMiddleLobeBronchus, GetBackgroundValue(), centroids); m_OriginOfRightMiddleLobeBronchusZ = centroids[1][2]; // add one slice to include carina ? - m_OriginOfRightMiddleLobeBronchusZ += LeftBronchus->GetSpacing()[2]; + m_OriginOfRightMiddleLobeBronchusZ += RightBronchus->GetSpacing()[2]; // We dont need Carina structure from now OriginOfRightMiddleLobeBronchus->Delete(); - LeftBronchus = - clitk::CropImageBelow(LeftBronchus, 2, + RightBronchus = + clitk::CropImageBelow(RightBronchus, 2, m_OriginOfRightMiddleLobeBronchusZ, true, // AutoCrop GetBackgroundValue()); - RightBronchus = - clitk::CropImageBelow(RightBronchus, 2, + LeftBronchus = + clitk::CropImageBelow(LeftBronchus, 2, m_OriginOfRightMiddleLobeBronchusZ, true, // AutoCrop GetBackgroundValue()); @@ -331,37 +335,42 @@ ExtractStation_8_Ant_Limits() // Search for points that are the most left/post/ant and most // right/post/ant of the left and right bronchus // 15 = not 15 mm more distance than the middle point. - FindExtremaPointsInBronchus(LeftBronchus, 0, 10, m_RightMostInLeftBronchus, - m_AntMostInLeftBronchus, m_PostMostInLeftBronchus); - FindExtremaPointsInBronchus(RightBronchus, 1, 10, m_LeftMostInRightBronchus, + FindExtremaPointsInBronchus(RightBronchus, 0, 10, m_LeftMostInRightBronchus, m_AntMostInRightBronchus, m_PostMostInRightBronchus); + FindExtremaPointsInBronchus(LeftBronchus, 1, 10, m_RightMostInLeftBronchus, + m_AntMostInLeftBronchus, m_PostMostInLeftBronchus); + // DEBUG : write the list of points ListOfPointsType v; + v.reserve(m_LeftMostInRightBronchus.size()+m_AntMostInRightBronchus.size()+ + m_PostMostInRightBronchus.size()); + v.insert(v.end(), m_LeftMostInRightBronchus.begin(), m_LeftMostInRightBronchus.end() ); + v.insert(v.end(), m_AntMostInRightBronchus.begin(), m_AntMostInRightBronchus.end() ); + v.insert(v.end(), m_PostMostInRightBronchus.begin(), m_PostMostInRightBronchus.end() ); + clitk::WriteListOfLandmarks(v, "S8-RightBronchus-points.txt"); + + v.clear(); v.reserve(m_RightMostInLeftBronchus.size()+m_AntMostInLeftBronchus.size()+ m_PostMostInLeftBronchus.size()); v.insert(v.end(), m_RightMostInLeftBronchus.begin(), m_RightMostInLeftBronchus.end() ); v.insert(v.end(), m_AntMostInLeftBronchus.begin(), m_AntMostInLeftBronchus.end() ); v.insert(v.end(), m_PostMostInLeftBronchus.begin(), m_PostMostInLeftBronchus.end() ); - clitk::WriteListOfLandmarks(v, "LeftBronchusPoints.txt"); + clitk::WriteListOfLandmarks(v, "S8-LeftBronchus-points.txt"); v.clear(); - v.reserve(m_LeftMostInRightBronchus.size()+m_AntMostInRightBronchus.size()+ - m_PostMostInRightBronchus.size()); - v.insert(v.end(), m_LeftMostInRightBronchus.begin(), m_LeftMostInRightBronchus.end() ); - v.insert(v.end(), m_AntMostInRightBronchus.begin(), m_AntMostInRightBronchus.end() ); + v.reserve(m_PostMostInLeftBronchus.size()+m_PostMostInRightBronchus.size()); + v.insert(v.end(), m_PostMostInLeftBronchus.begin(), m_PostMostInLeftBronchus.end() ); v.insert(v.end(), m_PostMostInRightBronchus.begin(), m_PostMostInRightBronchus.end() ); - clitk::WriteListOfLandmarks(v, "RightBronchusPoints.txt"); + clitk::WriteListOfLandmarks(v, "S8-RightLeftBronchus-points.txt"); // Now uses these points to limit, slice by slice // line is mainly horizontal, so mainDirection=1 - writeImage(m_Working_Support, "before.mhd"); clitk::SliceBySliceSetBackgroundFromLineSeparation(m_Working_Support, - m_PostMostInLeftBronchus, m_PostMostInRightBronchus, + m_PostMostInLeftBronchus, GetBackgroundValue(), 1, 10); - writeImage(m_Working_Support, "after.mhd"); // Keep main 3D CCL : m_Working_Support = Labelize(m_Working_Support, 0, false, 10); @@ -382,11 +391,11 @@ ExtractStation_8_Ant_Limits() template void clitk::ExtractLymphStationsFilter:: -ExtractStation_8_LR_Limits() +ExtractStation_8_Ant_Inf_Limits() { //-------------------------------------------------------------------- - StartNewStep("[Station8] Ant part (not post to Esophagus)"); + StartNewStep("[Station8] Ant part: not post to Esophagus"); /* Consider Esophagus, dilate it and remove ant part. It remains part on L & R, than can be partly removed by cutting what remains at @@ -394,7 +403,7 @@ ExtractStation_8_LR_Limits() */ // Get Esophagus - MaskImagePointer Esophagus = GetAFDB()->template GetImage("Esophagus"); + m_Esophagus = GetAFDB()->template GetImage("Esophagus"); clitk::PrintMemory(GetVerboseMemoryFlag(), "after read Esophagus"); // In images from the original article, Atlas – UM, the oesophagus @@ -406,30 +415,36 @@ ExtractStation_8_LR_Limits() //Dr. Aletta Frasier exclude this structure. From an oncological //prospective, the oesophagus should be excluded from these nodal //stations. - typedef clitk::BooleanOperatorLabelImageFilter BoolFilterType; - typename BoolFilterType::Pointer boolFilter = BoolFilterType::New(); - boolFilter->InPlaceOn(); - boolFilter->SetInput1(m_Working_Support); - boolFilter->SetInput2(Esophagus); - boolFilter->SetOperationType(BoolFilterType::AndNot); - boolFilter->Update(); - m_Working_Support = boolFilter->GetOutput(); + + /* NOT YET !! DO IT LATER + + typedef clitk::BooleanOperatorLabelImageFilter BoolFilterType; + typename BoolFilterType::Pointer boolFilter = BoolFilterType::New(); + boolFilter->InPlaceOn(); + boolFilter->SetInput1(m_Working_Support); + boolFilter->SetInput2(m_Esophagus); + boolFilter->SetOperationType(BoolFilterType::AndNot); + boolFilter->Update(); + m_Working_Support = boolFilter->GetOutput(); + + */ // Crop Esophagus : keep only below the OriginOfRightMiddleLobeBronchusZ - Esophagus = - clitk::CropImageAbove(Esophagus, 2, + m_Esophagus = + clitk::CropImageAbove(m_Esophagus, 2, m_OriginOfRightMiddleLobeBronchusZ, true, // AutoCrop GetBackgroundValue()); // Dilate to keep only not Anterior positions MaskImagePointType radiusInMM = GetEsophagusDiltationForAnt(); - Esophagus = EnlargeEsophagusDilatationRadiusInferiorly(Esophagus); - Esophagus = clitk::Dilate(Esophagus, - radiusInMM, - GetBackgroundValue(), - GetForegroundValue(), true); - writeImage(Esophagus, "enlarged-eso.mhd"); + + // m_Esophagus = EnlargeEsophagusDilatationRadiusInferiorly(m_Esophagus); + + m_Esophagus = clitk::Dilate(m_Esophagus, + radiusInMM, + GetBackgroundValue(), + GetForegroundValue(), true); // Remove Anterior part according to this dilatated esophagus. Note: // because we crop Esophagus with ORML, the support will also be @@ -444,16 +459,19 @@ ExtractStation_8_LR_Limits() typedef clitk::SliceBySliceRelativePositionFilter RelPosFilterType; typename RelPosFilterType::Pointer relPosFilter = RelPosFilterType::New(); relPosFilter->VerboseStepFlagOff(); + relPosFilter->WriteStepFlagOff(); relPosFilter->SetInput(m_Working_Support); - relPosFilter->SetInputObject(Esophagus); + relPosFilter->SetInputObject(m_Esophagus); relPosFilter->AddOrientationTypeString("P"); relPosFilter->InverseOrientationFlagOff(); relPosFilter->SetDirection(2); // Z axis relPosFilter->UniqueConnectedComponentBySliceOff(); relPosFilter->SetIntermediateSpacing(3); - relPosFilter->ResampleBeforeRelativePositionFilterOn(); + relPosFilter->IntermediateSpacingFlagOn(); relPosFilter->SetFuzzyThreshold(GetFuzzyThresholdForS8()); relPosFilter->RemoveObjectFlagOff(); // Do not exclude here because it is dilated + relPosFilter->CombineWithOrFlagOff(); // NO ! + relPosFilter->IgnoreEmptySliceObjectFlagOn(); relPosFilter->Update(); m_Working_Support = relPosFilter->GetOutput(); @@ -461,6 +479,216 @@ ExtractStation_8_LR_Limits() m_Working_Support = clitk::AutoCrop(m_Working_Support, GetBackgroundValue()); StopCurrentStep(m_Working_Support); + m_ListOfStations["8"] = m_Working_Support; +} +//-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +ExtractStation_8_Ant_Injected_Limits() +{ + + //-------------------------------------------------------------------- + StartNewStep("[Station8] Ant part (remove high density, injected part)"); + + // Consider initial image, crop to current support + ImagePointer working_input = + clitk::ResizeImageLike(m_Input, + m_Working_Support, + (short)GetBackgroundValue()); + + // Threshold + typedef itk::BinaryThresholdImageFilter BinarizeFilterType; + typename BinarizeFilterType::Pointer binarizeFilter = BinarizeFilterType::New(); + binarizeFilter->SetInput(working_input); + binarizeFilter->SetLowerThreshold(GetInjectedThresholdForS8()); + binarizeFilter->SetInsideValue(GetForegroundValue()); + binarizeFilter->SetOutsideValue(GetBackgroundValue()); + binarizeFilter->Update(); + MaskImagePointer injected = binarizeFilter->GetOutput(); + + // Combine with current support + clitk::AndNot(m_Working_Support, injected, GetBackgroundValue()); + + StopCurrentStep(m_Working_Support); + m_ListOfStations["8"] = m_Working_Support; +} +//-------------------------------------------------------------------- + + + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +ExtractStation_8_LR_1_Limits() +{ + //-------------------------------------------------------------------- + StartNewStep("[Station8] Left and Right (from Carina to PulmonaryTrunk): Right to LeftPulmonaryArtery"); + + /* + We remove LeftPulmonaryArtery structure and what is at Left to + this structure. + */ + MaskImagePointer LeftPulmonaryArtery = GetAFDB()->template GetImage("LeftPulmonaryArtery"); + + // Relative Position : not at Left + typedef clitk::AddRelativePositionConstraintToLabelImageFilter RelPosFilterType; + typename RelPosFilterType::Pointer relPosFilter = RelPosFilterType::New(); + relPosFilter->VerboseStepFlagOff(); + relPosFilter->WriteStepFlagOff(); + relPosFilter->SetBackgroundValue(GetBackgroundValue()); + relPosFilter->SetInput(m_Working_Support); + relPosFilter->SetInputObject(LeftPulmonaryArtery); + relPosFilter->RemoveObjectFlagOn(); // remove the object too + relPosFilter->AddOrientationTypeString("L"); + relPosFilter->InverseOrientationFlagOn(); // Not at Left + relPosFilter->SetIntermediateSpacing(3); + relPosFilter->IntermediateSpacingFlagOn(); + relPosFilter->SetFuzzyThreshold(0.7); + relPosFilter->AutoCropFlagOn(); + relPosFilter->Update(); + m_Working_Support = relPosFilter->GetOutput(); + + // Release LeftPulmonaryArtery + GetAFDB()->template ReleaseImage("LeftPulmonaryArtery"); + + StopCurrentStep(m_Working_Support); + m_ListOfStations["8"] = m_Working_Support; +} +//-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +ExtractStation_8_LR_2_Limits() +{ + //-------------------------------------------------------------------- + StartNewStep("[Station8] Left and Right (from PulmTrunk to OriginMiddleLobeBronchus) Right to line from Aorta to PulmonaryTrunk"); + + /* + We consider a line from Left part of Aorta to left part of + PulmonaryTrunk and remove what is at Left. + */ + + // Load the structures + MaskImagePointer Aorta = GetAFDB()->template GetImage("Aorta"); + MaskImagePointer PulmonaryTrunk = GetAFDB()->template GetImage("PulmonaryTrunk"); + + // Resize like the PT and define the slices + MaskImagePointType min, max; + clitk::GetMinMaxPointPosition(PulmonaryTrunk, min, max); + Aorta = clitk::CropImageAlongOneAxis(Aorta, 2, min[2], max[2], false, GetBackgroundValue()); + std::vector slices_aorta; + clitk::ExtractSlices(Aorta, 2, slices_aorta); + std::vector slices_PT; + clitk::ExtractSlices(PulmonaryTrunk, 2, slices_PT); + + // Find the points at left + std::vector p_A; + std::vector p_B; + MaskImagePointType pA; + MaskImagePointType pB; + for(uint i=0; i(slices_aorta[i], GetBackgroundValue(), 0, false, ps); + clitk::PointsUtils::Convert2DTo3D(ps, Aorta, i, pA); + + if (found) { + // In PT : generally 2 CCL, we keep the most at left + found = + clitk::FindExtremaPointInAGivenDirection(slices_PT[i], GetBackgroundValue(), 0, false, ps); + clitk::PointsUtils::Convert2DTo3D(ps, PulmonaryTrunk, i, pB); + } + + if (found) { + p_A.push_back(pA); + p_B.push_back(pB); + } + } + clitk::WriteListOfLandmarks(p_A, "S8-Aorta-Left-points.txt"); + clitk::WriteListOfLandmarks(p_B, "S8-PT-Left-points.txt"); + + // Remove part at Left + clitk::SliceBySliceSetBackgroundFromLineSeparation(m_Working_Support, + p_A, p_B, + GetBackgroundValue(), + 0, -10); + + StopCurrentStep(m_Working_Support); + m_ListOfStations["8"] = m_Working_Support; +} +//-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +ExtractStation_8_Single_CCL_Limits() +{ + //-------------------------------------------------------------------- + StartNewStep("[Station8] Single CCL"); + + // Consider slices + std::vector slices; + std::vector slices_vb; + clitk::ExtractSlices(m_Working_Support, 2, slices); + MaskImagePointer VertebralBody = + GetAFDB()->template GetImage ("VertebralBody"); + VertebralBody = clitk::ResizeImageLike(VertebralBody, m_Working_Support, GetBackgroundValue()); + clitk::ExtractSlices(VertebralBody, 2, slices_vb); + + for(uint i=0; i(slices[i], 0, true, 100); + // Compute centroids coordinate + std::vector centroids; + 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(slices[i], GetBackgroundValue(), + GetForegroundValue(), index, index, true); + } + } + m_Working_Support = clitk::JoinSlices(slices, m_Working_Support, 2); + + StopCurrentStep(m_Working_Support); + m_ListOfStations["8"] = m_Working_Support; +} +//-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +ExtractStation_8_LR_Limits_old2() +{ + //-------------------------------------------------------------------- StartNewStep("[Station8] Left and Right limits arround esophagus (below Carina)"); @@ -468,8 +696,22 @@ ExtractStation_8_LR_Limits() std::vector slices; clitk::ExtractSlices(m_Working_Support, 2, slices); + // Dilate the Esophagus to consider a margins around + MaskImagePointType radiusInMM = GetEsophagusDiltationForAnt(); + m_Esophagus = clitk::Dilate(m_Esophagus, + radiusInMM, + GetBackgroundValue(), + GetForegroundValue(), true); + + // Remove what is outside the mediastinum in this enlarged Esophagus -> it allows to select + // 'better' extrema points (not too post). + MaskImagePointer Lungs = GetAFDB()->template GetImage("Lungs"); + clitk::AndNot(m_Esophagus, Lungs, GetBackgroundValue()); + GetAFDB()->template ReleaseImage("Lungs"); + // Estract slices of Esophagus (resize like support before to have the same set of slices) - MaskImagePointer EsophagusForSlice = clitk::ResizeImageLike(Esophagus, m_Working_Support, GetBackgroundValue()); + MaskImagePointer EsophagusForSlice = clitk::ResizeImageLike(m_Esophagus, m_Working_Support, GetBackgroundValue()); + std::vector eso_slices; clitk::ExtractSlices(EsophagusForSlice, 2, eso_slices); @@ -497,6 +739,8 @@ ExtractStation_8_LR_Limits() std::vector p_LeftMostAnt; std::vector p_LeftMostPost; std::vector p_AllPoints; + std::vector p_LeftAorta; + std::vector p_LeftEso; /* In the following, we search for the LeftRight limits. We search @@ -507,7 +751,7 @@ ExtractStation_8_LR_Limits() is concave between Eso/VertebralBody. Esophagus is a bit dilatated. On VertebralBody we go right (or left) until we reach the lung (but no more 20 mm). - */ + */ // Loop slices MaskImagePointType p; @@ -603,9 +847,26 @@ ExtractStation_8_LR_Limits() p_LeftMostPost.push_back(p); p_AllPoints.push_back(p); + // Find Eso slice centroid and do not consider what is post to + // this centroid. + std::vector 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] = + clitk::ResizeImageLike(eso_slices[i], aorta_slices[i], GetBackgroundValue()); + // writeImage(eso_slices[i], "eso-slice-"+toString(i)+".mhd"); + } + // Find right limit of Esophagus and Aorta - clitk::FindExtremaPointInAGivenDirection(eso_slices[i], GetBackgroundValue(), 0, true, sp_maxRight_Eso); - bool f = clitk::FindExtremaPointInAGivenDirection(aorta_slices[i], GetBackgroundValue(), 0, true, sp_maxRight_Aorta); + bool f = + clitk::FindExtremaPointInAGivenDirection(eso_slices[i], GetBackgroundValue(), + 0, true, sp_maxRight_Eso); + f = f && + clitk::FindExtremaPointInAGivenDirection(aorta_slices[i], GetBackgroundValue(), + 0, true, sp_maxRight_Aorta); clitk::PointsUtils::Convert2DTo3D(sp_maxRight_Eso, EsophagusForSlice, i, p); clitk::PointsUtils::Convert2DTo3D(sp_maxRight_Aorta, Aorta, i, pp); pp[0] -= 2; // Add a margin of 2 mm to include the Aorta 'wall' @@ -619,13 +880,14 @@ ExtractStation_8_LR_Limits() if (s>0) p_RightMostAnt.push_back(p); else p_RightMostAnt.push_back(pp); } - else { - p_RightMostAnt.push_back(p); + else { // No more Esophagus in this slice : do nothing + // p_RightMostAnt.push_back(p); + p_RightMostPost.pop_back(); } // -------------------------------------------------------------------------- // Find the limit on the Left: most left point between Eso and - // Vertebra. (Left is left on screen, coordinate increase) + // Eso. (Left is left on screen, coordinate increase) // Find left limit of Esophagus clitk::FindExtremaPointInAGivenDirection(eso_slices[i], GetBackgroundValue(), 0, false, sp_maxLeft_Eso); @@ -640,15 +902,26 @@ ExtractStation_8_LR_Limits() MaskImagePointType B = p; MaskImagePointType C = pp; double s = (B[0] - A[0]) * (C[1] - A[1]) - (B[1] - A[1]) * (C[0] - A[0]); - if (s<0) p_LeftMostAnt.push_back(p); // Insert point most at Left - else p_LeftMostAnt.push_back(pp); + if (s<0) { + p_LeftMostAnt.push_back(p); // Insert point most at Left, Eso + } + else { + // in this case -> two lines ! + p_LeftMostAnt.push_back(pp); // Insert point most at Left, Aorta (Vert to Aorta) + // but also consider Aorta to Eso + p_LeftAorta.push_back(pp); + p_LeftEso.push_back(p); + } } - else { - p_LeftMostAnt.push_back(p); + else { // No more Esophagus in this slice : do nothing + p_LeftMostPost.pop_back(); + //p_LeftMostAnt.push_back(p); } } // End of slice loop - clitk::WriteListOfLandmarks(p_AllPoints, "LR-Eso-Vert.txt"); + clitk::WriteListOfLandmarks(p_AllPoints, "S8-LR-Eso-Vert.txt"); + clitk::WriteListOfLandmarks(p_LeftEso, "S8-Left-Eso.txt"); + clitk::WriteListOfLandmarks(p_LeftAorta, "S8-Left-Aorta.txt"); // Now uses these points to limit, slice by slice // Line is mainly vertical, so mainDirection=0 @@ -658,6 +931,222 @@ ExtractStation_8_LR_Limits() clitk::SliceBySliceSetBackgroundFromLineSeparation(m_Working_Support, p_LeftMostAnt, p_LeftMostPost, GetBackgroundValue(), 0, -10); + clitk::SliceBySliceSetBackgroundFromLineSeparation(m_Working_Support, + p_LeftEso,p_LeftAorta, + GetBackgroundValue(), 0, -10); + // END + StopCurrentStep(m_Working_Support); + m_ListOfStations["8"] = m_Working_Support; + return; +} +//-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +ExtractStation_8_LR_Limits() +{ + + //-------------------------------------------------------------------- + StartNewStep("[Station8] Left and Right limits arround esophagus with lines from VertebralBody-Aorta-Esophagus"); + + // Estract slices for current support for slice by slice processing + std::vector slices; + clitk::ExtractSlices(m_Working_Support, 2, slices); + + // Dilate the Esophagus to consider a margins around + MaskImagePointType radiusInMM = GetEsophagusDiltationForAnt(); + m_Esophagus = clitk::Dilate(m_Esophagus, radiusInMM, + GetBackgroundValue(), GetForegroundValue(), true); + // m_Esophagus = EnlargeEsophagusDilatationRadiusInferiorly(m_Esophagus); + + // Remove what is outside the mediastinum in this enlarged Esophagus -> it allows to select + // 'better' extrema points (not too post). + MaskImagePointer Lungs = GetAFDB()->template GetImage("Lungs"); + clitk::AndNot(m_Esophagus, Lungs, GetBackgroundValue()); + GetAFDB()->template ReleaseImage("Lungs"); + + // Estract slices of Esophagus (resize like support before to have the same set of slices) + MaskImagePointer EsophagusForSlice = clitk::ResizeImageLike(m_Esophagus, m_Working_Support, GetBackgroundValue()); + std::vector eso_slices; + clitk::ExtractSlices(EsophagusForSlice, 2, eso_slices); + + // Estract slices of Vertebral (resize like support before to have the same set of slices) + MaskImagePointer VertebralBody = GetAFDB()->template GetImage("VertebralBody"); + VertebralBody = clitk::ResizeImageLike(VertebralBody, m_Working_Support, GetBackgroundValue()); + std::vector vert_slices; + clitk::ExtractSlices(VertebralBody, 2, vert_slices); + + // Estract slices of Aorta (resize like support before to have the same set of slices) + MaskImagePointer Aorta = GetAFDB()->template GetImage("Aorta"); + Aorta = clitk::ResizeImageLike(Aorta, m_Working_Support, GetBackgroundValue()); + std::vector aorta_slices; + clitk::ExtractSlices(Aorta, 2, aorta_slices); + + // Extract slices of Mediastinum (resize like support before to have the same set of slices) + m_Mediastinum = GetAFDB()->template GetImage("Mediastinum"); + m_Mediastinum = clitk::ResizeImageLike(m_Mediastinum, m_Working_Support, GetBackgroundValue()); + std::vector mediast_slices; + clitk::ExtractSlices(m_Mediastinum, 2, mediast_slices); + + // List of points + std::vector p_MostLeftVertebralBody; + std::vector p_MostRightVertebralBody; + std::vector p_MostLeftAorta; + std::vector p_MostLeftEsophagus; + + /* + In the following, we search for the LeftRight limits. We search + for the most Right points in Esophagus and in VertebralBody and + consider a line between those to most right points. All points in + the support which are most right to this line are discarded. Same + for the left part. The underlying assumption is that the support + is concave between Eso/VertebralBody. Esophagus is a bit + dilatated. On VertebralBody we go right (or left) until we reach + the lung (but no more 20 mm). + */ + + // Temporary 3D point + MaskImagePointType p; + + // Loop slices + for(uint i=0; i(vert_slices[j], GetBackgroundValue(), 1, true, sp_MostAntVertebralBody); + if (!found) { + // It should not happen ! But sometimes, a contour is missing or + // the VertebralBody is not delineated enough inferiorly ... in + // those cases, we consider the first found slice. + std::cerr << "No foreground pixels in this VertebralBody slices !?? I try with the previous/next slice" << std::endl; + j++; + } + } + // DD(j); + sp_MostAntVertebralBody[1] += GetDistanceMaxToAnteriorPartOfTheSpine(); // Consider offset + + // Crop the vertebralbody below this most post line + vert_slices[j] = + clitk::CropImageAbove(vert_slices[j], 1, sp_MostAntVertebralBody[1], false, GetBackgroundValue()); + vert_slices[j] = + clitk::ResizeImageLike(vert_slices[j], aorta_slices[i], GetBackgroundValue()); + // writeImage(vert_slices[i], "vert-slice-"+toString(i)+".mhd"); + + // Find first point not in mediastinum + clitk::FindExtremaPointInAGivenDirection(vert_slices[j], GetBackgroundValue(), 0, false, sp_MostLeftVertebralBody); + sp_MostLeftVertebralBody = clitk::FindExtremaPointInAGivenLine(mediast_slices[i], 0, false, sp_MostLeftVertebralBody, GetBackgroundValue(), 30); + sp_MostLeftVertebralBody[0] += 2; // 2mm margin + clitk::FindExtremaPointInAGivenDirection(vert_slices[j], GetBackgroundValue(), 0, true, sp_MostRightVertebralBody); + sp_MostRightVertebralBody = clitk::FindExtremaPointInAGivenLine(mediast_slices[i], 0, true, sp_MostRightVertebralBody, GetBackgroundValue(),30); + sp_MostRightVertebralBody[0] -= 2; // 2 mm margin + + // Convert 2D points into 3D + clitk::PointsUtils::Convert2DTo3D(sp_MostRightVertebralBody, VertebralBody, i, p); + p_MostRightVertebralBody.push_back(p); + clitk::PointsUtils::Convert2DTo3D(sp_MostLeftVertebralBody, VertebralBody, i, p); + p_MostLeftVertebralBody.push_back(p); + + // DD(p_MostLeftVertebralBody.back()); + + /* ------------------------------------------------------------------------------------- + Find most Left point in Esophagus + */ + found = clitk::FindExtremaPointInAGivenDirection(eso_slices[i], GetBackgroundValue(), 0, false, sp_MostLeftEsophagus); + if (!found) { // No more Esophagus, I remove the previous point + + // if (p_MostLeftEsophagus.size() < 1) { + p_MostLeftVertebralBody.pop_back(); + // } + // else { + // // Consider the previous point + // p = p_MostLeftEsophagus.back(); + // 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()); + } + + /* ------------------------------------------------------------------------------------- + Find most Left point in Aorta + */ + if (found) { + clitk::FindExtremaPointInAGivenDirection(aorta_slices[i], GetBackgroundValue(), 0, false, sp_MostLeftAorta); + sp_MostLeftAorta = clitk::FindExtremaPointInAGivenLine(mediast_slices[i], 0, false, sp_MostLeftAorta, GetBackgroundValue(), 10); + typename MaskSliceType::PointType temp=sp_MostLeftEsophagus; + temp[0] -= 10; + if (clitk::IsOnTheSameLineSide(sp_MostLeftAorta, sp_MostLeftVertebralBody, sp_MostLeftEsophagus, temp)) { + // sp_MostLeftAorta is on the same side than temp (at Right) -> ignore Aorta + sp_MostLeftAorta = sp_MostLeftEsophagus; + } + clitk::PointsUtils::Convert2DTo3D(sp_MostLeftAorta, Aorta, i, p); + p_MostLeftAorta.push_back(p); + } + + } // End of slice loop + + clitk::WriteListOfLandmarks(p_MostLeftVertebralBody, "S8-MostLeft-VB-points.txt"); + clitk::WriteListOfLandmarks(p_MostRightVertebralBody, "S8-MostRight-VB-points.txt"); + clitk::WriteListOfLandmarks(p_MostLeftAorta, "S8-MostLeft-Aorta-points.txt"); + clitk::WriteListOfLandmarks(p_MostLeftEsophagus, "S8-MostLeft-eso-points.txt"); + + clitk::SliceBySliceSetBackgroundFromLineSeparation(m_Working_Support, + p_MostLeftVertebralBody, p_MostLeftAorta, + GetBackgroundValue(), 0, -10); + + clitk::SliceBySliceSetBackgroundFromLineSeparation(m_Working_Support, + p_MostLeftAorta, p_MostLeftEsophagus, + GetBackgroundValue(), 0, -10); + + // END + StopCurrentStep(m_Working_Support); + m_ListOfStations["8"] = m_Working_Support; + return; +} +//-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +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()); + // END StopCurrentStep(m_Working_Support); m_ListOfStations["8"] = m_Working_Support; @@ -672,12 +1161,25 @@ typename clitk::ExtractLymphStationsFilter::MaskImagePointer clitk::ExtractLymphStationsFilter:: EnlargeEsophagusDilatationRadiusInferiorly(MaskImagePointer & Esophagus) { - // Check if Esophagus is delineated at least until GOjunction. Now, + // Check if Esophagus is delineated at least until Diaphragm. Now, // because we use AutoCrop, Origin[2] gives this max inferior // position. - // double GOjunctionZ = GetAFDB()->GetPoint3D("GOjunction", 2); + + DD("BUGGY DONT USE"); + exit(0); if (Esophagus->GetOrigin()[2] > m_DiaphragmInferiorLimit) { + // crop first slice without mask + MaskImagePointType pt; + clitk::FindExtremaPointInAGivenDirection(Esophagus, GetBackgroundValue(), 2, true, pt); + DD(pt); + Esophagus = + clitk::CropImageBelow(Esophagus, 2, + pt[2], + false, // AutoCrop + GetBackgroundValue()); + writeImage(Esophagus, "crop-eso.mhd"); + std::cout << "Warning Esophagus is not delineated until Diaphragm. I mirror-pad it." << std::endl; double extraSize = Esophagus->GetOrigin()[2]-m_DiaphragmInferiorLimit; @@ -746,7 +1248,7 @@ ExtractStation_8_LR_Limits_old() relPosFilter->SetDirection(2); // Z axis relPosFilter->UniqueConnectedComponentBySliceOff(); // important relPosFilter->SetIntermediateSpacing(4); - relPosFilter->ResampleBeforeRelativePositionFilterOn(); + relPosFilter->IntermediateSpacingFlagOn(); relPosFilter->SetFuzzyThreshold(0.9); // remove few part only relPosFilter->RemoveObjectFlagOff(); relPosFilter->Update(); @@ -776,7 +1278,7 @@ ExtractStation_8_LR_Limits_old() relPosFilter->SetDirection(2); // Z axis relPosFilter->UniqueConnectedComponentBySliceOff(); // important relPosFilter->SetIntermediateSpacing(4); - relPosFilter->ResampleBeforeRelativePositionFilterOn(); + relPosFilter->IntermediateSpacingFlagOn(); relPosFilter->CombineWithOrFlagOn(); relPosFilter->SetFuzzyThreshold(0.9); // remove few part only relPosFilter->RemoveObjectFlagOff(); @@ -788,3 +1290,73 @@ ExtractStation_8_LR_Limits_old() m_ListOfStations["8"] = m_Working_Support; } //-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +FindExtremaPointsInBronchus(MaskImagePointer input, + int direction, + double distance_max_from_center_point, + ListOfPointsType & LR, + ListOfPointsType & Ant, + ListOfPointsType & Post) +{ + + // Other solution ==> with auto bounding box ! (but pb to prevent to + // be too distant from the center point + + // Extract slices + std::vector slices; + clitk::ExtractSlices(input, 2, slices); + + // Loop on slices + bool found; + for(uint i=0; i(slices[i], 0, true, 10); + slices[i] = KeepLabels(slices[i], + GetBackgroundValue(), + GetForegroundValue(), 1, 1, true); + */ + + // ------- Find rightmost or leftmost point ------- + MaskSliceType::PointType LRMost; + found = + clitk::FindExtremaPointInAGivenDirection(slices[i], + GetBackgroundValue(), + 0, // axis XY + (direction==0?false:true), // right or left according to direction + LRMost); + // ------- Find postmost point ------- + MaskSliceType::PointType postMost; + found = + clitk::FindExtremaPointInAGivenDirection(slices[i], + GetBackgroundValue(), + 1, false, LRMost, + distance_max_from_center_point, + postMost); + // ------- Find antmost point ------- + MaskSliceType::PointType antMost; + found = + clitk::FindExtremaPointInAGivenDirection(slices[i], + GetBackgroundValue(), + 1, true, LRMost, + distance_max_from_center_point, + antMost); + // Only add point if found + if (found) { + // ------- Convert 2D to 3D points -------- + MaskImageType::PointType p; + clitk::PointsUtils::Convert2DTo3D(LRMost, input, i, p); + LR.push_back(p); + clitk::PointsUtils::Convert2DTo3D(antMost, input, i, p); + Ant.push_back(p); + clitk::PointsUtils::Convert2DTo3D(postMost, input, i, p); + Post.push_back(p); + } + } +} +//-------------------------------------------------------------------- diff --git a/segmentation/clitkExtractLymphStations.ggo b/segmentation/clitkExtractLymphStations.ggo index 4b39015..c19c607 100644 --- a/segmentation/clitkExtractLymphStations.ggo +++ b/segmentation/clitkExtractLymphStations.ggo @@ -24,5 +24,5 @@ option "maxAntSpine" - "Distance max to anterior part of the spine in mm" doubl option "esophagusDilatationForAnt" - "Dilatation of esophagus, in mm, for 'anterior' limits (default=15,2,1)" double no multiple option "esophagusDilatationForRight" - "Dilatation of esophagus, in mm, for 'right' limits (default=5,10,1)" double no multiple 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 diff --git a/segmentation/clitkExtractLymphStationsFilter.h b/segmentation/clitkExtractLymphStationsFilter.h index 8c7b329..773e3bd 100644 --- a/segmentation/clitkExtractLymphStationsFilter.h +++ b/segmentation/clitkExtractLymphStationsFilter.h @@ -71,6 +71,7 @@ namespace clitk { typedef typename MaskImageType::PointType MaskImagePointType; typedef itk::Image MaskSliceType; + typedef typename MaskSliceType::Pointer MaskSlicePointer; /** ImageDimension constants */ itkStaticConstMacro(ImageDimension, unsigned int, ImageType::ImageDimension); @@ -91,6 +92,9 @@ namespace clitk { itkSetMacro(FuzzyThresholdForS8, double); itkGetConstMacro(FuzzyThresholdForS8, double); + itkSetMacro(InjectedThresholdForS8, double); + itkGetConstMacro(InjectedThresholdForS8, double); + // Station 7 itkSetMacro(FuzzyThreshold, double); itkGetConstMacro(FuzzyThreshold, double); @@ -124,6 +128,8 @@ namespace clitk { double m_CarinaZ; double m_OriginOfRightMiddleLobeBronchusZ; double m_FuzzyThresholdForS8; + double m_InjectedThresholdForS8; + MaskImagePointer m_Esophagus; MaskImagePointType m_EsophagusDiltationForAnt; MaskImagePointType m_EsophagusDiltationForRight; MaskImagePointer EnlargeEsophagusDilatationRadiusInferiorly(MaskImagePointer & eso); @@ -131,15 +137,27 @@ namespace clitk { void ExtractStation_8_SetDefaultValues(); void ExtractStation_8_SI_Limits(); void ExtractStation_8_Post_Limits(); - void ExtractStation_8_Ant_Limits(); + void ExtractStation_8_Ant_Sup_Limits(); + void ExtractStation_8_Ant_Inf_Limits(); + void ExtractStation_8_Ant_Injected_Limits(); + void ExtractStation_8_LR_1_Limits(); + void ExtractStation_8_LR_2_Limits(); + void ExtractStation_8_Single_CCL_Limits(); void ExtractStation_8_LR_Limits(); + void ExtractStation_8_Remove_Structures(); void ExtractStation_8_LR_Limits_old(); + void ExtractStation_8_LR_Limits_old2(); // Station 3P void ExtractStation_3P(); void ExtractStation_3P_SetDefaultValues(); void ExtractStation_3P_SI_Limits(); - + void ExtractStation_3P_Remove_Structures(); + void ExtractStation_3P_Ant_Limits(); + void ExtractStation_3P_Post_Limits(); + void ExtractStation_3P_LR_sup_Limits(); + void ExtractStation_3P_LR_inf_Limits(); + // Station 7 void ExtractStation_7(); void ExtractStation_7_SI_Limits(); diff --git a/segmentation/clitkExtractLymphStationsFilter.txx b/segmentation/clitkExtractLymphStationsFilter.txx index 70d540b..ea9e2b6 100644 --- a/segmentation/clitkExtractLymphStationsFilter.txx +++ b/segmentation/clitkExtractLymphStationsFilter.txx @@ -79,8 +79,6 @@ GenerateOutputInformation() { ExtractStation_8(); StopSubStep(); - DD(GetCurrentStepNumber()); - // Extract Station3P StartNewStep("Station 3P"); StartSubStep(); @@ -123,7 +121,7 @@ template void clitk::ExtractLymphStationsFilter:: GenerateInputRequestedRegion() { - DD("GenerateInputRequestedRegion (nothing?)"); + //DD("GenerateInputRequestedRegion (nothing?)"); } //-------------------------------------------------------------------- @@ -159,10 +157,10 @@ CheckForStation(std::string station) } // Define the starting support - if (found && GetComputeStation("8")) { + if (found && GetComputeStation(station)) { std::cout << "Station " << station << " already exists, but re-computation forced." << std::endl; } - if (!found || GetComputeStation("8")) { + if (!found || GetComputeStation(station)) { m_Working_Support = m_Mediastinum = GetAFDB()->template GetImage("Mediastinum", true); return true; } @@ -205,8 +203,14 @@ ExtractStation_8() if (CheckForStation("8")) { ExtractStation_8_SI_Limits(); ExtractStation_8_Post_Limits(); - ExtractStation_8_Ant_Limits(); + ExtractStation_8_Ant_Sup_Limits(); + ExtractStation_8_Ant_Inf_Limits(); + ExtractStation_8_LR_1_Limits(); + ExtractStation_8_LR_2_Limits(); ExtractStation_8_LR_Limits(); + ExtractStation_8_Ant_Injected_Limits(); + ExtractStation_8_Single_CCL_Limits(); + ExtractStation_8_Remove_Structures(); // Store image filenames into AFDB writeImage(m_ListOfStations["8"], "seg/Station8.mhd"); GetAFDB()->SetImageFilename("Station8", "seg/Station8.mhd"); @@ -224,9 +228,15 @@ 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_inf_Limits(); // Store image filenames into AFDB writeImage(m_ListOfStations["3P"], "seg/Station3P.mhd"); - GetAFDB()->SetImageFilename("Station3P", "seg/Station3P.mhd"); + GetAFDB()->SetImageFilename("Station3P", "seg/Station3P.mhd"); + WriteAFDB(); } } //-------------------------------------------------------------------- @@ -267,73 +277,5 @@ ExtractStation_4RL() { //-------------------------------------------------------------------- -//-------------------------------------------------------------------- -template -void -clitk::ExtractLymphStationsFilter:: -FindExtremaPointsInBronchus(MaskImagePointer input, - int direction, - double distance_max_from_center_point, - ListOfPointsType & LR, - ListOfPointsType & Ant, - ListOfPointsType & Post) -{ - - // Other solution ==> with auto bounding box ! (but pb to prevent to - // be too distant from the center point - - // Extract slices - std::vector slices; - clitk::ExtractSlices(input, 2, slices); - - // Loop on slices - bool found; - for(uint i=0; i(slices[i], 0, true, 10); - slices[i] = KeepLabels(slices[i], - GetBackgroundValue(), - GetForegroundValue(), 1, 1, true); - */ - - // ------- Find rightmost or leftmost point ------- - MaskSliceType::PointType LRMost; - found = - clitk::FindExtremaPointInAGivenDirection(slices[i], - GetBackgroundValue(), - 0, // axis XY - (direction==0?false:true), // right or left according to direction - LRMost); - // ------- Find postmost point ------- - MaskSliceType::PointType postMost; - found = - clitk::FindExtremaPointInAGivenDirection(slices[i], - GetBackgroundValue(), - 1, false, LRMost, - distance_max_from_center_point, - postMost); - // ------- Find antmost point ------- - MaskSliceType::PointType antMost; - found = - clitk::FindExtremaPointInAGivenDirection(slices[i], - GetBackgroundValue(), - 1, true, LRMost, - distance_max_from_center_point, - antMost); - // Only add point if found - if (found) { - // ------- Convert 2D to 3D points -------- - MaskImageType::PointType p; - clitk::PointsUtils::Convert2DTo3D(LRMost, input, i, p); - LR.push_back(p); - clitk::PointsUtils::Convert2DTo3D(antMost, input, i, p); - Ant.push_back(p); - clitk::PointsUtils::Convert2DTo3D(postMost, input, i, p); - Post.push_back(p); - } - } -} -//-------------------------------------------------------------------- #endif //#define CLITKBOOLEANOPERATORLABELIMAGEFILTER_TXX diff --git a/segmentation/clitkExtractLymphStationsGenericFilter.txx b/segmentation/clitkExtractLymphStationsGenericFilter.txx index ae1106a..c5c6378 100644 --- a/segmentation/clitkExtractLymphStationsGenericFilter.txx +++ b/segmentation/clitkExtractLymphStationsGenericFilter.txx @@ -70,6 +70,7 @@ SetOptionsFromArgsInfoToFilter(FilterType * f) f->SetAFDBFilename(mArgsInfo.afdb_arg); f->SetDistanceMaxToAnteriorPartOfTheSpine(mArgsInfo.maxAntSpine_arg); f->SetFuzzyThresholdForS8(mArgsInfo.fuzzyThresholdForS8_arg); + f->SetInjectedThresholdForS8(mArgsInfo.injectedThresholdForS8_arg); // Check multiple options for radius dilatation /* diff --git a/segmentation/clitkExtractMediastinumFilter.txx b/segmentation/clitkExtractMediastinumFilter.txx index 276c388..f9cd006 100644 --- a/segmentation/clitkExtractMediastinumFilter.txx +++ b/segmentation/clitkExtractMediastinumFilter.txx @@ -246,7 +246,7 @@ GenerateOutputInformation() { relPosFilter->SetInput(output); relPosFilter->SetInputObject(left_lung); // relPosFilter->SetInputObject(lung); - relPosFilter->AddOrientationType(RelPosFilterType::LeftTo); // warning left lung is at right ;) + relPosFilter->AddOrientationType(RelPosFilterType::AtRightTo); // warning left lung is at right ;) relPosFilter->SetIntermediateSpacing(GetIntermediateSpacing()); relPosFilter->SetFuzzyThreshold(GetFuzzyThreshold1()); relPosFilter->Update(); @@ -261,7 +261,7 @@ GenerateOutputInformation() { relPosFilter->SetInput(output); relPosFilter->SetInputObject(right_lung); //relPosFilter->SetInputObject(lung); - relPosFilter->AddOrientationType(RelPosFilterType::RightTo); + relPosFilter->AddOrientationType(RelPosFilterType::AtLeftTo); relPosFilter->SetIntermediateSpacing(GetIntermediateSpacing()); relPosFilter->SetFuzzyThreshold(GetFuzzyThreshold1()); relPosFilter->Update(); @@ -414,7 +414,7 @@ GenerateOutputInformation() { //-------------------------------------------------------------------- // Step 8: Lower limits from lung (need separate lung ?) - if (1) { + if (0) { // StartNewStep("Trial : minus segmented struct"); // MaskImagePointer heart = GetAFDB()->template GetImage ("heart"); // boolFilter = BoolFilterType::New();