[Insight-users] itkBilateralImageFilter only works for Origin 0,0,0

[jean-michel.rouet at philips.com]

I have the feeling that the itkBilateralImageFilter assumes that the input 
origin is 0,0,0.

As a demo, I tried the following code (inspired from 
itkBilateralImageFilterTest2.cxx):

#include "itkBilateralImageFilter.h"
#include "itkMinimumMaximumImageCalculator.h"
#include "itkImageFileReader.h"
#include "itkImageFileWriter.h"

int main(int argc, char* argv[] )
{
  typedef unsigned char PixelType;
  const unsigned int Dimension = 3;

  typedef itk::Image<PixelType,Dimension> myImageType;
  typedef itk::BilateralImageFilter<myImageType,myImageType> myFilterType;

  // read input
  itk::ImageFileReader<myImageType>::Pointer
      input = itk::ImageFileReader<myImageType>::New();
      input->SetFileName(argv[1]);
      input->Update();
 
  // Create the filter
  myFilterType::Pointer filt = myFilterType::New();
      filt->SetInput(input->GetOutput());
      filt->SetDomainSigma( 4.0 );
      filt->SetRangeSigma( 50.0 );
      filt->SetFilterDimensionality( Dimension );
      filt->Update();
 
  itk::MinimumMaximumImageCalculator<myImageType>::Pointer
      minmax = itk::MinimumMaximumImageCalculator<myImageType>::New();
  minmax->SetImage( filt->GetOutput() );
  minmax->ComputeMinimum();
  minmax->ComputeMaximum();

  int min = minmax->GetMinimum();
  int max = minmax->GetMaximum();
 
  std::cout << "minimum = " << min << std::endl;
  std::cout << "maximum = " << max << std::endl;

  // Generate test image
  itk::ImageFileWriter<myImageType>::Pointer
      writer = itk::ImageFileWriter<myImageType>::New();
      writer->SetInput( filt->GetOutput() );
      writer->SetFileName( "filtered.mhd" );
      writer->Update();
 
  return 0;
}

and as input I took the HeadrMRVolume.mhd from the Testing/Data/Input 
directory.
It works correctly (at least it outputs something!!! :-) )

Then I tried to modify the Origin of the input volume either by using 
SetOrigin() on the input->GetOutput(), or by adding an Offset line in the 
mhd header as 
follows (origin set to -100 -200 -300) :

rouet at joebar /tmp>cat HeadMRVolume.mhd
NDims = 3
DimSize = 48 62 42
Offset = -100 -200 -300
ElementSize = 4.000000e+000 4.000000e+000 4.000000e+000
ElementSpacing = 4.000000e+000 4.000000e+000 4.000000e+000
ElementType = MET_UCHAR
ElementByteOrderMSB = False
ElementDataFile = HeadMRVolume.raw

then I ran again the test program, and it generated an empty image (filled 
with 0), (maximum = minimum = 0)

I think this is not an expected result since the filter has no reason to 
be origin dependant !
Am I right ?

Regards,
Jean-Michel

PS: 

I guess that the error comes from line 147 in itkBilateralImageFilter.txx 
where the center of the gaussian should take into account the origin 
    mean[i] = this->GetInput()->GetSpacing()[i]*radius[i]; // center pixel 
pos
should become
    mean[i] = this->GetInput()->GetSpacing()[i]*radius[i] + 
this->GetInput()->GetOrigin()[i]; // center pixel pos

Other solution could be to remove the line 141:
          gaussianImage->SetOrigin( this->GetInput()->GetOrigin() );
so that the gaussianImage kernel is independent of the origin.