[Insight-developers] image information reset by gradient magnitude filter

[Richard Beare]

Hi,
I've come across a strange resetting of image information when using
the gradient magnitude filter (code below). All spacing and
orientation information seems to have been reset to default values.
This is using itk 3.10.1 retrieved using the cvs repository, and the
cvs head from about 2 weeks ago under linux. I haven't yet checked any
other filters to see if the problem is more general.

nifti file report from fslhd of input image:

sizeof_hdr     348
data_type      INT16
dim0           3
dim1           256
dim2           120
dim3           256
dim4           1
dim5           1
dim6           1
dim7           1
vox_units      mm
time_units     Unknown
datatype       4
nbyper         2
bitpix         16
pixdim0        0.0000000000
pixdim1        0.8203120232
pixdim2        1.1999969482
pixdim3        0.8203120232
pixdim4        0.0000000000
pixdim5        0.0000000000
pixdim6        0.0000000000
pixdim7        0.0000000000
vox_offset     352
cal_max        0.0000
cal_min        0.0000
scl_slope      1.000000
scl_inter      0.000000
phase_dim      0
freq_dim       0
slice_dim      0
slice_name     Unknown
slice_code     0
slice_start    0
slice_end      0
slice_duration 0.000000
time_offset    0.000000
intent         Unknown
intent_code    0
intent_name
intent_p1      0.000000
intent_p2      0.000000
intent_p3      0.000000
qform_name     Aligned Anat
qform_code     2
qto_xyz:1      -0.820312  -0.000000  -0.000000  105.000000
qto_xyz:2      0.000000  1.199997  -0.000000  -74.299637
qto_xyz:3      -0.000000  0.000000  0.820312  -114.179565
qto_xyz:4      0.000000  0.000000  0.000000  1.000000
qform_xorient  Right-to-Left
qform_yorient  Posterior-to-Anterior
qform_zorient  Inferior-to-Superior
sform_name     Scanner Anat
sform_code     1
sto_xyz:1      -0.820312  -0.000000  0.000000  105.000000
sto_xyz:2      -0.000000  1.199997  -0.000000  -74.299637
sto_xyz:3      0.000000  0.000000  0.820312  -114.179565
sto_xyz:4      0.000000  0.000000  0.000000  1.000000
sform_xorient  Right-to-Left
sform_yorient  Posterior-to-Anterior
sform_zorient  Inferior-to-Superior
file_type      NIFTI-1+
file_code      1
descrip
aux_file

Result image:
sizeof_hdr     348
data_type      INT16
dim0           3
dim1           256
dim2           120
dim3           256
dim4           1
dim5           1
dim6           1
dim7           1
vox_units      mm
time_units     Unknown
datatype       4
nbyper         2
bitpix         16
pixdim0        0.0000000000
pixdim1        1.0000000000
pixdim2        1.0000000000
pixdim3        1.0000000000
pixdim4        0.0000000000
pixdim5        0.0000000000
pixdim6        0.0000000000
pixdim7        0.0000000000
vox_offset     352
cal_max        0.0000
cal_min        0.0000
scl_slope      1.000000
scl_inter      0.000000
phase_dim      0
freq_dim       0
slice_dim      0
slice_name     Unknown
slice_code     0
slice_start    0
slice_end      0
slice_duration 0.000000
time_offset    0.000000
intent         Unknown
intent_code    0
intent_name
intent_p1      0.000000
intent_p2      0.000000
intent_p3      0.000000
qform_name     Aligned Anat
qform_code     2
qto_xyz:1      -1.000000  0.000000  0.000000  -0.000000
qto_xyz:2      0.000000  -1.000000  0.000000  -0.000000
qto_xyz:3      0.000000  0.000000  1.000000  0.000000
qto_xyz:4      0.000000  0.000000  0.000000  1.000000
qform_xorient  Right-to-Left
qform_yorient  Anterior-to-Posterior
qform_zorient  Inferior-to-Superior
sform_name     Scanner Anat
sform_code     1
sto_xyz:1      -1.000000  0.000000  0.000000  -0.000000
sto_xyz:2      0.000000  -1.000000  0.000000  -0.000000
sto_xyz:3      0.000000  0.000000  1.000000  0.000000
sto_xyz:4      0.000000  0.000000  0.000000  1.000000
sform_xorient  Right-to-Left
sform_yorient  Anterior-to-Posterior
sform_zorient  Inferior-to-Superior
file_type      NIFTI-1+
file_code      1
descrip
aux_file

Code:

#include <itkGradientMagnitudeRecursiveGaussianImageFilter.h>
#include <itkImageFileWriter.h>
#include <itkImageFileReader.h>


template <class TImage>
void writeIm(typename TImage::Pointer Im, std::string filename)
{
  typedef typename itk::ImageFileWriter<TImage> WriterType;
  typename WriterType::Pointer writer = WriterType::New();
  writer->SetInput(Im);
  writer->SetFileName(filename.c_str());
  writer->Update();
}

template <class TImage>
typename TImage::Pointer readIm(std::string filename)
{
  typedef typename itk::ImageFileReader<TImage> ReaderType;
  typename ReaderType::Pointer reader = ReaderType::New();
  reader->SetFileName(filename.c_str());
  typename TImage::Pointer result = reader->GetOutput();
  try
    {
    result->Update();
    }
  catch(itk::ExceptionObject &ex)
    {
    std::cout << ex << std::endl;
    std::cout << filename << std::endl;
    return 0;
    }
    result->DisconnectPipeline();
    return(result);
}

int main(int argc, char * argv[])
{

  const int dim = 3;
  typedef short PixType;
  typedef itk::Image<PixType, dim> RawImType;

  typedef itk::GradientMagnitudeRecursiveGaussianImageFilter<RawImType,RawImType>
GradFiltType;

  RawImType::Pointer raw;
  RawImType::Pointer grad;

  raw = readIm<RawImType>(argv[1]);
  GradFiltType::Pointer gradfilt = GradFiltType::New();

  gradfilt->SetInput(raw);
  gradfilt->SetSigma(0.01);
  grad = gradfilt->GetOutput();
  grad->Update();
  grad->DisconnectPipeline();
  std::cout << raw->GetSpacing() << grad->GetSpacing() << std::endl;
  writeIm<RawImType>(grad, "grad.nii.gz");

  return (EXIT_SUCCESS);
}