template <class ImageType>
clitk::RelativePositionAnalyzerFilter<ImageType>::
RelativePositionAnalyzerFilter():
- clitk::FilterBase(),
- clitk::FilterWithAnatomicalFeatureDatabaseManagement(),
itk::ImageToImageFilter<ImageType, ImageType>()
{
- this->SetNumberOfRequiredInputs(3); // support, object, target
- VerboseFlagOff();
+ this->SetNumberOfRequiredInputs(3); // Input : support, object, target
SetBackgroundValue(0);
SetForegroundValue(1);
+ SetNumberOfBins(100);
+ SetAreaLossTolerance(0.01);
+ SetSupportSize(0);
+ SetTargetSize(0);
+ SetSizeWithThreshold(0);
+ SetSizeWithReverseThreshold(0);
}
//--------------------------------------------------------------------
//--------------------------------------------------------------------
-//--------------------------------------------------------------------
-template <class ImageType>
-void
-clitk::RelativePositionAnalyzerFilter<ImageType>::
-GenerateInputRequestedRegion()
-{
- // Call default
- itk::ImageToImageFilter<ImageType, ImageType>::GenerateInputRequestedRegion();
- // Get input pointers and set requested region to common region
- ImagePointer input1 = dynamic_cast<ImageType*>(itk::ProcessObject::GetInput(0));
- ImagePointer input2 = dynamic_cast<ImageType*>(itk::ProcessObject::GetInput(1));
- ImagePointer input3 = dynamic_cast<ImageType*>(itk::ProcessObject::GetInput(2));
- input1->SetRequestedRegion(input1->GetLargestPossibleRegion());
- input2->SetRequestedRegion(input2->GetLargestPossibleRegion());
- input3->SetRequestedRegion(input3->GetLargestPossibleRegion());
-}
-//--------------------------------------------------------------------
-
//--------------------------------------------------------------------
template <class ImageType>
void
clitk::RelativePositionAnalyzerFilter<ImageType>::
GenerateData()
{
- this->LoadAFDB();
-
- // Get input pointer
- m_Support = dynamic_cast<ImageType*>(itk::ProcessObject::GetInput(0));
+ ImagePointer temp = dynamic_cast<ImageType*>(itk::ProcessObject::GetInput(0));
m_Object = dynamic_cast<ImageType*>(itk::ProcessObject::GetInput(1));
m_Target = dynamic_cast<ImageType*>(itk::ProcessObject::GetInput(2));
- static const unsigned int dim = ImageType::ImageDimension;
- // Remove object from support
+ // Remove object from support (keep initial image)
+ m_Support = clitk::Clone<ImageType>(temp);
clitk::AndNot<ImageType>(m_Support, m_Object, GetBackgroundValue());
- // Resize object like target (to enable substraction later)
- ImagePointer objectLikeTarget = clitk::ResizeImageLike<ImageType>(m_Object, m_Target, GetBackgroundValue());
-
// Define filter to compute statics on mask image
typedef itk::LabelStatisticsImageFilter<ImageType, ImageType> StatFilterType;
typename StatFilterType::Pointer statFilter = StatFilterType::New();
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<ImageType>(m_Support, m_Target, GetBackgroundValue());
statFilter->SetInput(s);
statFilter->SetLabelInput(m_Target);
statFilter->Update();
- int m_TargetSize = statFilter->GetCount(GetForegroundValue());
- // DD(m_TargetSize);
-
- // Build the list of tested orientations
- std::vector<double> m_ListOfAngles;
- int m_NumberOfAngles = this->GetAFDB()->GetDouble("NumberOfAngles");
- for(uint i=0; i<m_NumberOfAngles; i++) {
- double a = i*360.0/m_NumberOfAngles;
- if (a>180) a = 180-a;
- m_ListOfAngles.push_back(clitk::deg2rad(a));
+ SetTargetSize(statFilter->GetCount(GetForegroundValue()));
+ // DD(GetTargetSize());
+
+ //
+ int bins = GetNumberOfBins();
+ double tolerance = GetAreaLossTolerance();
+
+ // Compute Fuzzy map
+ double angle = GetDirection().angle1;
+ typename FloatImageType::Pointer map = ComputeFuzzyMap(m_Object, m_Target, m_Support, angle);
+ writeImage<FloatImageType>(map, "fuzzy_"+toString(clitk::rad2deg(angle))+".mha");
+
+ // Compute the optimal thresholds (direct and inverse)
+ 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 = GetSupportSize();
+ if (mThreshold > 0.0) {
+ ImagePointer support1 =
+ clitk::SliceBySliceRelativePosition<ImageType>(m_Support, m_Object, 2,
+ mThreshold,
+ angle,false, // inverseFlag
+ false, // uniqueConnectedComponent
+ -1, true,
+ false);//singleObjectCCL
+ // Compute the new support size
+ statFilter->SetInput(support1);
+ statFilter->SetLabelInput(support1);
+ statFilter->Update();
+ s1 = statFilter->GetCount(GetForegroundValue());
}
-
- // Loop on all orientations
- for(int i=0; i<m_ListOfAngles.size(); i++) {
-
- // Compute Fuzzy map
- typename FloatImageType::Pointer map = ComputeFuzzyMap(objectLikeTarget, m_Target, m_ListOfAngles[i]);
- writeImage<FloatImageType>(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");
- ComputeOptimalThresholds(map, m_Target, bins, tolerance, mThreshold, mReverseThreshold);
-
- // Use the threshold to compute new support
- int s1;
- if (mThreshold > 0.0) {
- ImagePointer support1 =
- clitk::SliceBySliceRelativePosition<ImageType>(m_Support, m_Object, 2,
- mThreshold,
- m_ListOfAngles[i],false,
- false, -1, true, false);
- writeImage<ImageType>(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;
- if (mReverseThreshold < 1.0) {
- ImagePointer support2 =
- clitk::SliceBySliceRelativePosition<ImageType>(m_Support, m_Object, 2,
- mReverseThreshold,
- m_ListOfAngles[i],true,
- false, -1, true, false);
- writeImage<ImageType>(support2, "sup_rev_"+toString(i)+".mha");
- // Compute the new support size
- statFilter = StatFilterType::New();
- statFilter->SetInput(support2);
- statFilter->SetLabelInput(support2);
- 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);
+ int s2 = GetSupportSize();
+ if (mReverseThreshold < 1.0) {
+ ImagePointer support2 =
+ clitk::SliceBySliceRelativePosition<ImageType>(m_Support, m_Object, 2,
+ mReverseThreshold,
+ angle,true,// inverseFlag
+ false, // uniqueConnectedComponent
+ -1, true,
+ false); //singleObjectCCL
+ // Compute the new support size
+ statFilter = StatFilterType::New();
+ statFilter->SetInput(support2);
+ statFilter->SetLabelInput(support2);
+ statFilter->Update();
+ s2 = statFilter->GetCount(GetForegroundValue());
+ }
+
+ // Set results values
+ m_Info.threshold = mThreshold;
+ m_Info.sizeAfterThreshold = s1;
+ m_Info.sizeBeforeThreshold = GetSupportSize();
+ m_Info.sizeReference = GetTargetSize();
+ m_InfoReverse.threshold = mReverseThreshold;
+ m_InfoReverse.sizeAfterThreshold = s2;
+ m_InfoReverse.sizeBeforeThreshold = GetSupportSize();
+ m_InfoReverse.sizeReference = GetTargetSize();
}
//--------------------------------------------------------------------
template <class ImageType>
typename clitk::RelativePositionAnalyzerFilter<ImageType>::FloatImageType::Pointer
clitk::RelativePositionAnalyzerFilter<ImageType>::
-ComputeFuzzyMap(ImageType * object, ImageType * target, double angle)
+ComputeFuzzyMap(ImageType * object, ImageType * target, ImageType * support, double angle)
{
typedef clitk::SliceBySliceRelativePositionFilter<ImageType> SliceRelPosFilterType;
typedef typename SliceRelPosFilterType::FloatImageType FloatImageType;
typename SliceRelPosFilterType::Pointer sliceRelPosFilter = SliceRelPosFilterType::New();
sliceRelPosFilter->VerboseStepFlagOff();
sliceRelPosFilter->WriteStepFlagOff();
- sliceRelPosFilter->SetInput(target);
+ sliceRelPosFilter->SetInput(support);
sliceRelPosFilter->SetInputObject(object);
sliceRelPosFilter->SetDirection(2);
sliceRelPosFilter->SetIntermediateSpacingFlag(false);
//sliceRelPosFilter->AddOrientationTypeString(orientation);
- sliceRelPosFilter->AddAngles(angle, 0.0);
+ sliceRelPosFilter->AddAnglesInRad(angle, 0.0);
sliceRelPosFilter->FuzzyMapOnlyFlagOn(); // do not threshold, only compute the fuzzy map
// sliceRelPosFilter->PrintOptions();
sliceRelPosFilter->Update();
typename FloatImageType::Pointer map = sliceRelPosFilter->GetFuzzyMap();
+ // Resize map like object to allow SetBackground
+ map = clitk::ResizeImageLike<FloatImageType>(map, object, GetBackgroundValue());
+
// Remove initial object from the fuzzy map
map = clitk::SetBackground<FloatImageType, ImageType>(map, object, GetForegroundValue(), 0.0, true);
for(int j=0; j<bins; j++) {
sum += ((double)h->GetFrequency(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;
}
}
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;
}
}
git://git.creatis.insa-lyon.fr / clitk.git / blobdiff