script to compress MHDs

[clitk.git] / itk / clitkAddRelativePositionConstraintToLabelImageFilter.txx

diff --git a/itk/clitkAddRelativePositionConstraintToLabelImageFilter.txx b/itk/clitkAddRelativePositionConstraintToLabelImageFilter.txx
index e25f98b829fa3c65ed497f5fb5ffd73de863a066..6d536186fbf5084da969d0750e05a4cc8cd3c8fb 100644 (file)
--- a/itk/clitkAddRelativePositionConstraintToLabelImageFilter.txx
+++ b/itk/clitkAddRelativePositionConstraintToLabelImageFilter.txx
@@ -63,6 +63,8 @@ AddRelativePositionConstraintToLabelImageFilter():
   VerboseStepFlagOff();
   WriteStepFlagOff();
   FuzzyMapOnlyFlagOff();
+  FastFlagOff();
+  SetRadius(2.0);
 }
 //--------------------------------------------------------------------
 
@@ -380,12 +382,16 @@ GenerateData()
   for(int i=0; i<GetNumberOfAngles(); i++) {
     // Compute fuzzy map
     relPosFilter = RelPosFilterType::New();
+    relPosFilter->SetFast(GetFastFlag());
+    relPosFilter->SetRadius(GetRadius());
     relPosFilter->SetInput(working_image);
     relPosFilter->SetAlpha1(m_Angle1[i]); // xy plane
     relPosFilter->SetAlpha2(m_Angle2[i]);
     relPosFilter->SetK1(M_PI/2.0); // Opening parameter, default = pi/2
-    relPosFilter->SetFast(true);
-    relPosFilter->SetRadius(1); // seems sufficient in this case
+
+    // relPosFilter->SetFast(true);
+    // relPosFilter->SetRadius(1); // seems sufficient in this case
+
     // relPosFilter->SetVerboseProgress(true);
     relPosFilter->Update();
     relPos = relPosFilter->GetOutput();
@@ -426,7 +432,43 @@ GenerateData()
 
   relPos = m_FuzzyMap;
   StopCurrentStep<FloatImageType>(relPos);
-  if (GetFuzzyMapOnlyFlag()) return;
+  if (GetFuzzyMapOnlyFlag()) {
+    // If the spacing is used, recompute fuzzy map with the input size/spacing
+    if (m_IntermediateSpacingFlag) {
+      StartNewStep("Resample FuzzyMap to come back to the same sampling than input");
+      typedef clitk::ResampleImageWithOptionsFilter<FloatImageType> ResampleFilterType;
+      typename ResampleFilterType::Pointer resampleFilter = ResampleFilterType::New();
+      resampleFilter->SetDefaultPixelValue(0.0); //Default fuzzymap value FIXME
+      resampleFilter->SetInput(relPos);
+      resampleFilter->SetOutputSpacing(input->GetSpacing());
+      resampleFilter->SetGaussianFilteringEnabled(false);
+      resampleFilter->Update();
+      relPos = m_FuzzyMap = resampleFilter->GetOutput();
+
+      // Need to put exactly the same size
+      if (relPos->GetLargestPossibleRegion() != input->GetLargestPossibleRegion()) {
+        StartNewStep("Pad to get the same size than input");
+        typename FloatImageType::Pointer temp = FloatImageType::New();
+        temp->CopyInformation(input);
+        temp->SetRegions(input->GetLargestPossibleRegion()); // Do not forget !!
+        temp->Allocate();
+        temp->FillBuffer(0.0); // Default fuzzymap value FIXME
+        typename PasteFloatFilterType::Pointer padFilter2 = PasteFloatFilterType::New();
+        padFilter2->SetSourceImage(relPos);
+        padFilter2->SetDestinationImage(temp);
+        padFilter2->SetDestinationIndex(input->GetLargestPossibleRegion().GetIndex());
+        padFilter2->SetSourceRegion(relPos->GetLargestPossibleRegion());
+        padFilter2->Update();
+        relPos = padFilter2->GetOutput();
+        StopCurrentStep<FloatImageType>(relPos);
+        m_FuzzyMap = relPos;
+      }     
+    }
+    else {
+      StopCurrentStep<FloatImageType>(relPos);
+    }
+    return;
+  }
                
   //--------------------------------------------------------------------
   //--------------------------------------------------------------------