Transform a point index into its respective one after registration

Algorithms
1

I did a registration in two steps, Affine -> B-Spline, and now I want to know what is the correspondent index in the new space for a given index in the original space. Here is the code I’m using:

#Receiving the voxel coordinates
print("Please, insert the voxel coordinates (X, Y, Z): ")
voxel = [int(x) for x in input().split()]

physical_voxel = moving_volume.TransformIndexToPhysicalPoint(voxel) 
print(physical_voxel)

#Applying the transformation
transformed_point = BS_transformation.TransformPoint(Affine_transformation.TransformPoint(physical_voxel))

print(transformed_point)
fixed_index = Fixed_volume.TransformPhysicalPointToIndex(transformed_point)
print(fixed_index)
label = fixed_atlas.GetPixel(fixed_index)

Voxel is the index in the original space.

The problem is: I do not get the correct index in the registered space but for me everything is right. So, How to solve this?

2

Hello Everaldo,

Welcome to the ITK community! :sunny:

A few items to note:

physical_voxel = moving_volume.TransformIndexToPhysicalPoint(voxel)

Usually, the registered transform goes from fixed to moving space, so you may need to reverse the volumes used.

transformed_point = BS_transformation.TransformPoint(Affine_transformation.TransformPoint(physical_voxel))

A CompositeTransform could optionally be used to encapsulate these steps.

1 Like
3

Thanks for your reply Matt. What do you mean by “Reverse the volumes”? Do you mean assume that the coordinates VOXEL are in the fixed_volume space and change Fixed_volume by Moving_volume? Or do you mean by creating an inverse transformation?

4

Yes, in the code snippet, relpace moving_volume with Fixed_volume and Fixed_volume with moving_volume.

1 Like
5

Matt, I’ve made the changes you said, but I’m getting some errors. Here is the code I’m using:

#Receiving the voxel coordinates
print("Please, insert the voxel coordinates (X, Y, Z): ")
voxel = [int(x) for x in input().split()]

physical_voxel = fixed_volume.TransformIndexToPhysicalPoint(voxel) 
print(physical_voxel)

#Applying the transformation
transformed_point = BS_result.TransformPoint(Affine_result.TransformPoint(physical_voxel))

print(transformed_point)
fixed_index = moving_volume.TransformPhysicalPointToIndex(transformed_point)
print(fixed_index)
label = fixed_volume.GetPixel(fixed_index)

print(label)

When I test with the entry (192, 160, 58), I get the fixed_index (81, 127, -13) which cause an error. Is my code right or I didn’t understand what you said?

6

Look at this issue in the bug tracker, and its soon-to-be-proposed fix. It should be instructive on what you might be doing wrong in the registration, and how to fix it.

Side note: voxel coordinates entered are IJK, not XYZ as your message indicates.

7

Should be:

moving_index = moving_volume.TransformPhysicalPointToIndex(transformed_point)
print(moving_index)
label = moving_volume.GetPixel(moving_index)

Or, invert the transform if the label from the fixed volume is desired.

1 Like