2

I am using image registration toolbox to register two 3D images. I register the moving image to the fixed image. I use "imregtform" to save tform. 

tform = imregtform(moving,fixed,transformType,optimizer,metric)

for example this is my 'tform':

    1.0235    0.0022   -0.0607         0
   -0.0276    1.0002    0.0089         0
   -0.0170   -0.0141    1.1685         0
   12.8777    5.0311  -70.0325    1.0000

then I use 'Imwarp' to transfer moving image to the fixed image. In my code it is,

nii=load_untouch_nii(['mypath image.nii' ]);
I = nii.img; 
dii=nii.hdr.dime.pixdim(2:4);
Rfixed=imref3d(size(I),dii(2),dii(1),dii(3)); 
timg= imwarp(I, Rfixed, tform, 'OutputView', Rfixed);

'timg' is transfered image of source image.I checked it and it is working fine, but I need control points or displacement field of this transformation now. In another word, I need to know each voxel(3d pixel) moved to which position. If I know this I can draw vector field. In deformable image registartion methods such as NiftyReg package, control point command is provided to make it easy for users but I do not know how to do it in simple 3d affine in Matlab. Any help will be appreciated too much

Ehsan
  • 517
  • 1
  • 7
  • 32

1 Answers1

3

You can apply your geometric transform to a single point very easily.

You just need to have your point defined as p=[x;y;z;1]

and then obtain pt by pt=p*tform; pt=pt(1:3)./pt(4); (search for homogeneous coordinates for explanation of this last division). This is what imwarp does in the inside. It performs the said multiplication for each pixel p.

To then obtain the displacement, you would need just to disp=pt-p;. Note that the displacements are most likely not integer.

Note that in the field of deformable image registration, generally you will have a different tform for each control point (that is why it is called deformable and not rigid)

EDIT: as @Ashish Uthama suggests in the comments, you can also do it with Matlab inbuilt function transformpointsforward(). I will always be a promoter of "writing your own code" philosophy, specially the first time, so you do understand what are you actually doing.

http://uk.mathworks.com/help/images/ref/affine3d.transformpointsforward.html

Ander Biguri
  • 35,140
  • 11
  • 74
  • 120
  • Thanks a lot for quick reply, it seems working. I want to apply it to simple grid. `[U, V, W] = ndgrid(336, 166,336); p = [100;70;110;1]; pt=pt(1:3)./pt(4); ` the result is ` -0.0158 -0.0113 -0.0208`, so this is displacement or control point? – Ehsan May 14 '15 at 15:27
  • I called it `pt` as "point transformed". The displacement would be `disp=pt-p;` @Ehsan – Ander Biguri May 14 '15 at 15:47
  • 1
    Also: http://www.mathworks.com/help/images/ref/affine3d.transformpointsforward.html :) – Ashish Uthama May 14 '15 at 16:10
  • Thanks for both of you, i know what is going behind and what is matlab function now. Regards – Ehsan May 14 '15 at 17:05
  • @AnderBiguri I have implemented based on what you said, but the result of my code is not same as imwarp. Can you look at [my code](http://stackoverflow.com/questions/30437467/how-imwarp-transfer-points-in-matlab) ? Thanks. – Ehsan May 25 '15 at 11:38
  • @AshishUthama I have problem with `transformpointsforward` because in my case my result is not match with imwarp result. any help ? Thanks – Ehsan May 25 '15 at 11:40