From bf4928c59a1d39f53fe03deb4b73ecb7e1cf214b Mon Sep 17 00:00:00 2001 From: dsarrut Date: Tue, 22 Feb 2011 07:14:00 +0000 Subject: [PATCH] Station8 almost ok. --- .../clitkAnatomicalFeatureDatabase.cxx | 78 +++++---- segmentation/clitkAnatomicalFeatureDatabase.h | 3 +- .../clitkAnatomicalFeatureDatabase.txx | 27 +-- .../clitkConnectedComponentLabeling.cxx | 6 + ...onnectedComponentLabelingGenericFilter.txx | 82 +++++---- segmentation/clitkExtractLymphStation_3P.txx | 55 ++++++ segmentation/clitkExtractLymphStation_8.txx | 161 ++++++++++-------- segmentation/clitkExtractLymphStations.ggo | 1 + .../clitkExtractLymphStationsFilter.h | 15 +- .../clitkExtractLymphStationsFilter.txx | 121 ++++++++++--- ...clitkExtractLymphStationsGenericFilter.txx | 3 + ...rWithAnatomicalFeatureDatabaseManagement.h | 1 + .../clitkRegionGrowingGenericFilter.txx | 10 -- 13 files changed, 371 insertions(+), 192 deletions(-) create mode 100644 segmentation/clitkExtractLymphStation_3P.txx diff --git a/segmentation/clitkAnatomicalFeatureDatabase.cxx b/segmentation/clitkAnatomicalFeatureDatabase.cxx index 8bf17fb..a23deb4 100644 --- a/segmentation/clitkAnatomicalFeatureDatabase.cxx +++ b/segmentation/clitkAnatomicalFeatureDatabase.cxx @@ -114,46 +114,29 @@ double clitk::AnatomicalFeatureDatabase::GetPoint3D(std::string tag, int dim) //-------------------------------------------------------------------- void clitk::AnatomicalFeatureDatabase::GetPoint3D(std::string tag, PointType3D & p) { - if (m_MapOfTag.find(tag) == m_MapOfTag.end()) { + if (!TagExist(tag)) { clitkExceptionMacro("Could not find the tag <" << tag << "> of type Point3D in the DB"); + return; } - else { - std::string s = m_MapOfTag[tag]; - - // construct a stream from the string - std::stringstream strstr(s); - - // use stream iterators to copy the stream to the vector as - // whitespace separated strings - std::istream_iterator it(strstr); - std::istream_iterator end; - std::vector results(it, end); - // parse the string into 3 doubles - for(int i=0; i<3; i++) { + std::string s = m_MapOfTag[tag]; + + // construct a stream from the string + std::stringstream strstr(s); - if (!clitk::fromString(p[i], results[i].c_str())) { - clitkExceptionMacro("Error while reading Point3D, could not convert '" - << results[i].c_str() << "' into double."); - } + // use stream iterators to copy the stream to the vector as + // whitespace separated strings + std::istream_iterator it(strstr); + std::istream_iterator end; + std::vector results(it, end); - // p[i] = atof(results[i].c_str()); - } + // parse the string into 3 doubles + for(int i=0; i<3; i++) { - /* - // boost - #include - #include - // parse the string into 3 doubles - boost::char_separator sep(", "); - boost::tokenizer > tokens(s, sep); - int i=0; - BOOST_FOREACH(std::string t, tokens) { - std::cout << t << "." << std::endl; - p[i] = atof(t.c_str()); - i++; + if (!clitk::fromString(p[i], results[i].c_str())) { + clitkExceptionMacro("Error while reading Point3D, could not convert '" + << results[i].c_str() << "' into double."); } - */ } } //-------------------------------------------------------------------- @@ -167,4 +150,33 @@ void clitk::AnatomicalFeatureDatabase::SetImageFilename(std::string tag, std::st //-------------------------------------------------------------------- +//-------------------------------------------------------------------- +bool clitk::AnatomicalFeatureDatabase::TagExist(std::string tag) +{ + return (m_MapOfTag.find(tag) != m_MapOfTag.end()); +} +//-------------------------------------------------------------------- + +//-------------------------------------------------------------------- +void clitk::AnatomicalFeatureDatabase::SetDouble(std::string tag, double value) +{ + m_MapOfTag[tag] = clitk::toString(value); +} +//-------------------------------------------------------------------- + +//-------------------------------------------------------------------- +double clitk::AnatomicalFeatureDatabase::GetDouble(std::string tag) +{ + if (!TagExist(tag)) { + clitkExceptionMacro("Could not find the tag <" << tag << "> of type Double in the DB"); + return -1; + } + double a; + if (!clitk::fromString(a, m_MapOfTag[tag])) { + clitkExceptionMacro("Error while reading Double (tag='" << tag << "'), could not convert '" + << m_MapOfTag[tag] << "' into double."); + } + return a; +} +//-------------------------------------------------------------------- diff --git a/segmentation/clitkAnatomicalFeatureDatabase.h b/segmentation/clitkAnatomicalFeatureDatabase.h index e2850e1..e359192 100644 --- a/segmentation/clitkAnatomicalFeatureDatabase.h +++ b/segmentation/clitkAnatomicalFeatureDatabase.h @@ -51,11 +51,12 @@ namespace clitk { void SetPoint3D(TagType tag, PointType3D & p); void GetPoint3D(TagType tag, PointType3D & p); double GetPoint3D(std::string tag, int dim); + bool TagExist(std::string tag); // Set Get image void SetImageFilename(TagType tag, std::string f); template - typename ImageType::Pointer GetImage(TagType tag); + typename ImageType::Pointer GetImage(TagType tag, bool reload=false); template void SetImage(TagType tag, std::string f, diff --git a/segmentation/clitkAnatomicalFeatureDatabase.txx b/segmentation/clitkAnatomicalFeatureDatabase.txx index f1ee25e..7ec7f97 100644 --- a/segmentation/clitkAnatomicalFeatureDatabase.txx +++ b/segmentation/clitkAnatomicalFeatureDatabase.txx @@ -20,7 +20,7 @@ //-------------------------------------------------------------------- template typename ImageType::Pointer AnatomicalFeatureDatabase:: -GetImage(std::string tag) +GetImage(std::string tag, bool reload) { if (m_MapOfTag.find(tag) == m_MapOfTag.end()) { clitkExceptionMacro("Could not find the tag <" << tag << "> of type 'Image' in the DB ('" @@ -28,7 +28,7 @@ GetImage(std::string tag) } else { typename ImageType::Pointer image; - if (m_MapOfImage[tag]) { + if ((!reload) && (m_MapOfImage[tag])) { image = static_cast(m_MapOfImage[tag]); } else { @@ -37,7 +37,8 @@ GetImage(std::string tag) image = readImage(s); // I add a reference count because the cache is not a smartpointer image->SetReferenceCount(image->GetReferenceCount()+1); - // Insert into the cache + // Insert into the cache + m_MapOfImage.erase(tag); m_MapOfImage[tag] = &(*image); // pointer } return image; @@ -69,22 +70,10 @@ ReleaseImage(std::string tag) clitkExceptionMacro("Could not find the tag <" << tag << "> of type Image Filename in the DB"); } else { - DD("TODO"); - exit(0); - if (m_MapOfImage[tag]) { - DD(m_MapOfImage[tag]->GetReferenceCount()); - ImageType * image = static_cast(m_MapOfImage[tag]); - image->SetReferenceCount(image->GetReferenceCount()-1); - m_MapOfImage.erase(tag); - /* - DD(image->GetReferenceCount()); - image->Delete(); - */ - // DD(image->GetReferenceCount()); - } - else { - // Do nothing in this case (image not loaded) - } + typename ImageType::Pointer image = GetImage(tag); + DD(image->GetReferenceCount()); + image->SetReferenceCount(image->GetReferenceCount()-1); + m_MapOfImage.erase(tag); } } //-------------------------------------------------------------------- diff --git a/segmentation/clitkConnectedComponentLabeling.cxx b/segmentation/clitkConnectedComponentLabeling.cxx index bc2fcb1..31d63f4 100644 --- a/segmentation/clitkConnectedComponentLabeling.cxx +++ b/segmentation/clitkConnectedComponentLabeling.cxx @@ -19,26 +19,32 @@ // clitk #include "clitkConnectedComponentLabeling_ggo.h" #include "clitkConnectedComponentLabelingGenericFilter.h" +#include "clitkMemoryUsage.h" //-------------------------------------------------------------------- int main(int argc, char * argv[]) { + clitk::PrintMemory(true, "start"); // Init command line GGO(clitkConnectedComponentLabeling, args_info); CLITK_INIT; // Filter + clitk::PrintMemory(true, "before filter"); typedef clitk::ConnectedComponentLabelingGenericFilter FilterType; FilterType::Pointer filter = FilterType::New(); filter->SetArgsInfo(args_info); + clitk::PrintMemory(true, "before update"); try { filter->Update(); } catch(std::runtime_error e) { std::cout << e.what() << std::endl; } + clitk::PrintMemory(true, "after filter"); + return EXIT_SUCCESS; } // This is the end, my friend //-------------------------------------------------------------------- diff --git a/segmentation/clitkConnectedComponentLabelingGenericFilter.txx b/segmentation/clitkConnectedComponentLabelingGenericFilter.txx index caf5cdf..33f4f3b 100644 --- a/segmentation/clitkConnectedComponentLabelingGenericFilter.txx +++ b/segmentation/clitkConnectedComponentLabelingGenericFilter.txx @@ -21,6 +21,7 @@ // clitk #include "clitkImageCommon.h" +#include "clitkSegmentationUtils.h" // itk #include "itkConnectedComponentImageFilter.h" @@ -78,42 +79,57 @@ void clitk::ConnectedComponentLabelingGenericFilter::UpdateWithInp // Output image type typedef itk::Image OutputImageType; - - // Create CCL filter - DD("CCL"); - typedef itk::ConnectedComponentImageFilter ConnectFilterType; - typename ConnectFilterType::Pointer connectFilter = ConnectFilterType::New(); - connectFilter->SetInput(input); - connectFilter->SetBackgroundValue(mArgsInfo.inputBG_arg); - connectFilter->SetFullyConnected(mArgsInfo.full_flag); - - // TODO SetBackgroud to zero forr relabel ? - - - // Sort by size and remove too small area. - typedef itk::RelabelComponentImageFilter RelabelFilterType; - typename RelabelFilterType::Pointer relabelFilter = RelabelFilterType::New(); - // relabelFilter->InPlaceOn(); - relabelFilter->SetInput(connectFilter->GetOutput()); - relabelFilter->SetMinimumObjectSize(mArgsInfo.minSize_arg); - relabelFilter->Update(); - - DD(mArgsInfo.inputBG_arg); - DD(mArgsInfo.full_flag); - DD(mArgsInfo.minSize_arg); - - // Set information - const std::vector & a = relabelFilter->GetSizeOfObjectsInPixels(); - m_SizeOfObjectsInPixels.resize(a.size()); - for(unsigned int i=0; iGetSizeOfObjectsInPhysicalUnits(); - m_OriginalNumberOfObjects = relabelFilter->GetOriginalNumberOfObjects(); - DD(m_OriginalNumberOfObjects); - DD(m_SizeOfObjectsInPhysicalUnits.size()); + PrintMemory(true, "initial"); + + typename OutputImageType::Pointer output; + { + typename OutputImageType::Pointer temp; + { + // Create CCL filter + typedef itk::ConnectedComponentImageFilter ConnectFilterType; + typename ConnectFilterType::Pointer connectFilter = ConnectFilterType::New(); + // connectFilter->ReleaseDataFlagOn(); // release earlier + connectFilter->SetInput(input); + connectFilter->SetBackgroundValue(mArgsInfo.inputBG_arg); + connectFilter->SetFullyConnected(mArgsInfo.full_flag); + // connectFilter->SetNumberOfThreads(8); + connectFilter->Update(); + temp = connectFilter->GetOutput(); + PrintMemory(true, "after udpate"); + } + PrintMemory(true, "after CCL block"); + DD(input->GetReferenceCount()); + DD(temp->GetReferenceCount()); + + // Sort by size and remove too small area. + typedef itk::RelabelComponentImageFilter RelabelFilterType; + typename RelabelFilterType::Pointer relabelFilter = RelabelFilterType::New(); + // relabelFilter->SetInput(connectFilter->GetOutput()); + relabelFilter->SetInput(temp); + relabelFilter->SetMinimumObjectSize(mArgsInfo.minSize_arg); + relabelFilter->Update(); + + DD(mArgsInfo.inputBG_arg); + DD(mArgsInfo.full_flag); + DD(mArgsInfo.minSize_arg); + + // Set information + const std::vector & a + = relabelFilter->GetSizeOfObjectsInPixels(); + m_SizeOfObjectsInPixels.resize(a.size()); + for(unsigned int i=0; iGetSizeOfObjectsInPhysicalUnits(); + m_OriginalNumberOfObjects = relabelFilter->GetOriginalNumberOfObjects(); + DD(m_OriginalNumberOfObjects); + DD(m_SizeOfObjectsInPhysicalUnits.size()); + + output = relabelFilter->GetOutput(); + } + PrintMemory(true, "after block"); // Write/Save results - typename OutputImageType::Pointer output = relabelFilter->GetOutput(); this->template SetNextOutput(output); + PrintMemory(true, "end filter "); } //-------------------------------------------------------------------- diff --git a/segmentation/clitkExtractLymphStation_3P.txx b/segmentation/clitkExtractLymphStation_3P.txx new file mode 100644 index 0000000..effcd4f --- /dev/null +++ b/segmentation/clitkExtractLymphStation_3P.txx @@ -0,0 +1,55 @@ + +#include +#include + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +ExtractStation_3P_SetDefaultValues() +{ +} +//-------------------------------------------------------------------- + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +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 + MaskImagePointer Lungs = GetAFDB()->template GetImage("Lungs"); + DD("lung ok"); + MaskImagePointType p; + bool found = clitk::FindExtremaPointInAGivenDirection(Lungs, + GetBackgroundValue(), + 1, true, p); + DD(found); + DD(p); + 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, + GetBackgroundValue()); + + StopCurrentStep(m_Working_Support); + m_ListOfStations["3P"] = m_Working_Support; +} +//-------------------------------------------------------------------- diff --git a/segmentation/clitkExtractLymphStation_8.txx b/segmentation/clitkExtractLymphStation_8.txx index 69e24c9..8349261 100644 --- a/segmentation/clitkExtractLymphStation_8.txx +++ b/segmentation/clitkExtractLymphStation_8.txx @@ -2,6 +2,22 @@ #include #include +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +ExtractStation_8_SetDefaultValues() +{ + SetDistanceMaxToAnteriorPartOfTheSpine(10); + MaskImagePointType p; + p[0] = 15; p[1] = 2; p[2] = 1; + SetEsophagusDiltationForAnt(p); + p[0] = 5; p[1] = 10; p[2] = 1; + SetEsophagusDiltationForRight(p); + SetFuzzyThresholdForS8(0.5); +} +//-------------------------------------------------------------------- + //-------------------------------------------------------------------- template void @@ -29,18 +45,17 @@ ExtractStation_8_SI_Limits() // Get Carina Z position MaskImagePointer Carina = GetAFDB()->template GetImage("Carina"); - clitk::PrintMemory(GetVerboseMemoryFlag(), "after read Carina"); + + // CHANGE TO -> OriginOfRightMiddleLobeBronchus ??? + std::vector centroids; clitk::ComputeCentroids(Carina, GetBackgroundValue(), centroids); - DD(centroids[1]); m_CarinaZ = centroids[1][2]; - DD(m_CarinaZ); // add one slice to include carina ? m_CarinaZ += m_Mediastinum->GetSpacing()[2]; - DD(m_CarinaZ); // We dont need Carina structure from now Carina->Delete(); - clitk::PrintMemory(GetVerboseMemoryFlag(), "after delete Carina"); + GetAFDB()->SetDouble("CarinaZ", m_CarinaZ); // Get left lower lobe bronchus (L), take the upper border // Get right bronchus intermedius (RML), take the lower border @@ -55,10 +70,10 @@ ExtractStation_8_SI_Limits() MaskImagePointer Lungs = GetAFDB()->template GetImage("Lungs"); // It should be already croped, so I took the origin and add 10mm above m_DiaphragmInferiorLimit = Lungs->GetOrigin()[2]+10; - DD(m_DiaphragmInferiorLimit); - clitk::PrintMemory(GetVerboseMemoryFlag(), "after read Lung"); - Lungs->Delete(); // we don't need it, release memory - clitk::PrintMemory(GetVerboseMemoryFlag(), "after release Lung"); + // Lungs->Delete(); // we don't need it, release memory -> it we want to release, also free in AFDB + clitk::PrintMemory(GetVerboseMemoryFlag(), "after reading lungs"); + GetAFDB()->template ReleaseImage("Lungs"); + clitk::PrintMemory(GetVerboseMemoryFlag(), "after release lungs"); /* Crop support : Superior limit = carina @@ -228,12 +243,10 @@ ExtractStation_8_Ant_Limits() // Ant limit from carina (start) to end of S7 = originRMLB // Get Trachea MaskImagePointer Trachea = GetAFDB()->template GetImage("Trachea"); - clitk::PrintMemory(GetVerboseMemoryFlag(), "after read Trachea"); MaskImagePointer m_Working_Trachea = clitk::CropImageAbove(Trachea, 2, m_CarinaZ, true, // AutoCrop GetBackgroundValue()); - clitk::PrintMemory(GetVerboseMemoryFlag(), "after crop"); // Seprate into two main bronchi MaskImagePointer LeftBronchus; @@ -241,7 +254,6 @@ ExtractStation_8_Ant_Limits() // Labelize and consider the two first (main) labels m_Working_Trachea = Labelize(m_Working_Trachea, 0, true, 1); - clitk::PrintMemory(GetVerboseMemoryFlag(), "after labelize"); // Carina position must at the first slice that separate the two // main bronchus (not superiorly). We thus first check that the @@ -297,21 +309,13 @@ ExtractStation_8_Ant_Limits() MaskImagePointer OriginOfRightMiddleLobeBronchus = GetAFDB()->template GetImage("OriginOfRightMiddleLobeBronchus"); - clitk::PrintMemory(GetVerboseMemoryFlag(), "after read OriginOfRightMiddleLobeBronchus"); std::vector centroids; clitk::ComputeCentroids(OriginOfRightMiddleLobeBronchus, GetBackgroundValue(), centroids); - DD(centroids.size()); - DD(centroids[0]); // BG - DD(centroids[1]); m_OriginOfRightMiddleLobeBronchusZ = centroids[1][2]; - DD(m_OriginOfRightMiddleLobeBronchusZ); // add one slice to include carina ? m_OriginOfRightMiddleLobeBronchusZ += LeftBronchus->GetSpacing()[2]; - DD(m_OriginOfRightMiddleLobeBronchusZ); - DD(OriginOfRightMiddleLobeBronchus->GetReferenceCount()); // We dont need Carina structure from now OriginOfRightMiddleLobeBronchus->Delete(); - clitk::PrintMemory(GetVerboseMemoryFlag(), "after delete OriginOfRightMiddleLobeBronchus"); LeftBronchus = clitk::CropImageBelow(LeftBronchus, 2, @@ -359,8 +363,6 @@ ExtractStation_8_Ant_Limits() GetBackgroundValue(), 1, 10); writeImage(m_Working_Support, "after.mhd"); -HERE - // Keep main 3D CCL : m_Working_Support = Labelize(m_Working_Support, 0, false, 10); m_Working_Support = KeepLabels(m_Working_Support, @@ -384,7 +386,7 @@ ExtractStation_8_LR_Limits() { //-------------------------------------------------------------------- - StartNewStep("[Station8] Left and Right limits arround esophagus (below Carina)"); + 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 @@ -395,6 +397,24 @@ ExtractStation_8_LR_Limits() MaskImagePointer Esophagus = GetAFDB()->template GetImage("Esophagus"); clitk::PrintMemory(GetVerboseMemoryFlag(), "after read Esophagus"); + // In images from the original article, Atlas – UM, the oesophagus + //was included in nodal stations 3p and 8. Having said that, in the + //description for station 8, it indicates that “the delineation of + //station 8 is limited to the soft tissues surrounding the + //oesophagus”. In the recent article, The IASLC Lung Cancer Staging + //Project (J Thorac Oncol 4:5, 568-77), the images drawn by + //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(); + // Crop Esophagus : keep only below the OriginOfRightMiddleLobeBronchusZ Esophagus = clitk::CropImageAbove(Esophagus, 2, @@ -409,10 +429,18 @@ ExtractStation_8_LR_Limits() radiusInMM, GetBackgroundValue(), GetForegroundValue(), true); - clitk::PrintMemory(GetVerboseMemoryFlag(), "after dilate Esophagus"); writeImage(Esophagus, "enlarged-eso.mhd"); - // Remove Anterior part according to this dilatated esophagus + // Remove Anterior part according to this dilatated esophagus. Note: + // because we crop Esophagus with ORML, the support will also be + // croped in the same way. Here it is a desired feature. If we dont + // want, use SetIgnoreEmptySliceObject(true) + + // In the new IASCL definition, it is not clear if sup limits is + // around carina or On the right, it is “the lower border of the + // bronchus intermedius”, indicated on the image set as a point + // (“lower border of the bronchus intermedius”) + typedef clitk::SliceBySliceRelativePositionFilter RelPosFilterType; typename RelPosFilterType::Pointer relPosFilter = RelPosFilterType::New(); relPosFilter->VerboseStepFlagOff(); @@ -425,63 +453,43 @@ ExtractStation_8_LR_Limits() relPosFilter->SetIntermediateSpacing(3); relPosFilter->ResampleBeforeRelativePositionFilterOn(); relPosFilter->SetFuzzyThreshold(GetFuzzyThresholdForS8()); - relPosFilter->RemoveObjectFlagOff(); + relPosFilter->RemoveObjectFlagOff(); // Do not exclude here because it is dilated relPosFilter->Update(); m_Working_Support = relPosFilter->GetOutput(); - clitk::PrintMemory(GetVerboseMemoryFlag(), "after SbS Rel pos Post"); // AutoCrop (OR SbS ?) m_Working_Support = clitk::AutoCrop(m_Working_Support, GetBackgroundValue()); - writeImage(m_Working_Support, "step1.4.1.mhd"); - - clitk::PrintMemory(GetVerboseMemoryFlag(), "after autocrop"); + StopCurrentStep(m_Working_Support); + //-------------------------------------------------------------------- + StartNewStep("[Station8] Left and Right limits arround esophagus (below Carina)"); // Estract slices for current support for slice by slice processing std::vector slices; clitk::ExtractSlices(m_Working_Support, 2, slices); - clitk::PrintMemory(GetVerboseMemoryFlag(), "after support slices"); // Estract slices of Esophagus (resize like support before to have the same set of slices) MaskImagePointer EsophagusForSlice = clitk::ResizeImageLike(Esophagus, m_Working_Support, GetBackgroundValue()); - clitk::PrintMemory(GetVerboseMemoryFlag(), "after Eso resize"); std::vector eso_slices; clitk::ExtractSlices(EsophagusForSlice, 2, eso_slices); - clitk::PrintMemory(GetVerboseMemoryFlag(), "after Eso slices"); // Estract slices of Vertebral (resize like support before to have the same set of slices) MaskImagePointer VertebralBody = GetAFDB()->template GetImage("VertebralBody"); - clitk::PrintMemory(GetVerboseMemoryFlag(), "after Read VertebralBody"); VertebralBody = clitk::ResizeImageLike(VertebralBody, m_Working_Support, GetBackgroundValue()); - clitk::PrintMemory(GetVerboseMemoryFlag(), "after VertebralBody Resize"); std::vector vert_slices; clitk::ExtractSlices(VertebralBody, 2, vert_slices); - clitk::PrintMemory(GetVerboseMemoryFlag(), "after VertebralBody slices"); // Estract slices of Aorta (resize like support before to have the same set of slices) MaskImagePointer Aorta = GetAFDB()->template GetImage("Aorta"); - clitk::PrintMemory(GetVerboseMemoryFlag(), "after Read Aorta"); Aorta = clitk::ResizeImageLike(Aorta, m_Working_Support, GetBackgroundValue()); - clitk::PrintMemory(GetVerboseMemoryFlag(), "after Aorta Resize"); std::vector aorta_slices; clitk::ExtractSlices(Aorta, 2, aorta_slices); - clitk::PrintMemory(GetVerboseMemoryFlag(), "after Aorta slices"); // Extract slices of Mediastinum (resize like support before to have the same set of slices) m_Mediastinum = GetAFDB()->template GetImage("Mediastinum"); - clitk::PrintMemory(GetVerboseMemoryFlag(), "after read Mediastinum"); m_Mediastinum = clitk::ResizeImageLike(m_Mediastinum, m_Working_Support, GetBackgroundValue()); - clitk::PrintMemory(GetVerboseMemoryFlag(), "after resize Mediastinum"); std::vector mediast_slices; clitk::ExtractSlices(m_Mediastinum, 2, mediast_slices); - clitk::PrintMemory(GetVerboseMemoryFlag(), "after Mediastinum slices"); - - writeImage(EsophagusForSlice, "slices_eso.mhd"); - writeImage(VertebralBody, "slices_vert.mhd"); - writeImage(Aorta, "slices_aorta.mhd"); - writeImage(m_Mediastinum, "slices_medias.mhd"); - writeImage(m_Working_Support, "slices_support.mhd"); - // List of points std::vector p_RightMostAnt; @@ -516,17 +524,6 @@ ExtractStation_8_LR_Limits() // Right is at left on screen, coordinate decrease // Left is at right on screen, coordinate increase - // Find right limit of Esophagus and Aorta - clitk::FindExtremaPointInAGivenDirection(eso_slices[i], GetBackgroundValue(), 0, true, sp_maxRight_Eso); - 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 'wall' - p_AllPoints.push_back(p); - p_AllPoints.push_back(pp); - if (p[0] only at most Post part of current // slice support. First found most ant point in VertebralBody typedef MaskSliceType SliceType; @@ -537,7 +534,6 @@ ExtractStation_8_LR_Limits() // 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; - DD(i); int j=i++; bool found = false; while (!found) { @@ -545,7 +541,6 @@ ExtractStation_8_LR_Limits() //clitkExceptionMacro("No foreground pixels in this VertebralBody slices ??"); j++; } - DD(j); } p_slice_ant[1] += GetDistanceMaxToAnteriorPartOfTheSpine(); // Consider offset @@ -584,6 +579,7 @@ ExtractStation_8_LR_Limits() else { // in that case, we are probably below the diaphragm, so we // add aribtrarly few mm sp_maxRight_Vertebra[0] -= 2; // Leave 2 mm around the VertebralBody + clitk::PointsUtils::Convert2DTo3D(sp_maxRight_Vertebra, m_Mediastinum, i, p); } p_RightMostPost.push_back(p); p_AllPoints.push_back(p); @@ -602,24 +598,54 @@ ExtractStation_8_LR_Limits() else { // in that case, we are probably below the diaphragm, so we // add aribtrarly few mm sp_maxLeft_Vertebra[0] += 2; // Leave 2 mm around the VertebralBody + clitk::PointsUtils::Convert2DTo3D(sp_maxLeft_Vertebra, m_Mediastinum, i, p); } p_LeftMostPost.push_back(p); p_AllPoints.push_back(p); + // 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); + 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' + p_AllPoints.push_back(p); + if (f) { + p_AllPoints.push_back(pp); + MaskImagePointType A = p_RightMostPost.back(); + 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_RightMostAnt.push_back(p); + else p_RightMostAnt.push_back(pp); + } + else { + p_RightMostAnt.push_back(p); + } + // -------------------------------------------------------------------------- // Find the limit on the Left: most left point between Eso and // Vertebra. (Left is left on screen, coordinate increase) // Find left limit of Esophagus clitk::FindExtremaPointInAGivenDirection(eso_slices[i], GetBackgroundValue(), 0, false, sp_maxLeft_Eso); - clitk::FindExtremaPointInAGivenDirection(aorta_slices[i], GetBackgroundValue(), 0, false, sp_maxLeft_Aorta); + f = clitk::FindExtremaPointInAGivenDirection(aorta_slices[i], GetBackgroundValue(), 0, false, sp_maxLeft_Aorta); clitk::PointsUtils::Convert2DTo3D(sp_maxLeft_Eso, EsophagusForSlice, i, p); clitk::PointsUtils::Convert2DTo3D(sp_maxLeft_Aorta, Aorta, i, pp); p_AllPoints.push_back(p); pp[0] += 2; // Add a margin of 2 mm to include the 'wall' - p_AllPoints.push_back(pp); - if (p[0]>pp[0]) p_LeftMostAnt.push_back(p); // Insert point most at right - else p_LeftMostAnt.push_back(pp); + if (f) { // not below Aorta + p_AllPoints.push_back(pp); + MaskImagePointType A = p_LeftMostPost.back(); + 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); + } + else { + p_LeftMostAnt.push_back(p); + } } // End of slice loop clitk::WriteListOfLandmarks(p_AllPoints, "LR-Eso-Vert.txt"); @@ -632,7 +658,8 @@ ExtractStation_8_LR_Limits() clitk::SliceBySliceSetBackgroundFromLineSeparation(m_Working_Support, p_LeftMostAnt, p_LeftMostPost, GetBackgroundValue(), 0, -10); - // DEBUG + // END + StopCurrentStep(m_Working_Support); m_ListOfStations["8"] = m_Working_Support; return; } diff --git a/segmentation/clitkExtractLymphStations.ggo b/segmentation/clitkExtractLymphStations.ggo index bbaa460..4b39015 100644 --- a/segmentation/clitkExtractLymphStations.ggo +++ b/segmentation/clitkExtractLymphStations.ggo @@ -17,6 +17,7 @@ section "I/O" option "afdb" a "Input Anatomical Feature DB" string no option "input" i "Input filename" string no option "output" o "Output lungs mask filename" string no +option "station" - "Force to compute station even if already exist in the DB" string no multiple section "Options for Station 8" option "maxAntSpine" - "Distance max to anterior part of the spine in mm" double no default="10" diff --git a/segmentation/clitkExtractLymphStationsFilter.h b/segmentation/clitkExtractLymphStationsFilter.h index 0a9099b..8c7b329 100644 --- a/segmentation/clitkExtractLymphStationsFilter.h +++ b/segmentation/clitkExtractLymphStationsFilter.h @@ -96,6 +96,9 @@ namespace clitk { itkGetConstMacro(FuzzyThreshold, double); itkSetMacro(Station7Filename, std::string); itkGetConstMacro(Station7Filename, std::string); + + bool GetComputeStation(std::string s); + void AddComputeStation(std::string station) ; protected: ExtractLymphStationsFilter(); @@ -110,7 +113,10 @@ namespace clitk { MaskImagePointer m_Working_Support; std::map m_ListOfStations; MaskImagePixelType m_BackgroundValue; - MaskImagePixelType m_ForegroundValue; + MaskImagePixelType m_ForegroundValue; + std::map m_ComputeStationMap; + + bool CheckForStation(std::string station); // Station 8 double m_DistanceMaxToAnteriorPartOfTheSpine; @@ -122,11 +128,17 @@ namespace clitk { MaskImagePointType m_EsophagusDiltationForRight; MaskImagePointer EnlargeEsophagusDilatationRadiusInferiorly(MaskImagePointer & eso); void ExtractStation_8(); + void ExtractStation_8_SetDefaultValues(); void ExtractStation_8_SI_Limits(); void ExtractStation_8_Post_Limits(); void ExtractStation_8_Ant_Limits(); void ExtractStation_8_LR_Limits(); void ExtractStation_8_LR_Limits_old(); + + // Station 3P + void ExtractStation_3P(); + void ExtractStation_3P_SetDefaultValues(); + void ExtractStation_3P_SI_Limits(); // Station 7 void ExtractStation_7(); @@ -174,6 +186,7 @@ namespace clitk { #ifndef ITK_MANUAL_INSTANTIATION #include "clitkExtractLymphStationsFilter.txx" #include "clitkExtractLymphStation_8.txx" +#include "clitkExtractLymphStation_3P.txx" #include "clitkExtractLymphStation_7.txx" #include "clitkExtractLymphStation_4RL.txx" #endif diff --git a/segmentation/clitkExtractLymphStationsFilter.txx b/segmentation/clitkExtractLymphStationsFilter.txx index 2172dae..70d540b 100644 --- a/segmentation/clitkExtractLymphStationsFilter.txx +++ b/segmentation/clitkExtractLymphStationsFilter.txx @@ -52,14 +52,9 @@ ExtractLymphStationsFilter(): SetBackgroundValue(0); SetForegroundValue(1); - // Station 8 - SetDistanceMaxToAnteriorPartOfTheSpine(10); - MaskImagePointType p; - p[0] = 15; p[1] = 2; p[2] = 1; - SetEsophagusDiltationForAnt(p); - p[0] = 5; p[1] = 10; p[2] = 1; - SetEsophagusDiltationForRight(p); - SetFuzzyThresholdForS8(0.5); + // Default values + ExtractStation_8_SetDefaultValues(); + ExtractStation_3P_SetDefaultValues(); // Station 7 SetFuzzyThreshold(0.5); @@ -84,14 +79,13 @@ GenerateOutputInformation() { ExtractStation_8(); StopSubStep(); - // Compute some interesting points in trachea - // ( ALTERNATIVE -> SKELETON ANALYSIS ? - // Pb : not sufficient for mostXX points ... ) + DD(GetCurrentStepNumber()); - /* ==> todo (but why ???) - ComputeTracheaCentroidsAboveCarina(); - ComputeBronchusExtremaPointsBelowCarina(); - */ + // Extract Station3P + StartNewStep("Station 3P"); + StartSubStep(); + ExtractStation_3P(); + StopSubStep(); if (0) { // temporary suppress // Extract Station7 @@ -149,23 +143,91 @@ GenerateData() { //-------------------------------------------------------------------- template -void +bool clitk::ExtractLymphStationsFilter:: -ExtractStation_8() { +CheckForStation(std::string station) +{ + // Compute Station name + std::string s = "Station"+station; + - // Check if m_ListOfStations["8"] exist. If yes -> use it as initial - // support instead of m_Mediastinum - if (m_ListOfStations["8"]) { - DD("Station 8 support already exist -> use it"); - m_Working_Support = m_ListOfStations["8"]; + // Check if station already exist in DB + bool found = false; + if (GetAFDB()->TagExist(s)) { + m_ListOfStations[station] = GetAFDB()->template GetImage(s); + found = true; } - else m_Working_Support = m_Mediastinum; - ExtractStation_8_SI_Limits(); - ExtractStation_8_Post_Limits(); - ExtractStation_8_Ant_Limits(); - ExtractStation_8_LR_Limits(); - // ExtractStation_8_LR_Limits(); + // Define the starting support + if (found && GetComputeStation("8")) { + std::cout << "Station " << station << " already exists, but re-computation forced." << std::endl; + } + if (!found || GetComputeStation("8")) { + m_Working_Support = m_Mediastinum = GetAFDB()->template GetImage("Mediastinum", true); + return true; + } + else { + std::cout << "Station " << station << " found. I used it" << std::endl; + return false; + } +} +//-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template +bool +clitk::ExtractLymphStationsFilter:: +GetComputeStation(std::string station) +{ + return (m_ComputeStationMap.find(station) != m_ComputeStationMap.end()); +} +//-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +AddComputeStation(std::string station) +{ + m_ComputeStationMap[station] = true; +} +//-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +ExtractStation_8() +{ + if (CheckForStation("8")) { + ExtractStation_8_SI_Limits(); + ExtractStation_8_Post_Limits(); + ExtractStation_8_Ant_Limits(); + ExtractStation_8_LR_Limits(); + // Store image filenames into AFDB + writeImage(m_ListOfStations["8"], "seg/Station8.mhd"); + GetAFDB()->SetImageFilename("Station8", "seg/Station8.mhd"); + WriteAFDB(); + } +} +//-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template +void +clitk::ExtractLymphStationsFilter:: +ExtractStation_3P() +{ + if (CheckForStation("3P")) { + ExtractStation_3P_SI_Limits(); + // Store image filenames into AFDB + writeImage(m_ListOfStations["3P"], "seg/Station3P.mhd"); + GetAFDB()->SetImageFilename("Station3P", "seg/Station3P.mhd"); + } } //-------------------------------------------------------------------- @@ -202,6 +264,9 @@ ExtractStation_4RL() { //-------------------------------------------------------------------- +//-------------------------------------------------------------------- + + //-------------------------------------------------------------------- template void diff --git a/segmentation/clitkExtractLymphStationsGenericFilter.txx b/segmentation/clitkExtractLymphStationsGenericFilter.txx index bc5c936..ae1106a 100644 --- a/segmentation/clitkExtractLymphStationsGenericFilter.txx +++ b/segmentation/clitkExtractLymphStationsGenericFilter.txx @@ -104,6 +104,9 @@ SetOptionsFromArgsInfoToFilter(FilterType * f) } } f->SetEsophagusDiltationForRight(p); + + for(uint i=0; iAddComputeStation(mArgsInfo.station_arg[i]); } //-------------------------------------------------------------------- diff --git a/segmentation/clitkFilterWithAnatomicalFeatureDatabaseManagement.h b/segmentation/clitkFilterWithAnatomicalFeatureDatabaseManagement.h index 9bcd9d9..0691059 100644 --- a/segmentation/clitkFilterWithAnatomicalFeatureDatabaseManagement.h +++ b/segmentation/clitkFilterWithAnatomicalFeatureDatabaseManagement.h @@ -53,6 +53,7 @@ namespace clitk { void WriteAFDB(); void LoadAFDB(); + AnatomicalFeatureDatabase * GetAFDB(); void SetAFDB(AnatomicalFeatureDatabase * a) { m_AFDB = a; } diff --git a/segmentation/clitkRegionGrowingGenericFilter.txx b/segmentation/clitkRegionGrowingGenericFilter.txx index 2ea7678..9ebfe18 100755 --- a/segmentation/clitkRegionGrowingGenericFilter.txx +++ b/segmentation/clitkRegionGrowingGenericFilter.txx @@ -18,16 +18,6 @@ #ifndef clitkRegionGrowingGenericFilter_txx #define clitkRegionGrowingGenericFilter_txx -/* ================================================= - * @file clitkRegionGrowingGenericFilter.txx - * @author - * @date - * - * @brief - * - ===================================================*/ - - namespace clitk { -- 2.45.1