Constructing 10-10 coordinates on segmented MRI scans

[1]:
import cedalion
import cedalion.io
import cedalion.geometry.segmentation
import cedalion.geometry.landmarks
from cedalion.datasets import get_colin27_segmentation
import os.path
import pyvista

#pyvista.set_jupyter_backend("html")
pyvista.set_jupyter_backend("static")

Load segmentation masks

This example constructs the 10-10 system on the Colin27 average brain.

[2]:
SEG_DATADIR, mask_files, landmarks_files = get_colin27_segmentation()
masks, t_ijk2ras = cedalion.io.read_segmentation_masks(SEG_DATADIR, mask_files)
Downloading file 'colin27_segmentation.zip' from 'https://doc.ml.tu-berlin.de/cedalion/datasets/colin27_segmentation.zip' to '/home/runner/.cache/cedalion'.
Unzipping contents of '/home/runner/.cache/cedalion/colin27_segmentation.zip' to '/home/runner/.cache/cedalion/colin27_segmentation.zip.unzip'

Load hand-picked landmarks

[3]:
initial_landmarks = cedalion.io.read_mrk_json(os.path.join(SEG_DATADIR, landmarks_files), crs="aligned")
initial_landmarks
[3]:
<xarray.DataArray (label: 4, aligned: 3)> Size: 96B
<Quantity([[-4.99239750e-02  7.97894669e+01 -3.61411209e+01]
 [ 1.75606728e+00 -1.00949997e+02 -5.48579865e+01]
 [-7.19499969e+01 -1.61194057e+01 -5.46569786e+01]
 [ 7.59499969e+01 -1.30380936e+01 -5.40778427e+01]], 'millimeter')>
Coordinates:
  * label    (label) <U3 48B 'Nz' 'Iz' 'LPA' 'RPA'
    type     (label) object 32B PointType.LANDMARK ... PointType.LANDMARK
Dimensions without coordinates: aligned

Construct scalp surface and transform to scanner RAS

[4]:
scalp_surface = cedalion.geometry.segmentation.surface_from_segmentation(
    masks,
    masks.segmentation_type.values,
    fill_holes_in_mask=True
)
scalp_surface = scalp_surface.apply_transform(t_ijk2ras)
scalp_surface
[4]:
TrimeshSurface(mesh=<trimesh.Trimesh(vertices.shape=(152529, 3), faces.shape=(304426, 3))>, crs='aligned', units=<Unit('millimeter')>)

Construct landmarks

[5]:
lmbuilder = cedalion.geometry.landmarks.LandmarksBuilder1010(scalp_surface, initial_landmarks)
all_landmarks = lmbuilder.build()
/home/runner/work/cedalion/cedalion/src/cedalion/geometry/landmarks.py:216: UserWarning: WIP: distance calculation around ears
  warnings.warn("WIP: distance calculation around ears")

Visualize

[6]:
lmbuilder.plot()
display(all_landmarks)
../_images/examples_1010_system_11_0.png
<xarray.DataArray (label: 73, aligned: 3)> Size: 2kB
<Quantity([[-4.99239750e-02  7.97894669e+01 -3.61411209e+01]
 [ 1.75606728e+00 -1.00949997e+02 -5.48579865e+01]
 [-7.19499969e+01 -1.61194057e+01 -5.46569786e+01]
 [ 7.59499969e+01 -1.30380936e+01 -5.40778427e+01]
 [ 1.48924732e+00 -2.80293236e+01  8.91000061e+01]
 [ 6.12858720e-02  8.11088104e+01  7.37103820e-01]
 [ 2.87298560e-01  7.13674164e+01  3.77344780e+01]
 [ 6.55276775e-01  4.53056488e+01  6.68061066e+01]
 [ 1.07111287e+00  1.00979099e+01  8.20730667e+01]
 [ 1.85324764e+00 -6.58003693e+01  8.12256699e+01]
 [ 2.05901766e+00 -9.44351578e+01  5.45259399e+01]
 [ 2.09403157e+00 -1.09099968e+02  2.00678768e+01]
 [ 1.99313700e+00 -1.12099998e+02 -1.88979530e+01]
 [-7.40999985e+01 -1.95488510e+01 -1.85347042e+01]
 [ 7.71044159e+01 -1.64808197e+01 -1.70659504e+01]
 [-2.78219814e+01  7.71388321e+01 -3.02536416e+00]
 [-4.96633530e+01  6.04367256e+01 -7.37024927e+00]
 [-6.30418015e+01  3.72775650e+01 -1.13610516e+01]
 [-6.90998383e+01  1.09453211e+01 -1.47975235e+01]
 [-7.24252396e+01 -4.38436012e+01 -2.08314400e+01]
...
 [ 3.38006325e+01  6.86656570e+01  2.46341743e+01]
 [ 4.37709770e+01  6.53327713e+01  9.80877018e+00]
 [-7.11000137e+01 -5.21063499e+01  1.62688789e+01]
 [-6.04117126e+01 -5.93472214e+01  4.92480507e+01]
 [-3.35214653e+01 -6.43834763e+01  7.31401901e+01]
 [ 3.70119171e+01 -6.36731415e+01  7.34150696e+01]
 [ 6.40872269e+01 -5.81273575e+01  4.98854256e+01]
 [ 7.49943771e+01 -5.05679817e+01  1.65376186e+01]
 [-6.02266083e+01 -7.88670959e+01  6.22620153e+00]
 [-4.79947624e+01 -8.71012192e+01  3.13694839e+01]
 [-2.63569469e+01 -9.27776489e+01  4.89596596e+01]
 [ 3.03763351e+01 -9.22147827e+01  4.83610687e+01]
 [ 5.25887566e+01 -8.63931656e+01  3.11849651e+01]
 [ 6.40652771e+01 -7.80682831e+01  6.22846365e+00]
 [-4.19541206e+01 -9.81257248e+01 -7.08912945e+00]
 [-3.21703415e+01 -1.03932159e+02  6.96908283e+00]
 [-1.80093327e+01 -1.08100220e+02  1.72244530e+01]
 [ 2.07628002e+01 -1.08123955e+02  1.81712074e+01]
 [ 3.61665001e+01 -1.03939384e+02  8.55530739e+00]
 [ 4.61439362e+01 -9.79618149e+01 -5.45690823e+00]], 'millimeter')>
Coordinates:
  * label    (label) <U3 876B 'Nz' 'Iz' 'LPA' 'RPA' ... 'PO1' 'PO2' 'PO4' 'PO6'
    type     (label) object 584B PointType.LANDMARK ... PointType.LANDMARK
Dimensions without coordinates: aligned