Quick Summary
It is possible to convert between Slicer registration transform matrices and FSL/FLIRT registration transform matrices. They differ because Slicer transforms describe transformations within RAS space, whereas FLIRT transforms describe transformations between scaled (and sometimes reflected) voxel spaces. Converting between them requires understanding exactly when FLIRT transforms performed a reflection/translation before registration, as well as what other elements of the images are included in the FLIRT transformation matrix which are not part of the Slicer transform matrix. This is not well documented and was not intuitive to me (hence the long explanation below). If you want to skip to something which might “just work” for you, start from the code in the gist linked at the bottom of the post.
Detailed Explanation
Slicer transformation matrices relate RAS spatial coordinates. That is, image data for each image is first situated in the RAS spatial coordinate system, and then the Slicer transformation matrix converts a vector of RAS coordinate (for example the location of some anatomical feature in the moving image) to the corresponding RAS coordinate relative to the fixed image (the location of the same anatomical feature in the fixed image in spatial RAS space)
FSL/FLIRT transformation matrices relate scaled (and sometimes reflected) ijk spatial coordinates. The units are millimeters, but the header information regarding the orientation and location of image volume axes in a world coordinate system (the qform or sform in NIFTI files, or the IJKToRASMatrix in Slicer) is almost entirely ignored. I say “almost entirely” because FSL/FLIRT does look at the IJKToRASMatrix (or equivalent), and if the determinant of that matrix is positive, before registration is carried out, the orientation of the i-axis is reversed and the origin is moved to the other corner voxel on the i-axis. This reflection-translation inverts the sign of the determinant so that it is now negative.
It is possible to obtain the Slicer transformation matrix from the FSL/FLIRT transformation matrix, or vice versa, as long as you also have the two image volumes which were used in the registration which generated the matrix you have, or at the absolute minimum, image volumes which you know to have identical header data (same IJKToRASDirectionMatrix, voxel size, image array dimensions; the voxel contents are not needed) to the image volumes which were registered.
This was not straightforward to figure out, but I think I have it all worked out now, and am going to try to lay it out clearly enough that hopefully this can be helpful to others. In the final understanding, I found the series of coordinate systems in the table below very helpful. Some nomenclature notes first:
I use V for a generic image volume, F for the image volume used as the fixed image in the registration pair (“-ref” in the FLIRT command), and M for the moving image volume in the registration pair (“-in” in the FLIRT command). I use ijk to refer to image voxel array coordinates. I use hjk to refer to the image voxel array coordinates after the i-axis has been optionally flipped and shifted. Note that hjk is sometimes identical to ijk (if the determinant of the IJKToRasDirectionMatrix is negative), that is, I use hjk notation after the flip could have been applied, even if it wasn’t actually applied (because the determinant was already negative). I use hjkMm to refer the coordinate system where hjk indices have been scaled by the voxel size and therefore the coordinate system has millimeter units rather than voxel units. The FLIRT transformation matrix transforms between the moving image’s hjkMm coordinate system to the fixed image’s hjkMm coordinate system. The next coordinate system has the HJKToRASDirectionMatrix rotation applied to hjkMm coordinates; I call this one hras0, where the “h” is a reminder that the voxel origin may have been switched, the “ras” is because we are in a RAS coordinate system, and the “0” is a reminder that the ”h” voxel origin has not yet been moved from [0,0,0]. That is, the transformations from hjk to hras0 involve only scaling and rotation around the origin, so [0,0,0] in hjk space is also [0,0,0] in hras0 space. The last coordinate system is called ras and represents the completion of the transformation from ijk voxel array space to a spatial location in RAS space as defined by the image header. ras differs from hras0 only by a translation transformation which moves the location of [0,0,0] in ijk space, as represented in hras0 space, to the origin location specified in the image header. This ras coordinate system is the same one used as the 3D space RAS coordinate system in 3D Slicer. It no longer has an “h” in the notation because no element remains tying it to the h-axis flip/translation. Between ijk and ras, if a flip/translation occurred, it happened between the ijk coordinate system and the hjk coordinate system, and then was undone by the combination of moving from the hjkMm space to the hras0 space and the translation from hras0 space to ras.
| Coord System Name | Units | Description |
|---|---|---|
| V_ijk | voxels | Voxel array coordinate system for volume V |
| V_hjk | voxels | Voxel array coordinate system for volume V, except that a h coordinate replaces the i coord, and where, Option A: h = i if the determinant of the IJKToRAS matrix for V is negative, or Option B: h = N - 1 - i if the determinant is negative. N = the number of image array voxels in the i dimension. This effectively switches the origin of the image from one corner of the array to the opposite corner (from the original [0,0,0] voxel to the original [N-1,0,0] voxel). This is done by FSL to enforce that all images have a fixed handedness (i.e. forcing a negative determinant). I’m not totally sure why anyone would do this; perhaps because if they didn’t they would have to allow reflections in their transform matrices? |
| V_hjkMm | mm | scaled voxel coordinate system. Origin has not moved relative to the V_hjk coordinate system, but [1,1,1] in the V_hjk coordinate system is located at [Sp_i, Sp_j, Sp_k] in the V_hjkMm coordinate system, where Sp_i is the voxel spacing in the i direction (note Sp_h and Sp_i are the same), and likewise for Sp_j, and Sp_k. More generally, the point [h,j,k] is located at [h*Sp_i, j*Sp_j, k*Sp_k] in the V_hjkMm coordinate system. |
| V_hras0 | mm | rotated coordinate system. Origin has not moved relative to the V_hjk coordinate system, but the axes which point in the h, j, and k directions in the V_hjk and V_hjkMm coordinate systems now point in the directions of the columns of the V_HJKDirs matrix. The V_HJKDirs matrix is the same as the IJKToRAS direction matrix of V if det(V_IJKToRASDirs) < 0, but if det(V_IJKToRASDirs) > 0, then V_HJKDirs = V_IJKToRASDirs @ diag([-1, 1, 1, 1]), because the h axis is reversed with respect to the i axis. Note that V_IJKToRASDirs is IJKToRAS **DIRECTIONS** matrix, not the IJKToRASMatrix (which includes scaling and origin translation). In contrast, the IJKToRASDirection matrix (obtained by volNode.GetIJKToRASDirectionMatrix() in Slicer) has unit length columns and no translation component. |
| V_ras | mm | spatial coordinates in the Slicer world coordinate system. The axis directions match those in V_hras0, but there is a translation-only difference between the two coordinate systems. In V_hras0, the V_hjk origin voxel is still at [0,0,0], but in V_ras, that voxel is located at V_hOri, and the V_ijk origin voxel is located at volNode.GetOrigin(). V_hOri is calculated to compensate for the switch of origin from i to h (if it happened; if det(V_IJKToRASDirs)<0, h=i, and hOri=ori). |
Here is a table of what each of the stepwise transformation matrices are between adjacent coordinate systems:
| Transformation | Calculation |
|---|---|
| V_ijk_to_hjk | Identity matrix if det(IJKToRas) < 0, or translationMatrix([N-1, 0 , 0]) @ diag([-1, 1,1,1]) if det(IJKToRas) > 0, where @ means matrix multiplication (as in numpy) |
| V_hjk_to_hjkMm | diag([Sp_i, Sp_j, Sp_k, 1]) if the voxel spacing is [Sp_i, Sp_j, Sp_k] in the i, j, and k directions, respectively |
| V_hjkMm_to_hras0 | IJKToRASDirectionMatrix if det(IJKToRas) < 0, or IJKToRASDirectionMatrix @ diag([-1, 1, 1, 1]) if det(IJKToRas) > 0 |
| V_hras0_to_ras | translationMatrix(ori - ijk_origin_in_hras0), where ori is the origin listed in the header information for image volume V, and ijk_origin_in_hras0 is the location of the point with ijk=[0,0,0] in the hras0 coordinate system. This can be calculated as V_hjkMm_to_hras0 @ V_hjk_to_hjkMm @ V_ijk_to_hjk @ [0,0,0,1]. |
If a registration has been carried out between a moving image volume M and a fixed image volume F, then one can find the corresponding ijk coordinate in F to an ijk coordinate in M through the following chain of coordinate system transformations:
M_ijk to M_hjk to M_hjkMm to M_hras0 to M_ras to F_ras to F_hras0 to F_hjkMm to F_hjk to F_ijk
Slicer and FLIRT registrations respresent different portions of this chain. A Slicer registration transfrom matrix would represent just the M_ras to F_ras transformation, whereas a FLIRT registration transform matrix represents the M_hjkMm to M_hras0 to M_ras to F_ras to F_hras0 to F_hjkMm transformation.
With this knowledge, we can convert between them!
FLIRT transformation matrix from Slicer transformation matrix:
FLIRT_Matrix = F_ras_to_F_hjkMm @ Slicer_Matrix @ M_hjkMm_to_M_ras
where
M_hjkMm_to_M_ras = M_hras0_to_ras @ M_hjkMm_to_hras0
F_ras_to_F_hjkMm = inv(F_hjkMm_to_hras0) @ inv(F_hras0_to_ras)
Slicer transformation matrix from FLIRT transformation matrix:
Slicer_matrix = F_hjkMm_to_ras @ FLIRT_Matrix @ M_ras_to_hjkMm
where
F_hjkMm_to_ras = F_hras0_to_ras @ F_hjkMm_to_hras0
M_ras_to_hjkMm = inv(M_hjkMm_to_hras0) @ inv(M_hras0_to_ras)