concatenating dense deformation fields

I was wondering if it is possible to concatenate two dense deformation fields in ITK python or SimpleITK.
Having

B = T2 ( T1 ( A(x) ) )

In which T1, and T2 are both dense deformation fields (both with size ofImage A and with vector type voxels), is it possible to derive

B = T ( A(x) )

which would be concatenation of T1 and T2.
I know there is a CompositeTransform class but does AddTransform, add the deformation fields or it applies the first one and then the second one.

My reason for concatenation of dense fields is that in my registration framework I predict a dense deformation field at each iteration, but I don’t want to resample the images to avoid interpolation artifacts as a result of multiple resampling.

That might be as simple as adding the deformation fields together using AddImageFilter. If that doesn’t work, report. Also, someone else with more experience with deformation fields might offer more/better adivce.

Problem is adding two deformation fields is not the same as concatenating them. AddImageFilter doesn’t work

You could iterate through all of the points of the final deformation field, and calculate its transform vector via CompositeTransform's TransformPoint method.

The SimpleITK/ITK composite transform does not combine the transform. It applies one transformation then the other.

I believe this can be accomplished with ResampleImageFilter with one transform’s displacement field image being the input image and the other being transform. Careful attention is needed to determine the proper order.

1 Like

Thanks for your comment
I suppose ResampleImageFilter can handle vector type input image and proper interpolation. right?

1 Like

Yes it can!

1 Like

Applying one transformation on the result of the other is the correct way of combining/concatenating two transforms.

ITK handles these composite transforms correctly as long as you only need to compute forward transform and you only try to compute the displacement in a region that is inside of all displacement fields in the transform. Inverse transform computation, especially near the region boundaries may not work reliably (and no inverse can be computed outside the displacement field regions).

Composite transforms are very well supported in 3D Slicer (ITK composite transforms can be read, written, visualized, applied to any kind of data types, resampled to a combined displacement field, etc.), so you can use 3D Slicer for testing, validation, and even batch processing from Python.

3 Likes

Thanks for your comment,
Is there any correct way to avoid multiple resampling of the image as a result of applying multiple transformation? So, instead of T2 ( T1( A(x) ) we go straight to T( A(x) ) since in the former we are resampling twice (and even more, if we have more than 2 dense deformations) and in the latter we are resampling once (no matter how many deformations)

Yes, ITK and Slicer do the resampling correctly - resampling the image only once, with the combined transform.

2 Likes

Here is a SimpleITK example which creates a composite transform of an affine transformed followed by a displacement field: https://simpleitk.readthedocs.io/en/master/Examples/ImageRegistrationMethodDisplacement1/Documentation.html

At the end resampling is done once, with the combination of the transform. So the pretense of needing to manually combine transforms into one to avoid multiple resampling is not correct.

But if you need to combine two displacement field transforms into one: T_2( T_1(\vec{p}))=T(\vec{p}), then resampling T_2's displacement field image with the T_1 transform will yield the same results on the defined region as using a composite transform. The same linear interpolation and sampling of T_2 with will occur when the fixed image’s domain ( the set of \vec{p} ) , and the displacement field images of T_1 and T_2 are the same. This resampling can be done to update T_2, so as the T_1 displacement field approaches zero, T_2 will change less and less.