reading surface ras2vox transformation using nibabel

[EhsanTadayon]

Hi guys,

I'm trying to get the ras2vox transformation matrix for a surface file generated by freesurfer, which can be extracted using mris_info. Is there any way to get surface ras2vox transformation matrix using nibabel?

thanks
Ehsan.

[effigies]

The header affine takes coordinates from voxel space to RAS, so the inverse should be RAS2VOX.

import numpy as np
import nibabel as nb

np.linalg.inv(nb.load(fname).affine)

FreeSurfer does some weird things with a center offset that may or may not come into play here. Let me know if the above doesn't work for you.

[EhsanTadayon]

Hi Chris,

Thanks. It didn't work. I believe that works as long as you input brain volume.
Let me go a little bit more into detail what I have done. I have created brain skin surfaces using freesurfer mri_watershed. I am trying to overlay bunch of electrode positions ( RAS) onto the created skin surface, which has the vertices in tkrRAS. For that, I need to find the transformation matrices ( vox2tkrRAS and vox2ras) stored in the generated surface file. I can do that for each single subject by calling mris_info to get these matrices, but I was wondering if something has been already implemented in nibabel to retrieve those matrices from the surface files. I couldn't find any on the nibabel's freesurfer documentation.

Thanks.

[effigies]

Ah, I misunderstood. So surface files should have associated volume geometry structures at the end, which nibabel calls "metadata". The can be accessed like so:

import numpy as np
import nibabel as nb

coords, faces, volume_info = nb.freesurfer.read_geometry(fname, read_metadata=True)

From this, you can build your VOX2RAS affine (I think I've gotten the shapes right, but haven't tested. You can fiddle with this until you get it right):

A = np.vstack((volume_info['xras'],
               volume_info['yras'],
               volume_info['zras'])).T * volume_info['voxelsize']
b = volume_info['cras'].reshape(3, 1) - A.dot(volume_info['volume'].reshape(3, 1)) / 2

affine = np.eye(4)
affine[:3, :3] = A
affine[:3, [3]] = b

And then you can go ahead and invert:

ras2vox = np.linalg.inv(affine)

I haven't tested this code, but the principle is that the x_ras .. z_ras are the direction cosine vectors, (the columns in a rotation matrix):

x_r y_r z_r
x_a y_a z_a
x_s y_s z_s

These need to be scaled by the size of the voxels in each direction (the voxelsize vector), resulting in A, above.

c_ras is the offset of the origin from the center of the volume. The center of the volume in RAS can be found by taking the dimensions of the volume, dividing by half, and applying the rotation and scaling matrix A. However, the translation vector b is relative to (0, 0, 0), so we need to subtract the centroid to make that work.

Hopefully that explanation is enough to help you work past any coding errors I made. (I'm drawing off of #565, which is work I started on over 6 months ago.)

[EhsanTadayon]

Hi Chris,

Thanks for this. It was very informative. It worked for me.