X-Git-Url: https://git.creatis.insa-lyon.fr/pubgit/?a=blobdiff_plain;f=segmentation%2FclitkExtractLymphStation_7.txx;h=d6df42895faa021bd9298d4dc95cc11bdbc47402;hb=5a7da4aedae5c204bc55c187717193e5950f9a44;hp=e1cbe1f7750f47ae1319c5f8d8ef5c28c1580386;hpb=5ff46b355b778ae675ca60085430c062ea6a3945;p=clitk.git diff --git a/segmentation/clitkExtractLymphStation_7.txx b/segmentation/clitkExtractLymphStation_7.txx index e1cbe1f..d6df428 100644 --- a/segmentation/clitkExtractLymphStation_7.txx +++ b/segmentation/clitkExtractLymphStation_7.txx @@ -1,150 +1,254 @@ + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +ExtractStation_7_SetDefaultValues() +{ + SetFuzzyThreshold("7", "Bronchi", 0.1); + SetFuzzyThreshold("7", "LeftSuperiorPulmonaryVein", 0.3); + SetFuzzyThreshold("7", "RightSuperiorPulmonaryVein", 0.2); + SetFuzzyThreshold("7", "RightPulmonaryArtery", 0.3); + SetFuzzyThreshold("7", "LeftPulmonaryArtery", 0.5); + SetFuzzyThreshold("7", "SVC", 0.2); + SetS7_UseMostInferiorPartOnlyFlag(false); +} +//-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +ExtractStation_7() { + if (CheckForStation("7")) { + ExtractStation_7_SI_Limits(); + ExtractStation_7_RL_Interior_Limits(); + + // ExtractStation_7_Posterior_Limits(); + ExtractStation_8_Single_CCL_Limits(); + 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_7_SI_Limits() { - // Local variables - double m_CarinaZPositionInMM; - double m_MiddleLobeBronchusZPositionInMM; - + StartNewStep("[Station7] Inf/Sup mediastinum limits with carina/LLL-RMLBronchus"); // Get Inputs - m_Trachea = GetAFDB()->template GetImage ("trachea"); - m_CarinaZPositionInMM = GetAFDB()->GetPoint3D("carina", 2); - DD(m_CarinaZPositionInMM); - m_MiddleLobeBronchusZPositionInMM = GetAFDB()->GetPoint3D("rightMiddleLobeBronchus", 2); - DD(m_MiddleLobeBronchusZPositionInMM); - - /* Crop support : - Superior limit = carina - Inferior limit = origin right middle lobe bronchus */ - StartNewStep("Inf/Sup mediastinum limits with carina/bronchus"); + MaskImagePointer Trachea = GetAFDB()->template GetImage ("Trachea"); + + // Create line from A to B with + // A = upper border of LLL at left + // B = lower border of bronchus intermedius (BI) or RightMiddleLobeBronchus + ImagePointType A; + ImagePointType B; + FindLineForS7S8Separation(A, B); + + // // if option -> replace A[0] with B[0] + // if (GetS7_UseMostInferiorPartOnlyFlag()) { + // A[0] = B[0]; + // } + + // Use line to remove the inferior part + m_Working_Support = + SliceBySliceSetBackgroundFromSingleLine(m_Working_Support, GetBackgroundValue(), + A, B, 2, 0, false); + + // Get the CarinaZ position + double m_CarinaZ = FindCarina(); + + // Crop support m_Working_Support = - clitk::CropImageAlongOneAxis(m_Support, 2, - m_MiddleLobeBronchusZPositionInMM, - m_CarinaZPositionInMM, true, + clitk::CropImageAlongOneAxis(m_Working_Support, 2, + A[2], m_CarinaZ, true, GetBackgroundValue()); + // Crop trachea + m_Working_Trachea = + clitk::CropImageAlongOneAxis(Trachea, 2, + A[2], m_CarinaZ, true, + GetBackgroundValue()); + StopCurrentStep(m_Working_Support); - - /* Crop trachea - Superior limit = carina - Inferior limit = origin right middle lobe bronchus*/ - StartNewStep("Inf/Sup trachea limits with carina/bronchus"); - m_working_trachea = - clitk::CropImageAlongOneAxis(m_Trachea, 2, - m_MiddleLobeBronchusZPositionInMM, - m_CarinaZPositionInMM, true, - GetBackgroundValue()); - StopCurrentStep(m_working_trachea); - - m_Station7 = m_Working_Support; + m_ListOfStations["7"] = m_Working_Support; } //-------------------------------------------------------------------- + //-------------------------------------------------------------------- template void clitk::ExtractLymphStationsFilter:: -ExtractStation_7_RL_Limits() +ExtractStation_7_RL_Interior_Limits() { // ---------------------------------------------------------------- - // Separate trachea in two CCL - StartNewStep("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()) { - // DD(iter.GetIndex()); - if (iter.Get() > maxLabel) maxLabel = iter.Get(); - --iter; + StartNewStep("[Station7] RL limits with bronchi"); + + /* + Slice by Slice, consider most Left point of the Right + bronchus. Remove the Ant/Right corner + */ + + // First consider bronchi + FindLeftAndRightBronchi(); + m_RightBronchus = GetAFDB()->template GetImage ("RightBronchus"); + m_LeftBronchus = GetAFDB()->template GetImage ("LeftBronchus"); + + // Resize like m_Working_Support + m_LeftBronchus = + clitk::ResizeImageLike(m_LeftBronchus, m_Working_Support, GetBackgroundValue()); + m_RightBronchus = + clitk::ResizeImageLike(m_RightBronchus, m_Working_Support, GetBackgroundValue()); + + // Extract slices, Label, compute centroid, keep most central connected component + std::vector slices_leftbronchus; + std::vector slices_rightbronchus; + std::vector slices_support; + clitk::ExtractSlices(m_LeftBronchus, 2, slices_leftbronchus); + clitk::ExtractSlices(m_RightBronchus, 2, slices_rightbronchus); + clitk::ExtractSlices(m_Working_Support, 2, slices_support); + + // Keep only the CCL of the bronchus with the closest to the center + // Loop on slices for left bronchus + 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 + } + + // Loop on slices for right bronchus + 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); + + // For Right bronchus, Find most Left point. Remove corner Ant/Right corner + for(uint i=0; i(slices_rightbronchus[i], GetBackgroundValue(), 0, false, p_left); + if (b) { + b = clitk::FindExtremaPointInAGivenDirection(slices_rightbronchus[i], GetBackgroundValue(), 1, false, p_post); + } + if (b) { + MaskSlicePointType p = p_left; + p[1] = p_post[1]; + MaskSliceIndexType pi; + slices_rightbronchus[i]->TransformPhysicalPointToIndex(p, pi); + + // Build region to remove + MaskSliceRegionType region = slices_rightbronchus[i]->GetLargestPossibleRegion(); + MaskSliceIndexType index = region.GetIndex(); + MaskSliceSizeType size = region.GetSize(); + size[0] = pi[0] - index[0]; + size[1] = pi[1] - index[1]; + region.SetSize(size); + + // Fill region with Background value + clitk::FillRegionWithValue(slices_support[i], GetBackgroundValue(), region); } - iter.PreviousLine(); } - if (maxLabel < 2) { - clitkExceptionMacro("First slice form Carina does not seems to seperate the two main bronchus. Abort"); + + // For Left bronchus, Find most Right point. Remove corner Ant/Left corner + for(uint i=0; i(slices_leftbronchus[i], GetBackgroundValue(), 0, true, p_right); + if (b) { + b = clitk::FindExtremaPointInAGivenDirection(slices_rightbronchus[i], GetBackgroundValue(), 1, false, p_post); + } + if (b) { + MaskSlicePointType p = p_right; + p[1] = p_post[1]; + MaskSliceIndexType pi; + slices_leftbronchus[i]->TransformPhysicalPointToIndex(p, pi); + + /* typedef itk::ImageRegionIterator IteratorType; + IteratorType iter(input, region); + iter.GoToBegin(); + while (!iter.IsAtEnd()) { + MaskSliceIndexType index = iter.GetIndex(); + if (index[0] > pi[0]) && (index[1] > pi[1]) iter.Set(GetBackgroundValue()); + + ++iter; + } + */ + + + // Build region to remove + MaskSliceRegionType region = slices_leftbronchus[i]->GetLargestPossibleRegion(); + MaskSliceIndexType index = region.GetIndex(); + MaskSliceSizeType size = region.GetSize(); + index[0] = pi[0]; + size[0] = slices_leftbronchus[i]->GetLargestPossibleRegion().GetSize()[0] - pi[0]; + size[1] = pi[1] - index[1]; + region.SetSize(size); + region.SetIndex(index); + + // Fill region with Background value + clitk::FillRegionWithValue(slices_support[i], GetBackgroundValue(), region); + } } - // Compute centroid of both parts to identify the left from the right bronchus - typedef long LabelType; - static const unsigned int Dim = ImageType::ImageDimension; - typedef itk::ShapeLabelObject< LabelType, Dim > LabelObjectType; - typedef itk::LabelMap< LabelObjectType > LabelMapType; - typedef itk::LabelImageToLabelMapFilter ImageToMapFilterType; - typename ImageToMapFilterType::Pointer imageToLabelFilter = ImageToMapFilterType::New(); - typedef itk::ShapeLabelMapFilter ShapeFilterType; - typename ShapeFilterType::Pointer statFilter = ShapeFilterType::New(); - imageToLabelFilter->SetBackgroundValue(GetBackgroundValue()); - imageToLabelFilter->SetInput(m_working_trachea); - statFilter->SetInput(imageToLabelFilter->GetOutput()); - statFilter->Update(); - typename LabelMapType::Pointer labelMap = statFilter->GetOutput(); - - ImagePointType C1 = labelMap->GetLabelObject(1)->GetCentroid(); - ImagePointType C2 = labelMap->GetLabelObject(2)->GetCentroid(); - - ImagePixelType leftLabel; - ImagePixelType rightLabel; - if (C1[0] < C2[0]) { leftLabel = 1; rightLabel = 2; } - else { leftLabel = 2; rightLabel = 1; } - DD(leftLabel); - DD(rightLabel); - - StopCurrentStep(m_working_trachea); - - //----------------------------------------------------- - StartNewStep("Left limits with bronchus (slice by slice)"); - // Select LeftLabel (set one label to Backgroundvalue) - m_LeftBronchus = - SetBackground(m_working_trachea, m_working_trachea, - rightLabel, GetBackgroundValue()); - m_RightBronchus = - SetBackground(m_working_trachea, m_working_trachea, - leftLabel, GetBackgroundValue()); - writeImage(m_LeftBronchus, "left.mhd"); - writeImage(m_RightBronchus, "right.mhd"); + m_Working_Support = clitk::JoinSlices(slices_support, m_Working_Support, 2); + // Also remove what is at right of the Right bronchus (left respectively) m_Working_Support = - clitk::SliceBySliceRelativePosition(m_Working_Support, - m_LeftBronchus, 2, - GetFuzzyThreshold(), "RightTo", - true, 4); - m_Working_Support = - clitk::SliceBySliceRelativePosition(m_Working_Support, - m_RightBronchus, - 2, GetFuzzyThreshold(), "LeftTo", - true, 4); - m_Working_Support = - clitk::SliceBySliceRelativePosition(m_Working_Support, - m_LeftBronchus, - 2, GetFuzzyThreshold(), "AntTo", - true, 4, true); // NOT - m_Working_Support = - clitk::SliceBySliceRelativePosition(m_Working_Support, - m_RightBronchus, - 2, GetFuzzyThreshold(), "AntTo", - true, 4, true); - m_Working_Support = - clitk::SliceBySliceRelativePosition(m_Working_Support, - m_LeftBronchus, - 2, GetFuzzyThreshold(), "PostTo", - true, 4, true); + clitk::SliceBySliceRelativePosition(m_Working_Support, m_LeftBronchus, 2, + GetFuzzyThreshold("7", "Bronchi"), "NotLeftTo", + false, 3, false); m_Working_Support = - clitk::SliceBySliceRelativePosition(m_Working_Support, - m_RightBronchus, - 2, GetFuzzyThreshold(), "PostTo", - true, 4, true); - m_Station7 = m_Working_Support; - StopCurrentStep(m_Station7); + clitk::SliceBySliceRelativePosition(m_Working_Support, m_RightBronchus, 2, + GetFuzzyThreshold("7", "Bronchi"), "NotRightTo", + false, 3, false); + + // SECOND PART + + StopCurrentStep(m_Working_Support); } //-------------------------------------------------------------------- @@ -153,161 +257,310 @@ ExtractStation_7_RL_Limits() template void clitk::ExtractLymphStationsFilter:: -ExtractStation_7_Posterior_Limits() +ExtractStation_7_RL_Limits_OLD() { - StartNewStep("Posterior limits"); - - // Left Bronchus slices - typedef clitk::ExtractSliceFilter ExtractSliceFilterType; - typedef typename ExtractSliceFilterType::SliceType SliceType; - typename ExtractSliceFilterType::Pointer - extractSliceFilter = ExtractSliceFilterType::New(); - extractSliceFilter->SetInput(m_LeftBronchus); - extractSliceFilter->SetDirection(2); - extractSliceFilter->Update(); - std::vector leftBronchusSlices; - extractSliceFilter->GetOutputSlices(leftBronchusSlices); + // ---------------------------------------------------------------- + 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); - // Right Bronchus slices - extractSliceFilter = ExtractSliceFilterType::New(); - extractSliceFilter->SetInput(m_RightBronchus); - extractSliceFilter->SetDirection(2); - extractSliceFilter->Update(); - std::vector rightBronchusSlices ; - extractSliceFilter->GetOutputSlices(rightBronchusSlices); + // Version #1 with limit = centroid of the bronchus (OUTDATED) + // Step1 = keep only the CCL of the bronchus with the closest to the center + // Step2 = SliceBySlice Rel pos to both bronchi + + // 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 + } - assert(leftBronchusSlices.size() == rightBronchusSlices.size()); + 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 + } - std::vector leftPoints; - std::vector rightPoints; - for(uint i=0; i(leftBronchusSlices[i], 0, true, 10); - leftBronchusSlices[i] = KeepLabels(leftBronchusSlices[i], - GetBackgroundValue(), - GetForegroundValue(), 1, 1, true); - rightBronchusSlices[i] = Labelize(rightBronchusSlices[i], 0, true, 10); - rightBronchusSlices[i] = KeepLabels(rightBronchusSlices[i], - GetBackgroundValue(), - GetForegroundValue(), 1, 1, true); - double distance_max_from_center_point = 15; - - // ------- Find point in left Bronchus ------- - // find right most point in left = rightMost - SliceType::PointType a; - SliceType::PointType rightMost = - clitk::FindExtremaPointInAGivenDirection(leftBronchusSlices[i], - GetBackgroundValue(), - 0, false, a, 0); - // find post most point in left, not far away from rightMost - SliceType::PointType p = - clitk::FindExtremaPointInAGivenDirection(leftBronchusSlices[i], - GetBackgroundValue(), - 1, false, rightMost, - distance_max_from_center_point); - MaskImageType::PointType pp; - pp[0] = p[0]; pp[1] = p[1]; - pp[2] = i*m_LeftBronchus->GetSpacing()[2] + m_LeftBronchus->GetOrigin()[2]; - leftPoints.push_back(pp); - - // ------- Find point in right Bronchus ------- - // find left most point in right = leftMost - SliceType::PointType leftMost = - clitk::FindExtremaPointInAGivenDirection(rightBronchusSlices[i], - GetBackgroundValue(), - 0, true, a, 0); - // find post most point in left, not far away from leftMost - p = clitk::FindExtremaPointInAGivenDirection(rightBronchusSlices[i], - GetBackgroundValue(), - 1, false, leftMost, - distance_max_from_center_point); - pp[0] = p[0]; pp[1] = p[1]; - pp[2] = i*m_RightBronchus->GetSpacing()[2] + m_RightBronchus->GetOrigin()[2]; - rightPoints.push_back(pp); + // Joint slices + m_LeftBronchus = clitk::JoinSlices(slices_leftbronchus, m_LeftBronchus, 2); + m_RightBronchus = clitk::JoinSlices(slices_rightbronchus, m_RightBronchus, 2); + + writeImage(m_LeftBronchus, "step-left.mhd"); + writeImage(m_RightBronchus, "step-right.mhd"); + + m_Working_Support = + clitk::SliceBySliceRelativePosition(m_Working_Support, m_LeftBronchus, 2, + GetFuzzyThreshold("7", "Bronchi"), "RightTo", + false, 3, false); + StopCurrentStep(m_Working_Support); + + + // ---------------------------------------------------------------- + StartNewStep("[Station7] Limits with bronchus : LeftTo the right bronchus"); + m_Working_Support = + clitk::SliceBySliceRelativePosition(m_Working_Support, m_RightBronchus, 2, + GetFuzzyThreshold("7", "Bronchi"), "LeftTo", + false, 3, false); + StopCurrentStep(m_Working_Support); + + + // ---------------------------------------------------------------- + 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(GetFuzzyThreshold("7", "LeftSuperiorPulmonaryVein")); + sliceRelPosFilter->AddOrientationTypeString("NotLeftTo"); + sliceRelPosFilter->AddOrientationTypeString("NotAntTo"); + sliceRelPosFilter->SetIntermediateSpacingFlag(true); + sliceRelPosFilter->SetIntermediateSpacing(3); + sliceRelPosFilter->SetUniqueConnectedComponentBySliceFlag(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; } - // DEBUG - std::ofstream osl; - openFileForWriting(osl, "left.txt"); - osl << "LANDMARKS1" << std::endl; - std::ofstream osr; - openFileForWriting(osr, "right.txt"); - osr << "LANDMARKS1" << std::endl; - for(uint i=0; itemplate GetImage("RightSuperiorPulmonaryVein"); + typedef SliceBySliceRelativePositionFilter SliceRelPosFilterType; + typename SliceRelPosFilterType::Pointer sliceRelPosFilter = SliceRelPosFilterType::New(); + sliceRelPosFilter->SetInput(m_Working_Support); + sliceRelPosFilter->SetInputObject(RightSuperiorPulmonaryVein); + sliceRelPosFilter->SetDirection(2); + sliceRelPosFilter->SetFuzzyThreshold(GetFuzzyThreshold("7", "RightSuperiorPulmonaryVein")); + sliceRelPosFilter->AddOrientationTypeString("NotRightTo"); + sliceRelPosFilter->AddOrientationTypeString("NotAntTo"); + sliceRelPosFilter->AddOrientationTypeString("NotPostTo"); + sliceRelPosFilter->SetIntermediateSpacingFlag(true); + sliceRelPosFilter->SetIntermediateSpacing(3); + sliceRelPosFilter->SetUniqueConnectedComponentBySliceFlag(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; } - osl.close(); - osr.close(); - // Now uses these points to limit, slice by slice - // http://www.gamedev.net/community/forums/topic.asp?topic_id=542870 - /* - Assuming the points are (Ax,Ay) (Bx,By) and (Cx,Cy), you need to compute: - (Bx - Ax) * (Cy - Ay) - (By - Ay) * (Cx - Ax) - This will equal zero if the point C is on the line formed by - points A and B, and will have a different sign depending on the - side. Which side this is depends on the orientation of your (x,y) - coordinates, but you can plug test values for A,B and C into this - formula to determine whether negative values are to the left or to - the right. - => to accelerate, start with formula, when change sign -> stop and fill + // ---------------------------------------------------------------- + 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(GetFuzzyThreshold("7", "RightPulmonaryArtery")); + sliceRelPosFilter->AddOrientationTypeString("NotAntTo"); + sliceRelPosFilter->SetIntermediateSpacingFlag(true); + sliceRelPosFilter->SetIntermediateSpacing(3); + sliceRelPosFilter->SetUniqueConnectedComponentBySliceFlag(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(GetFuzzyThreshold("7", "LeftPulmonaryArtery")); + sliceRelPosFilter->AddOrientationTypeString("NotAntTo"); + sliceRelPosFilter->SetIntermediateSpacingFlag(true); + sliceRelPosFilter->SetIntermediateSpacing(3); + sliceRelPosFilter->SetUniqueConnectedComponentBySliceFlag(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(GetFuzzyThreshold("7", "SVC")); + sliceRelPosFilter->AddOrientationTypeString("NotRightTo"); + sliceRelPosFilter->AddOrientationTypeString("NotAntTo"); + sliceRelPosFilter->SetIntermediateSpacingFlag(true); + sliceRelPosFilter->SetIntermediateSpacing(3); + sliceRelPosFilter->SetUniqueConnectedComponentBySliceFlag(false); + sliceRelPosFilter->SetAutoCropFlag(true); + sliceRelPosFilter->IgnoreEmptySliceObjectFlagOn(); + sliceRelPosFilter->Update(); + m_Working_Support = sliceRelPosFilter->GetOutput(); + StopCurrentStep(m_Working_Support); + + // End + m_ListOfStations["7"] = m_Working_Support; +} +//-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +ExtractStation_7_Posterior_Limits() +{ + 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); */ - 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 = leftPoints[i]; - B = rightPoints[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; + + // 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; + std::ofstream ospl; openFileForWriting(ospl, "ospl.txt"); ospl << "LANDMARKS1" << std::endl; + std::ofstream osrr; openFileForWriting(osrr, "osrr.txt"); osrr << "LANDMARKS1" << std::endl; + std::ofstream osar; openFileForWriting(osar, "osar.txt"); osar << "LANDMARKS1" << std::endl; + std::ofstream ospr; openFileForWriting(ospr, "ospr.txt"); ospr << "LANDMARKS1" << std::endl; + + for(uint i=0; i(m_Working_Support, + m_PostMostInRightBronchus, + m_PostMostInLeftBronchus, + GetBackgroundValue(), 1, -10); + // If needed -> can do the same with AntMost. + // 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("7", "AzygousVein"); + Remove_Structures("7", "Aorta"); + Remove_Structures("7", "RightPulmonaryArtery"); + Remove_Structures("7", "LeftPulmonaryArtery"); + Remove_Structures("7", "LeftSuperiorPulmonaryVein"); + Remove_Structures("7", "PulmonaryTrunk"); + Remove_Structures("7", "VertebralBody"); - // Set output - m_Station7 = m_Working_Support; + // Keep only one CCL by slice (before removing Esophagus) + // DD("SliceBySliceKeepMainCCL"); + + // TODO -> replace by keep the one that contains point at the middle of the line between the bronchus + // -> new function "keep/select" the ccl that contains this point (2D) + + //m_Working_Support = clitk::SliceBySliceKeepMainCCL(m_Working_Support, GetBackgroundValue(), GetForegroundValue()); + + Remove_Structures("7", "Esophagus"); + + // END + StopCurrentStep(m_Working_Support); + m_ListOfStations["7"] = m_Working_Support; + return; } //-------------------------------------------------------------------- + + +