ensuring correct data type of image

I am creating an image
gridImage = sitk.Image([2, 2, slices], sitk.sitkVectorFloat64, 3)

then modifying the contents:

displacements = sitk.GetArrayFromImage(gridImage)
# modify displacements contents
displacementImage = sitk.GetImageFromArray(displacements)
displacementImage.CopyInformation(gridImage)
print(displacementImage)

and the dimensions are still correct, but after creating a displacement field the type info changes to itk::Image<unsigned char, 2u> rather than itk::VectorImage<double, 3u>

dis_tx = sitk.DisplacementFieldTransform(im_vol.GetDimension())
dis_tx.SetInterpolator(sitk.sitkLinear)
dis_tx.SetDisplacementField(get_grid_image)
print(get_grid_image)

im_vol is itk::VectorImage<unsigned char, 3u>

What is causing the type change and How does one manage that?

Platform: OSX
Version: SimpleITK 2.0.0.dev68

Hello @kayarre,

This is a feature not a bug. When you provide the displacement image to the displacement transform it is “consumed” by the transform (the transform owns the memory now and the image has zero size). Please see the doxygen documentation for the SetDisplacementField method. This is a design choice that minimizes the large memory footprint associated with displacement fields (transform/image).

1 Like

I am totally cool with it being a feature. 95% of the time the issue is between chair and keyboard. Thank you for the reply. What is the proper way to do this then?

I think I want to create the inverse transform.
and when I do this following the above.

dis_tx_inv_gen = sitk.InvertDisplacementFieldImageFilter()
dis_tx_inv = dis_tx_inv_gen.Execute(get_grid_image)

it fails because get_grid_image is not the correct type anymore.

Here is what I am trying to do:

  dis_tx = sitk.DisplacementFieldTransform(im_vol.GetDimension())
  dis_tx.SetInterpolator(sitk.sitkLinear)
  dis_tx.SetDisplacementField(get_grid_image)
  inv_image = dis_tx.GetInverseDisplacementField()
  inv_trans = sitk.Transform(inv_image)
  dis_image = dis_tx.GetDisplacementField()
  dis_trans = sitk.Transform(dis_image)

  resampler = sitk.ResampleImageFilter()
  resampler.SetReferenceImage(im_vol_template)
  resampler.SetInterpolator(sitk.sitkLinear)
  resampler.SetTransform(dis_trans)
  im_trans = resampler.Execute(im_vol)

I get a solution using dis_trans, but the scaling appears to be compressing the stack of images rather than expanding them.

when I try the inverse inv_trans I get the following error

Traceback (most recent call last):
  File "build_volume.py", line 392, in <module>
    main()
  File "build_volume.py", line 184, in main
    im_trans = resampler.Execute(im_vol)
  File "/Users/sansomk/anaconda3/envs/reg/lib/python3.7/site-packages/SimpleITK/SimpleITK.py", line 59356, in Execute
    return _SimpleITK.ResampleImageFilter_Execute(self, *args)
RuntimeError: Exception thrown in SimpleITK ResampleImageFilter_Execute: Code/BasicFilters/src/sitkResampleImageFilter.cxx:210:
sitk::ERROR: Unexpected error converting transform! Possible miss matching dimensions!

So Here is the working example:

  dis_tx = sitk.DisplacementFieldTransform(im_vol.GetDimension())
  dis_tx.SetInterpolator(sitk.sitkLinear)
  dis_tx.SetDisplacementField(get_grid_image)
  dis_image = dis_tx.GetDisplacementField()
  invert = sitk.InverseDisplacementFieldImageFilter()
  invert.SetReferenceImage(dis_image)
  inverted = invert.Execute(dis_image)
  inv_trans = sitk.Transform(inverted)

  resampler = sitk.ResampleImageFilter()
  resampler.SetReferenceImage(im_vol_template)
  resampler.SetInterpolator(sitk.sitkLinear)
  resampler.SetTransform(inv_trans)
  im_trans = resampler.Execute(im_vol)

It now works without an error.