[Insight-developers] LevelsetMotionRegistrationImageFilter

[Kishore Mosaliganti]

Hi all,

I was using the LevelsetMotionRegistrationImageFilter to register two
3D image datasets whose voxels have anisotropic spacing:

  RegistrationFilterType::Pointer filter = RegistrationFilterType::New();
  CommandIterationUpdate::Pointer observer = CommandIterationUpdate::New();
  filter->AddObserver( itk::IterationEvent(), observer );
  filter->SetFixedImage( fixedImage );
  filter->SetMovingImage( movingImage );
  filter->SetNumberOfIterations( 50 );
  filter->SetGradientSmoothingStandardDeviations( 0.1 );
  filter->UseImageSpacingOn();
  filter->Update();

The problem is I get a Segmentation fault when I turn on
UseImageSpacingOn(). Otherwise, the filter executes normally.

I used the debugger and went into the code.
LevelsetMotionRegistrationImageFilter is a sub-class of
PDEDeformableRegistrationImageFilter which is a subclass of
FiniteDifferenceImageFilter. The code failed in the
FiniteDifferenceImageFilter::GenerateData() method in the following
lines:

....

if (this->GetState() == UNINITIALIZED)
    {
    // Set the coefficients for the deriviatives
    double coeffs[TInputImage::ImageDimension];
    if (m_UseImageSpacing)
      {
      for (unsigned int i = 0; i < TInputImage::ImageDimension; i++)
        {
        coeffs[i] = 1.0 / this->GetInput()->GetSpacing()[i];
        }
      }
    else
      {
      for (unsigned int i = 0; i < TInputImage::ImageDimension; i++)
        {
        coeffs[i] = 1.0;
        }
      }

.....

When I went to see the state of the input pointer, this->m_Input, it
is a vector of 3 image pointers. Basically, this->m_Input[1] refers to
the fixed image and this->m_Input[2] refers to the moving image.
this->m_Input[0] is NULL. this->m_Input[0] is supposed to contain the
initial deformation field which is not a prerequisite for running the
code.

The code in itkPDEDeformableRegistrationImageFilter::SetFixedImage()
and itkPDEDeformableRegistrationImageFilter::SetMovingImage() call
this->ProcessObject::SetNthInput( 1, const_cast< FixedImageType * >(
ptr ) ) and this->ProcessObject::SetNthInput( 2, const_cast<
FixedImageType * >( ptr ) ) respectively.

This filter derives from
ImageToImageFilter<TDeformationField,TDeformationField> but the vector
of inputs have been set as scalar-valued fixed and moving image
pointers. Also, for the case when the initial deformation field is
unknown, the FiniteDifferenceFilter is checking for the first input
which doesn't exist.

I have filed a bug report:
http://public.kitware.com/Bug/view.php?id=9709



Thanks,
Kishore