[Insight-users] Explicite representation of the registration matrix

[Maxime Toussaint]

Hi to you unknows reader,

I believe that it is my first time using that kind of help. Most of the
time, I only use answer that are already in the web. I will try my best to
make my question/situation as clear as possible.

Also, I am sorry for my mild english. I do know that it is not my strenght,
but I do try my best.


The question in one line :
Is there a class/function that produce the deformation matrix of an known
BSpline deformation?


The context :
I use ITK to register an image G1 to an image G0. The registration
parameter are the weight of an uniformly spaced BSpline. (Thus, my
deformation is non-rigid.) The interpolation function is the basic
one(linea interpolation).

What I want is the matrix W that wrap G0 to G1. In maths : G1 = W * G0
(Where G0 and G1 are column vector)

Why would I want to obtain its explicite form and not just use ITK to wrap
G0? Well, I study in optimisation and I need to be able to evaluate the
following : v * W (Where v is a line vector)

Which cannot, from my knowledge, be done directly by ITK.

The following are work around that I tried.


First attempt:
A first work around was to create an image with one voxel setted has 1 and
the other has 0. Wrap the image with ITK and extract the deformed image.
Lets say that the voxel with a 1 is the i voxel. Thus, the deformed image
give me the i column of W.
W * e_i = ( W_1 | W_2 | ... | W_i | ... | W_N) * e_i = W_i         (Where
e_i is a vector column of 0 except at i where it is equal to 1)

But this solution is sub sub sub optimal. Having to wrap the image N time
is unthinkable in computing time.


Second attempt :
Use the function TransformPoint( ) of the transformation class to obtain
the new position of the voxel and interpolate the W matrix from that.

The problem is I can't understand how the linear interpolation work.
Example :
Position of the non null value : [0, 0, 0]   (e_0)
Deformed position : [0.0231481, 0, 0]
Value of W * e_0 at [0, 0, 0] is 0.976852   (the rest is 0)
Which blow my mind. The voxel [0, 0, 0] got nearer from the voxel [1, 0, 0]
since he is now at [0.0231481, 0, 0], but its contribution to [1, 0, 0] is
null.


Well, I hope that my explanation was complete. I didn't post my code since
I believe that it is, out of context, incomprehensible.

Thank you in advance for your help,
One of the many users of ITK
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.itk.org/pipermail/insight-users/attachments/20130413/883c6f8e/attachment.htm>