Surfaces#

A Surface is a triangle mesh in world coordinates (RAS+, in mm), such as the white matter surface FreeSurfer makes, with named values per vertex, such as its curvature or thickness, in vertex_data. Larmor reads and writes FreeSurfer surfaces (lh.white, lh.pial, …) and GIfTI surfaces (.surf.gii).

A surface and its values#

The faces list the vertex indices of each triangle, counter-clockwise seen from outside, so that the normals point out. An octahedron of radius 20 mm, with a value per vertex, saved as a FreeSurfer surface with a FreeSurfer values file beside it, and as GIfTI:

import numpy as np

from larmor.io import save_surface
from larmor.io.formats import write_freesurfer_values
from larmor.types import Surface

vertices = 20 * np.array(
    [[1, 0, 0], [-1, 0, 0], [0, 1, 0], [0, -1, 0], [0, 0, 1], [0, 0, -1]], float
)
top = [[0, 2, 4], [2, 1, 4], [1, 3, 4], [3, 0, 4]]
bottom = [[2, 0, 5], [1, 2, 5], [3, 1, 5], [0, 3, 5]]
faces = np.array(top + bottom)
thickness = np.array([2.0, 2.5, 3.0, 2.2, 2.8, 2.4])
surface = Surface(vertices, faces, vertex_data={"thickness": thickness})

save_surface(surface, "lh.white")
write_freesurfer_values(thickness, "lh.thickness", num_faces=surface.num_faces)
save_surface(surface, "lh.white.surf.gii")

A FreeSurfer surface file holds the mesh only, so its values go into files of their own, such as lh.thickness; a GIfTI surface keeps them. load_surface() reads either, with values files given by values, each named after the part of its name after the hemisphere:

from larmor.io import load_surface

white = load_surface("lh.white", values="lh.thickness")
print(white.num_vertices, white.num_faces, sorted(white.vertex_data))

Normals, cutting and joining#

from larmor.surfaces import concatenate_surfaces, keep_vertices, vertex_normals

normals = vertex_normals(white)
upper = keep_vertices(white, white.vertices[:, 2] >= 0)
both = concatenate_surfaces([white, load_surface("lh.white.surf.gii")])
print(upper.num_vertices, both.num_vertices)

vertex_normals() gives a unit normal per vertex, pointing out. keep_vertices() keeps the vertices where a mask is True, with the faces between them, for example the cortex of a FreeSurfer label read with read_freesurfer_label(), which removes the medial wall. concatenate_surfaces() joins surfaces into one mesh, such as the two hemispheres.

Moving a surface into another space#

A registration computed by ANTs, FSL or MRtrix3, such as one from a T1 image to the diffusion image, brings a surface into the diffusion space. Here, an MRtrix3 linear transform, a 4x4 matrix in world coordinates, written for the example:

np.savetxt("t1_to_dwi.txt", [[1, 0, 0, 5], [0, 1, 0, 0], [0, 0, 1, 0], [0, 0, 0, 1]])

from larmor.io import load_surface, save_surface
from larmor.registration import apply_transform, load_registration

transform = load_registration(["t1_to_dwi.txt"])
moved = apply_transform(load_surface("lh.white"), transform)
save_surface(moved, "lh.white.dwi.surf.gii")

lar transform lh.white lh.white.dwi.surf.gii -t t1_to_dwi.txt

The transform files are given as for an image; the vertices move so that the surface lands where the registered image does. Bringing an atlas region into diffusion space shows a registration on images.

Viewing surfaces#

lar bv --surfaces lh.white

The viewer draws a surface as a shaded mesh in 3D and as the line where each slice cuts it in 2D. A FreeSurfer surface opens with its lh.curv and ../label/lh.cortex.label when they exist; values per vertex can colour it.

See also#