Hi everyone,
Does anybody know how to use MCFLIRT’s motion parameters to recreate the affine transformation matrix? I’m struggling in particular with the translation parameters.
Setting up the question
Let’s assume we use the following command to motion correct a functional image:
mcflirt -in func.nii.gz -mats -plots -meanvol
This will create a motion parameter file called func_mcf.par and a MAT_???? file (one for each volume) containing the transformation matrix that would motion correct a given volume.
The first line in my .par file contains the following parameters (3 rotation and 3 translation):
0.000388147 -0.000175142 0.000708284 -0.113317 0.265861 0.19274
And the affine transformation matrix for the first volume in MAT_0000 contains the following parameters:
1.000000 0.000708 0.000175 -0.193283
-0.000708 1.000000 0.000388 0.310994
-0.000175 -0.000388 1.000000 0.247711
0.000000 0.000000 0.000000 1.000000
Now, using the first three rotation parameters from the .par file, we can recreate the rotation aspect of MAT_0000. The details are also described here.
def rot_mat(theta):
R_x = np.matrix([[ 1, 0, 0, 0],
[ 0, cos(theta[0]), sin(theta[0]), 0],
[ 0, -sin(theta[0]), cos(theta[0]), 0],
[ 0, 0, 0, 1]
])
R_y = np.matrix([[ cos(theta[1]), 0, -sin(theta[1]), 0],
[ 0, 1, 0, 0],
[ sin(theta[1]), 0, cos(theta[1]), 0],
[ 0, 0, 0, 1]
])
R_z = np.matrix([[ cos(theta[2]), sin(theta[2]), 0, 0],
[-sin(theta[2]), cos(theta[2]), 0, 0],
[ 0, 0, 1, 0],
[ 0, 0, 0, 1]
])
return R_x*R_y*R_z
theta = [0.000388147, -0.000175142, 0.000708284]
R = rot_mat(theta)
>> R
1.000000 0.000708 0.000175 0.000000
-0.000708 1.000000 0.000388 0.000000
-0.000175 -0.000388 1.000000 0.000000
0.000000 0.000000 0.000000 1.000000
Problem
Now, using the translation parameters -0.113317 0.265861 0.19274 from the .par file, how do we get the translation parameters -0.193283 0.310994 0.247711 in the 4th column of the MAT_0000 file?
I have many sources that point to the right solution, but I’m not able to recreate it.
- Source 1 points out that I need to perform a shift to the center of mass of my image
-
Source 2 highlights the same issue and tells me to look into
src/mcflirt/mcflirt.ccandsrc/miscmaths/miscmaths.cc -
mcflirt.ccshows me how to compute the center of mass of an image -
miscmaths.ccdescribes how the transformation matrix is computed. I tried to recreate the code in python, but I wasn’t able to achieve it. I’m also not familiar with C++.
Question
FSL has also a function called avscale that helps to decompose the information inside an affine transformation matrix (stored in the MAT_0000 file).
For example, using the command avscale func_mcf.mat/MAT_0000 | head -n 5 will lead to the output:
Rotation & Translation Matrix:
1.000000 0.000708 0.000175 -0.193283
-0.000708 1.000000 0.000388 0.310994
-0.000175 -0.000388 1.000000 0.247711
0.000000 0.000000 0.000000 1.000000
This is what’s stored in the MAT_0000 file. Now, if you add the non-reference-volume parameter to the avscale command, you can recreate the translation parameters in the .par file.
I.e. using the code avscale func_mcf.mat/MAT_0000 func.nii.gz | head -n 5, will lead to the output:
Rotation & Translation Matrix:
1.000000 0.000708 0.000175 -0.113457
-0.000708 1.000000 0.000388 0.265843
-0.000175 -0.000388 1.000000 0.192804
0.000000 0.000000 0.000000 1.000000
The last column depicts now the same translation values, as stored in the .par file. But how did we get there? The information is somehow hidden here, but I’m not able to extract it.
Or in other words, how can we go from the translation parameters -0.113317 0.265861 0.19274 to the parameters -0.193283 0.310994 0.247711 stored in the MAT_0000 file?
Hints
Everything points to the fact that I need to consider the image center, and transform everything to mm (or back to voxel units?). But I’m unable to do this. Can anybody help me?
Here are also the relevant information from func.nii.gz's header:
$ fslhd func.nii.gz
filename short.nii.gz
dim0 4
dim1 106
dim2 106
dim3 64
dim4 50
vox_units mm
time_units s
datatype 4
pixdim0 0.000000
pixdim1 2.000000
pixdim2 2.000000
pixdim3 2.400000
pixdim4 1.000000
vox_offset 352
scl_slope 1.000000
qform_code 1
qto_xyz:1 -1.997255 -0.004184 0.125602 104.610275
qto_xyz:2 -0.065199 1.613788 -1.415511 -25.650452
qto_xyz:3 0.081988 1.181386 1.934051 -105.191925
qto_xyz:4 0.000000 0.000000 0.000000 1.000000
qform_xorient Right-to-Left
qform_yorient Posterior-to-Anterior
qform_zorient Inferior-to-Superior
sform_code 1
sto_xyz:1 -1.997255 -0.004183 0.125606 104.610275
sto_xyz:2 -0.065200 1.613788 -1.415511 -25.650452
sto_xyz:3 0.081992 1.181386 1.934050 -105.191925
sto_xyz:4 0.000000 0.000000 0.000000 1.000000
sform_xorient Right-to-Left
sform_yorient Posterior-to-Anterior
sform_zorient Inferior-to-Superior
file_type NIFTI-1+
file_code 1
I’m greatful for any help!
-Michael