segmentation/clitkExtractLymphStation_3P.txx

   1
   2 #include <itkBinaryDilateImageFilter.h>
   3 #include <itkMirrorPadImageFilter.h>
   4
   5 //--------------------------------------------------------------------
   6 template <class ImageType>
   7 void
   8 clitk::ExtractLymphStationsFilter<ImageType>::
   9 ExtractStation_3P_SetDefaultValues()
  10 {
  11 }
  12 //--------------------------------------------------------------------
  13
  14
  15 //--------------------------------------------------------------------
  16 template <class ImageType>
  17 void
  18 clitk::ExtractLymphStationsFilter<ImageType>::
  19 ExtractStation_3P_SI_Limits()
  20 {
  21   /*
  22     Apex of the chest & Carina.
  23   */
  24   StartNewStep("[Station 3P] Inf/Sup limits with apex of the chest and carina");
  25
  26   // Get Carina position (has been determined in Station8)
  27   m_CarinaZ = GetAFDB()->GetDouble("CarinaZ");
  28
  29   // Get Apex of the Chest. The "lungs" structure is autocroped so we
  30   // consider the most superior point.
  31   MaskImagePointer Lungs = GetAFDB()->template GetImage<MaskImageType>("Lungs");
  32   MaskImageIndexType index = Lungs->GetLargestPossibleRegion().GetIndex();
  33   index += Lungs->GetLargestPossibleRegion().GetSize();
  34   MaskImagePointType p;
  35   Lungs->TransformIndexToPhysicalPoint(index, p);
  36   p[2] -= 5; // 5 mm below because the autocrop is slightly above the apex
  37   double m_ApexOfTheChest = p[2];
  38
  39   /* Crop support :
  40      Superior limit = carina
  41      Inferior limit = Apex of the chest */
  42   m_Working_Support =
  43     clitk::CropImageAlongOneAxis<MaskImageType>(m_Working_Support, 2,
  44                                                 m_CarinaZ,
  45                                                 m_ApexOfTheChest, true,
  46                                                 GetBackgroundValue());
  47
  48   StopCurrentStep<MaskImageType>(m_Working_Support);
  49   m_ListOfStations["3P"] = m_Working_Support;
  50 }
  51 //--------------------------------------------------------------------
  52
  53
  54 //--------------------------------------------------------------------
  55 template <class ImageType>
  56 void
  57 clitk::ExtractLymphStationsFilter<ImageType>::
  58 ExtractStation_3P_Remove_Structures()
  59 {
  60   /*
  61     3 - (sup) remove Aorta, VB, subcl
  62     not LR subcl ! -> a séparer LR ?
  63     (inf) remove Eso, Aorta, Azygousvein
  64   */
  65
  66   StartNewStep("[Station 3P] Remove some structures.");
  67
  68   MaskImagePointer Aorta = GetAFDB()->template GetImage<MaskImageType>("Aorta");
  69   clitk::AndNot<MaskImageType>(m_Working_Support, Aorta, GetBackgroundValue());
  70
  71   MaskImagePointer VertebralBody = GetAFDB()->template GetImage<MaskImageType>("VertebralBody");
  72   clitk::AndNot<MaskImageType>(m_Working_Support, VertebralBody);
  73
  74   MaskImagePointer SubclavianArtery = GetAFDB()->template GetImage<MaskImageType>("SubclavianArtery");
  75   clitk::AndNot<MaskImageType>(m_Working_Support, SubclavianArtery);
  76
  77   MaskImagePointer Esophagus = GetAFDB()->template GetImage<MaskImageType>("Esophagus");
  78   clitk::AndNot<MaskImageType>(m_Working_Support, Esophagus);
  79
  80   MaskImagePointer AzygousVein = GetAFDB()->template GetImage<MaskImageType>("AzygousVein");
  81   clitk::AndNot<MaskImageType>(m_Working_Support, AzygousVein);
  82
  83   MaskImagePointer Thyroid = GetAFDB()->template GetImage<MaskImageType>("Thyroid");
  84   clitk::AndNot<MaskImageType>(m_Working_Support, Thyroid);
  85
  86   StopCurrentStep<MaskImageType>(m_Working_Support);
  87   m_ListOfStations["3P"] = m_Working_Support;
  88 }
  89 //--------------------------------------------------------------------
  90
  91
  92 //--------------------------------------------------------------------
  93 template <class ImageType>
  94 void
  95 clitk::ExtractLymphStationsFilter<ImageType>::
  96 ExtractStation_3P_Ant_Limits()
  97 {
  98   /*
  99     Ant Post limit :
 100
 101     Anteriorly, it is in contact with the posterior aspect of Stations
 102     1–2 superiorly (Fig. 2A–C) and with Stations 4R and 4L inferiorly
 103     (Fig. 2D–I and 3A–C). The anterior limit of Station 3P is kept
 104     posterior to the trachea, which is defined by an imaginary
 105     horizontal line running along the posterior wall of the trachea
 106     (Fig. 2B,E, red line). Posteriorly, it is delineated along the
 107     anterior and lateral borders of the vertebral body until an
 108     imaginary horizontal line running 1 cm posteriorly to the
 109     anterior border of the vertebral body (Fig. 2D).
 110
 111     1 - post to the trachea : define an imaginary line just post the
 112     Trachea and remove what is anterior.
 113
 114   */
 115   StartNewStep("[Station 3P] Ant limits with Trachea ");
 116
 117   // Load Trachea
 118   MaskImagePointer Trachea = GetAFDB()->template GetImage<MaskImageType>("Trachea");
 119
 120   // Crop like current support (need by SliceBySliceSetBackgroundFromLineSeparation after)
 121   Trachea =
 122     clitk::ResizeImageLike<MaskImageType>(Trachea, m_Working_Support, GetBackgroundValue());
 123
 124   // Slice by slice, determine the most post point of the trachea (A)
 125   // and consider a second point on the same horizontal line (B)
 126   std::vector<MaskImagePointType> p_A;
 127   std::vector<MaskImagePointType> p_B;
 128   std::vector<typename MaskSliceType::Pointer> slices;
 129   clitk::ExtractSlices<MaskImageType>(Trachea, 2, slices);
 130   MaskImagePointType p;
 131   typename MaskSliceType::PointType sp;
 132   for(uint i=0; i<slices.size(); i++) {
 133     // First only consider the main CCL (trachea, not bronchus)
 134     slices[i] = Labelize<MaskSliceType>(slices[i], 0, true, 100);
 135     slices[i] = KeepLabels<MaskSliceType>(slices[i], GetBackgroundValue(),
 136                                           GetForegroundValue(), 1, 1, true);
 137     // Find most posterior point
 138     clitk::FindExtremaPointInAGivenDirection<MaskSliceType>(slices[i], GetBackgroundValue(),
 139                                                             1, false, sp);
 140     clitk::PointsUtils<MaskImageType>::Convert2DTo3D(sp, Trachea, i, p);
 141     p_A.push_back(p);
 142     p[0] -= 20;
 143     p_B.push_back(p);
 144   }
 145   clitk::WriteListOfLandmarks<MaskImageType>(p_A, "trachea-post.txt");
 146
 147   // Remove Ant part above those lines
 148   clitk::SliceBySliceSetBackgroundFromLineSeparation<MaskImageType>(m_Working_Support,
 149                                                                     p_A, p_B,
 150                                                                     GetBackgroundValue(),
 151                                                                     1, 10);
 152   // End, release memory
 153   GetAFDB()->template ReleaseImage<MaskImageType>("Trachea");
 154   StopCurrentStep<MaskImageType>(m_Working_Support);
 155   m_ListOfStations["3P"] = m_Working_Support;
 156 }
 157 //--------------------------------------------------------------------
 158
 159
 160 //--------------------------------------------------------------------
 161 template <class ImageType>
 162 void
 163 clitk::ExtractLymphStationsFilter<ImageType>::
 164 ExtractStation_3P_Post_Limits()
 165 {
 166   /*
 167     Ant Post limit :
 168
 169     Anteriorly, it is in contact with the posterior aspect of Stations
 170     1–2 superiorly (Fig. 2A–C) and with Stations 4R and 4L inferiorly
 171     (Fig. 2D–I and 3A–C). The anterior limit of Station 3P is kept
 172     posterior to the trachea, which is defined by an imaginary
 173     horizontal line running along the posterior wall of the trachea
 174     (Fig. 2B,E, red line). Posteriorly, it is delineated along the
 175     anterior and lateral borders of the vertebral body until an
 176     imaginary horizontal line running 1 cm posteriorly to the
 177     anterior border of the vertebral body (Fig. 2D).
 178
 179     2 - post to the trachea : define an imaginary line just post the
 180     Trachea and remove what is anterior.
 181
 182   */
 183   StartNewStep("[Station 3P] Post limits with VertebralBody ");
 184
 185   // Load VertebralBody
 186   MaskImagePointer VertebralBody = GetAFDB()->template GetImage<MaskImageType>("VertebralBody");
 187
 188   // Crop like current support (need by SliceBySliceSetBackgroundFromLineSeparation after)
 189   VertebralBody = clitk::ResizeImageLike<MaskImageType>(VertebralBody, m_Working_Support, GetBackgroundValue());
 190
 191   writeImage<MaskImageType>(VertebralBody, "vb.mhd");
 192
 193   // Compute VertebralBody most Ant position (again because slices
 194   // changes). Slice by slice, determine the most post point of the
 195   // trachea (A) and consider a second point on the same horizontal
 196   // line (B)
 197   std::vector<MaskImagePointType> p_A;
 198   std::vector<MaskImagePointType> p_B;
 199   std::vector<typename MaskSliceType::Pointer> slices;
 200   clitk::ExtractSlices<MaskImageType>(VertebralBody, 2, slices);
 201   MaskImagePointType p;
 202   typename MaskSliceType::PointType sp;
 203   for(uint i=0; i<slices.size(); i++) {
 204     // Find most anterior point
 205     bool found = clitk::FindExtremaPointInAGivenDirection<MaskSliceType>(slices[i], GetBackgroundValue(),
 206                                                                          1, true, sp);
 207
 208     // If the VertebralBody stop superiorly before the end of
 209     // m_Working_Support, we consider the same previous point.
 210     if (!found) {
 211       p = p_A.back();
 212       p[2] += VertebralBody->GetSpacing()[2];
 213       p_A.push_back(p);
 214       p = p_B.back();
 215       p[2] += VertebralBody->GetSpacing()[2];
 216       p_B.push_back(p);
 217     }
 218     else {
 219       clitk::PointsUtils<MaskImageType>::Convert2DTo3D(sp, VertebralBody, i, p);
 220       p[1] += 10; // Consider 10 mm more post
 221       p_A.push_back(p);
 222       p[0] -= 20;
 223       p_B.push_back(p);
 224     }
 225   }
 226   clitk::WriteListOfLandmarks<MaskImageType>(p_A, "vb-ant.txt");
 227
 228
 229   // Remove Ant part above those lines
 230   clitk::SliceBySliceSetBackgroundFromLineSeparation<MaskImageType>(m_Working_Support,
 231                                                                     p_A, p_B,
 232                                                                     GetBackgroundValue(),
 233                                                                     1, -10);
 234   writeImage<MaskImageType>(m_Working_Support, "a.mhd");
 235
 236   StopCurrentStep<MaskImageType>(m_Working_Support);
 237   m_ListOfStations["3P"] = m_Working_Support;
 238 }
 239 //--------------------------------------------------------------------
 240
 241
 242 //--------------------------------------------------------------------
 243 template <class ImageType>
 244 void
 245 clitk::ExtractLymphStationsFilter<ImageType>::
 246 ExtractStation_3P_LR_sup_Limits()
 247 {
 248   /*
 249     "On the right side, the limit is deﬁned by the air–soft-tissue
 250     interface. On the left side, it is deﬁned by the air–tissue
 251     interface superiorly (Fig. 2A–C) and the aorta inferiorly
 252     (Figs. 2D–I and 3A–C)."
 253
 254     sup :
 255     Resizelike support : Trachea, SubclavianArtery
 256     Trachea, slice by slice, get centroid trachea
 257     SubclavianArtery, slice by slice, CCL
 258     prendre la 1ère à L et R, not at Left
 259
 260   */
 261   StartNewStep("[Station 3P] Left/Right limits (superior part) ");
 262
 263   // Load structures
 264   MaskImagePointer Trachea = GetAFDB()->template GetImage<MaskImageType>("Trachea");
 265   MaskImagePointer SubclavianArtery = GetAFDB()->template GetImage<MaskImageType>("SubclavianArtery");
 266
 267   // Crop like current support
 268   Trachea = clitk::ResizeImageLike<MaskImageType>(Trachea, m_Working_Support, GetBackgroundValue());
 269   SubclavianArtery = clitk::ResizeImageLike<MaskImageType>(SubclavianArtery, m_Working_Support, GetBackgroundValue());
 270
 271   writeImage<MaskImageType>(Trachea, "tr.mhd");
 272   writeImage<MaskImageType>(SubclavianArtery, "sca.mhd");
 273
 274   // Get list of slices
 275   std::vector<MaskSlicePointer> slices_support;
 276   std::vector<MaskSlicePointer> slices_trachea;
 277   std::vector<MaskSlicePointer> slices_subclavianartery;
 278   clitk::ExtractSlices<MaskImageType>(m_Working_Support, 2, slices_support);
 279   clitk::ExtractSlices<MaskImageType>(Trachea, 2, slices_trachea);
 280   clitk::ExtractSlices<MaskImageType>(SubclavianArtery, 2, slices_subclavianartery);
 281
 282   // Loop on slices
 283   std::vector<MaskImagePointType> points;
 284   MaskImagePointType p;
 285   for(uint i=0; i<slices_support.size(); i++) {
 286     // Get Trachea centroid
 287     std::vector<typename MaskSliceType::PointType> centroids;
 288     typename MaskSliceType::PointType c;
 289     ComputeCentroids<MaskSliceType>(slices_trachea[i], GetBackgroundValue(), centroids);
 290     c = centroids[1];
 291
 292     // [debug] Store point
 293     clitk::PointsUtils<MaskImageType>::Convert2DTo3D(centroids[1], Trachea, i, p);
 294     points.push_back(p);
 295
 296     // Get Right and Left CCL in SubclavianArtery
 297     slices_subclavianartery[i] = Labelize<MaskSliceType>(slices_subclavianartery[i], 0, true, 10);
 298     ComputeCentroids<MaskSliceType>(slices_subclavianartery[i], GetBackgroundValue(), centroids);
 299
 300     if (centroids.size() > 1) {
 301       // Determine the one at Right/Left -> first after Trachea
 302       // centroid
 303       typename MaskSliceType::PointType right;
 304       typename MaskSliceType::PointType left;
 305       int label_right=-1;
 306       int label_left=-1;
 307       right[0] = c[0]-100;
 308       left[0] = c[0]+100;
 309       for(uint j=1; j<centroids.size(); j++) {
 310         if (centroids[j][0] < c[0]) { // At Right of Trachea centroid
 311           if (centroids[j][0] >= right[0]) {
 312             right = centroids[j];
 313             label_right = j;
 314           }
 315         }
 316         if (centroids[j][0] > c[0]) { // At Left of Trachea centroid
 317           if (centroids[j][0] <= left[0]) {
 318             left = centroids[j];
 319             label_left = j;
 320           }
 321         }
 322       }
 323
 324       if (label_right != -1) {
 325
 326         // Debug points
 327         clitk::PointsUtils<MaskImageType>::Convert2DTo3D(centroids[label_right], SubclavianArtery, i, p);
 328         points.push_back(p);
 329
 330         // Set Background and ForegroundValue according to label_right
 331         MaskSlicePointer object =
 332           clitk::Binarize<MaskSliceType>(slices_subclavianartery[i], label_right, label_right,
 333                                          GetBackgroundValue(), GetForegroundValue());
 334
 335         // Relative Position : not at Right
 336         typedef clitk::AddRelativePositionConstraintToLabelImageFilter<MaskSliceType> RelPosFilterType;
 337         typename RelPosFilterType::Pointer relPosFilter = RelPosFilterType::New();
 338         relPosFilter->VerboseStepFlagOff();
 339         relPosFilter->WriteStepFlagOff();
 340         relPosFilter->SetBackgroundValue(GetBackgroundValue());
 341         relPosFilter->SetInput(slices_support[i]);
 342         relPosFilter->SetInputObject(object);
 343         relPosFilter->AddOrientationTypeString("R");
 344         relPosFilter->SetInverseOrientationFlag(true);
 345         //      relPosFilter->SetIntermediateSpacing(3);
 346         relPosFilter->SetIntermediateSpacingFlag(false);
 347         relPosFilter->SetFuzzyThreshold(0.7);
 348         relPosFilter->AutoCropFlagOff(); // important ! because we join the slices after this loop
 349         relPosFilter->Update();
 350         slices_support[i] = relPosFilter->GetOutput();
 351
 352         // Relative Position : not Anterior
 353         relPosFilter = RelPosFilterType::New();
 354         relPosFilter->VerboseStepFlagOff();
 355         relPosFilter->WriteStepFlagOff();
 356         relPosFilter->SetBackgroundValue(GetBackgroundValue());
 357         relPosFilter->SetInput(slices_support[i]);
 358         relPosFilter->SetInputObject(object);
 359         relPosFilter->AddOrientationTypeString("A");
 360         relPosFilter->SetInverseOrientationFlag(true);
 361         //      relPosFilter->SetIntermediateSpacing(3);
 362         relPosFilter->SetIntermediateSpacingFlag(false);
 363         relPosFilter->SetFuzzyThreshold(0.7);
 364         relPosFilter->AutoCropFlagOff(); // important ! because we join the slices after this loop
 365         relPosFilter->Update();
 366         slices_support[i] = relPosFilter->GetOutput();
 367
 368       } // End RelativePosition for Right
 369
 370
 371       if (label_left != -1) {
 372
 373         // Debug points
 374         clitk::PointsUtils<MaskImageType>::Convert2DTo3D(centroids[label_left], SubclavianArtery, i, p);
 375         points.push_back(p);
 376
 377         // Set Background and ForegroundValue according to label_right
 378         MaskSlicePointer object =
 379           clitk::Binarize<MaskSliceType>(slices_subclavianartery[i], label_left, label_left,
 380                                          GetBackgroundValue(), GetForegroundValue());
 381
 382         // Relative Position : not at Right
 383         typedef clitk::AddRelativePositionConstraintToLabelImageFilter<MaskSliceType> RelPosFilterType;
 384         typename RelPosFilterType::Pointer relPosFilter = RelPosFilterType::New();
 385         relPosFilter->VerboseStepFlagOff();
 386         relPosFilter->WriteStepFlagOff();
 387         relPosFilter->SetBackgroundValue(GetBackgroundValue());
 388         relPosFilter->SetInput(slices_support[i]);
 389         relPosFilter->SetInputObject(object);
 390         relPosFilter->AddOrientationTypeString("L");
 391         relPosFilter->SetInverseOrientationFlag(true);
 392         //      relPosFilter->SetIntermediateSpacing(3);
 393         relPosFilter->SetIntermediateSpacingFlag(false);
 394         relPosFilter->SetFuzzyThreshold(0.7);
 395         relPosFilter->AutoCropFlagOff(); // important ! because we join the slices after this loop
 396         relPosFilter->Update();
 397         slices_support[i] = relPosFilter->GetOutput();
 398
 399         // Relative Position : not Anterior
 400         relPosFilter = RelPosFilterType::New();
 401         relPosFilter->VerboseStepFlagOff();
 402         relPosFilter->WriteStepFlagOff();
 403         relPosFilter->SetBackgroundValue(GetBackgroundValue());
 404         relPosFilter->SetInput(slices_support[i]);
 405         relPosFilter->SetInputObject(object);
 406         relPosFilter->AddOrientationTypeString("A");
 407         relPosFilter->SetInverseOrientationFlag(true);
 408         //      relPosFilter->SetIntermediateSpacing(3);
 409         relPosFilter->SetIntermediateSpacingFlag(false);
 410         relPosFilter->SetFuzzyThreshold(0.7);
 411         relPosFilter->AutoCropFlagOff(); // important ! because we join the slices after this loop
 412         relPosFilter->Update();
 413         slices_support[i] = relPosFilter->GetOutput();
 414
 415       }
 416
 417
 418     } // if centroids.size > 1
 419   } // End loop slices
 420
 421   // Joint slices
 422   m_Working_Support = clitk::JoinSlices<MaskImageType>(slices_support, m_Working_Support, 2);
 423
 424   // Save list points
 425   clitk::WriteListOfLandmarks<MaskImageType>(points, "subcl-lr.txt");
 426
 427
 428   StopCurrentStep<MaskImageType>(m_Working_Support);
 429   m_ListOfStations["3P"] = m_Working_Support;
 430 }
 431 //--------------------------------------------------------------------
 432
 433 //--------------------------------------------------------------------
 434 template <class ImageType>
 435 void
 436 clitk::ExtractLymphStationsFilter<ImageType>::
 437 ExtractStation_3P_LR_inf_Limits()
 438 {
 439   /*
 440     "On the right side, the limit is deﬁned by the air–soft-tissue
 441     interface. On the left side, it is deﬁned by the air–tissue
 442     interface superiorly (Fig. 2A–C) and the aorta inferiorly
 443     (Figs. 2D–I and 3A–C)."
 444
 445     inf : not Right to Azygousvein
 446   */
 447   StartNewStep("[Station 3P] Left/Right limits (inferior part) ");
 448
 449   // Load structures
 450   MaskImagePointer AzygousVein = GetAFDB()->template GetImage<MaskImageType>("AzygousVein");
 451
 452   typedef clitk::AddRelativePositionConstraintToLabelImageFilter<MaskImageType> RelPosFilterType;
 453   typename RelPosFilterType::Pointer relPosFilter = RelPosFilterType::New();
 454   relPosFilter->VerboseStepFlagOff();
 455   relPosFilter->WriteStepFlagOff();
 456   relPosFilter->SetBackgroundValue(GetBackgroundValue());
 457   relPosFilter->SetInput(m_Working_Support);
 458   relPosFilter->SetInputObject(AzygousVein);
 459   relPosFilter->AddOrientationTypeString("R");
 460   relPosFilter->SetInverseOrientationFlag(true);
 461   //      relPosFilter->SetIntermediateSpacing(3);
 462   relPosFilter->SetIntermediateSpacingFlag(false);
 463   relPosFilter->SetFuzzyThreshold(0.7);
 464   relPosFilter->AutoCropFlagOn();
 465   relPosFilter->Update();
 466   m_Working_Support = relPosFilter->GetOutput();
 467
 468   StopCurrentStep<MaskImageType>(m_Working_Support);
 469   m_ListOfStations["3P"] = m_Working_Support;
 470 }
 471 //--------------------------------------------------------------------