From 5bff4d89c8ca42336267c36974e9ceb6fe2ac843 Mon Sep 17 00:00:00 2001 From: David Sarrut Date: Thu, 27 Oct 2011 15:52:29 +0200 Subject: [PATCH] Second version of RelativePositionAnalyzer --- itk/clitkRelativePositionAnalyzerFilter.h | 49 +++++++- itk/clitkRelativePositionAnalyzerFilter.txx | 119 +++++++++++--------- 2 files changed, 109 insertions(+), 59 deletions(-) diff --git a/itk/clitkRelativePositionAnalyzerFilter.h b/itk/clitkRelativePositionAnalyzerFilter.h index 7d70fe7..3554049 100644 --- a/itk/clitkRelativePositionAnalyzerFilter.h +++ b/itk/clitkRelativePositionAnalyzerFilter.h @@ -23,6 +23,7 @@ #include "clitkFilterWithAnatomicalFeatureDatabaseManagement.h" #include "clitkFilterBase.h" #include "clitkSliceBySliceRelativePositionFilter.h" +#include "clitkRelativePositionDataBase.h" // itk #include @@ -33,15 +34,17 @@ namespace clitk { //-------------------------------------------------------------------- /* Analyze the relative position of a Target (mask image) according - to some Object, in a given Support. Indicate the main important - position of this Target according the Object. + to some Object (mask image), in a given Support (mask + image). Compute the optimal threshold allowing to remove the + maximal area from the Support without removing area belonging to + the Target. */ //-------------------------------------------------------------------- template class RelativePositionAnalyzerFilter: - public virtual FilterBase, - public clitk::FilterWithAnatomicalFeatureDatabaseManagement, + // public virtual clitk::FilterBase, + public clitk::FilterWithAnatomicalFeatureDatabaseManagement, public itk::ImageToImageFilter { @@ -85,8 +88,28 @@ namespace clitk { itkGetConstMacro(ForegroundValue, PixelType); itkSetMacro(ForegroundValue, PixelType); + itkGetConstMacro(NumberOfBins, int); + itkSetMacro(NumberOfBins, int); + + itkGetConstMacro(NumberOfAngles, int); + itkSetMacro(NumberOfAngles, int); + + itkGetConstMacro(AreaLossTolerance, double); + itkSetMacro(AreaLossTolerance, double); + + itkGetConstMacro(SupportSize, int); + itkGetConstMacro(TargetSize, int); + itkGetConstMacro(SizeWithThreshold, int); + itkGetConstMacro(SizeWithReverseThreshold, int); + + std::vector & GetListOfInformation() { return m_ListOfInformation; } + std::vector & GetListOfOrientation() { return m_ListOfOrientation; } + // For debug void PrintOptions(); + + // Print output + void Print(std::string s=" ", std::ostream & os=std::cout); // I dont want to verify inputs information virtual void VerifyInputInformation() { } @@ -95,14 +118,28 @@ namespace clitk { RelativePositionAnalyzerFilter(); virtual ~RelativePositionAnalyzerFilter() {} + itkSetMacro(SupportSize, int); + itkSetMacro(TargetSize, int); + itkSetMacro(SizeWithThreshold, int); + itkSetMacro(SizeWithReverseThreshold, int); + PixelType m_BackgroundValue; PixelType m_ForegroundValue; ImagePointer m_Support; ImagePointer m_Object; ImagePointer m_Target; - + int m_NumberOfAngles; + int m_NumberOfBins; + double m_AreaLossTolerance; + int m_SupportSize; + int m_TargetSize; + int m_SizeWithReverseThreshold; + int m_SizeWithThreshold; + std::vector m_ListOfAngles; + std::vector m_ListOfInformation; + std::vector m_ListOfOrientation; + virtual void GenerateOutputInformation(); - virtual void GenerateInputRequestedRegion(); virtual void GenerateData(); typename FloatImageType::Pointer diff --git a/itk/clitkRelativePositionAnalyzerFilter.txx b/itk/clitkRelativePositionAnalyzerFilter.txx index 8b69c68..c8d126a 100644 --- a/itk/clitkRelativePositionAnalyzerFilter.txx +++ b/itk/clitkRelativePositionAnalyzerFilter.txx @@ -20,14 +20,22 @@ template clitk::RelativePositionAnalyzerFilter:: RelativePositionAnalyzerFilter(): - clitk::FilterBase(), - clitk::FilterWithAnatomicalFeatureDatabaseManagement(), + // clitk::FilterBase(), + clitk::FilterWithAnatomicalFeatureDatabaseManagement(), itk::ImageToImageFilter() { this->SetNumberOfRequiredInputs(3); // support, object, target VerboseFlagOff(); SetBackgroundValue(0); SetForegroundValue(1); + SetNumberOfBins(100); + SetNumberOfAngles(4); + SetAreaLossTolerance(0.01); + m_ListOfAngles.clear(); + SetSupportSize(0); + SetTargetSize(0); + SetSizeWithThreshold(0); + SetSizeWithReverseThreshold(0); } //-------------------------------------------------------------------- @@ -92,24 +100,6 @@ GenerateOutputInformation() //-------------------------------------------------------------------- -//-------------------------------------------------------------------- -template -void -clitk::RelativePositionAnalyzerFilter:: -GenerateInputRequestedRegion() -{ - // Call default - itk::ImageToImageFilter::GenerateInputRequestedRegion(); - // Get input pointers and set requested region to common region - ImagePointer input1 = dynamic_cast(itk::ProcessObject::GetInput(0)); - ImagePointer input2 = dynamic_cast(itk::ProcessObject::GetInput(1)); - ImagePointer input3 = dynamic_cast(itk::ProcessObject::GetInput(2)); - input1->SetRequestedRegion(input1->GetLargestPossibleRegion()); - input2->SetRequestedRegion(input2->GetLargestPossibleRegion()); - input3->SetRequestedRegion(input3->GetLargestPossibleRegion()); -} -//-------------------------------------------------------------------- - //-------------------------------------------------------------------- template void @@ -117,7 +107,7 @@ clitk::RelativePositionAnalyzerFilter:: GenerateData() { this->LoadAFDB(); - + // Get input pointer m_Support = dynamic_cast(itk::ProcessObject::GetInput(0)); m_Object = dynamic_cast(itk::ProcessObject::GetInput(1)); @@ -138,59 +128,66 @@ GenerateData() statFilter->SetInput(m_Support); statFilter->SetLabelInput(m_Support); statFilter->Update(); - int m_SupportSize = statFilter->GetCount(GetForegroundValue()); - // DD(m_SupportSize); + SetSupportSize(statFilter->GetCount(GetForegroundValue())); + // DD(GetSupportSize()); // Compute the initial target size ImagePointer s = clitk::ResizeImageLike(m_Support, m_Target, GetBackgroundValue()); statFilter->SetInput(s); statFilter->SetLabelInput(m_Target); statFilter->Update(); - int m_TargetSize = statFilter->GetCount(GetForegroundValue()); - // DD(m_TargetSize); + SetTargetSize(statFilter->GetCount(GetForegroundValue())); + // DD(GetTargetSize()); // Build the list of tested orientations - std::vector m_ListOfAngles; - int m_NumberOfAngles = this->GetAFDB()->GetDouble("NumberOfAngles"); - for(uint i=0; i180) a = 180-a; m_ListOfAngles.push_back(clitk::deg2rad(a)); + RelativePositionOrientationType r; + r.angle1 = clitk::deg2rad(a); + r.angle2 = 0; + r.notFlag = false; + m_ListOfOrientation.push_back(r); + r.notFlag = true; + m_ListOfOrientation.push_back(r); } // Loop on all orientations + int bins = GetNumberOfBins(); + double tolerance = GetAreaLossTolerance(); for(int i=0; i(map, "fuzzy_"+toString(i)+".mha"); // Compute the optimal thresholds (direct and inverse) - double mThreshold; - double mReverseThreshold; - int bins = this->GetAFDB()->GetDouble("bins"); - double tolerance = this->GetAFDB()->GetDouble("TargetAreaLossTolerance"); + double mThreshold=0.0; + double mReverseThreshold=1.0; ComputeOptimalThresholds(map, m_Target, bins, tolerance, mThreshold, mReverseThreshold); // Use the threshold to compute new support - int s1; + int s1 = GetSupportSize(); + // DD(mThreshold); + // DD(mReverseThreshold); if (mThreshold > 0.0) { ImagePointer support1 = clitk::SliceBySliceRelativePosition(m_Support, m_Object, 2, mThreshold, m_ListOfAngles[i],false, false, -1, true, false); - writeImage(support1, "sup_"+toString(i)+".mha"); + // writeImage(support1, "sup_"+toString(i)+".mha"); // Compute the new support size statFilter->SetInput(support1); statFilter->SetLabelInput(support1); statFilter->Update(); s1 = statFilter->GetCount(GetForegroundValue()); } - else s1 = m_SupportSize; - int s2; + int s2 = GetSupportSize(); if (mReverseThreshold < 1.0) { + // DD(m_ListOfAngles[1]); ImagePointer support2 = clitk::SliceBySliceRelativePosition(m_Support, m_Object, 2, mReverseThreshold, @@ -204,20 +201,20 @@ GenerateData() statFilter->Update(); s2 = statFilter->GetCount(GetForegroundValue()); } - else s2 =m_SupportSize; - - // Print results - std::cout << i << " " << clitk::rad2deg(m_ListOfAngles[i]) << "\t" - << m_SupportSize << " " << m_TargetSize << "\t" - << s1/(double)m_SupportSize << " " << s2/(double)m_SupportSize << "\t" - << mThreshold << " " << mReverseThreshold << std::endl; - - } // end loop on orientations - - // Final Step -> set output TODO - // this->SetNthOutput(0, working_image); - // this->GraftOutput(working_image); + // Set results values + RelativePositionInformationType r; + r.threshold = mThreshold; + r.sizeAfterThreshold = s1; // DD(s1); + r.sizeBeforeThreshold = GetSupportSize(); + r.sizeReference = GetTargetSize(); + m_ListOfInformation.push_back(r); + + r.threshold = mReverseThreshold; + r.sizeAfterThreshold = s2; // DD(s2); + m_ListOfInformation.push_back(r); + // Print(); + } // end loop on orientations } //-------------------------------------------------------------------- @@ -292,7 +289,8 @@ ComputeOptimalThresholds(FloatImageType * map, ImageType * target, int bins, dou for(int j=0; jGetFrequency(j)/(double)count); if ((!found) && (sum > tolerance)) { - threshold = h->GetBinMin(0,j); + if (j==0) threshold = h->GetBinMin(0,j); + else threshold = h->GetBinMin(0,j-1); // the last before reaching the threshold found = true; } } @@ -304,10 +302,25 @@ ComputeOptimalThresholds(FloatImageType * map, ImageType * target, int bins, dou for(int j=bins-1; j>=0; j--) { sum += ((double)h->GetFrequency(j)/(double)count); if ((!found) && (sum > tolerance)) { - reverseThreshold = h->GetBinMax(0,j); + if (j==bins-1) reverseThreshold = h->GetBinMax(0,j); + else reverseThreshold = h->GetBinMax(0,j+1); found = true; } } } //-------------------------------------------------------------------- + +//-------------------------------------------------------------------- +template +void +clitk::RelativePositionAnalyzerFilter:: +Print(std::string s, std::ostream & os) +{ + for(int i=0; i