Loading and saving data#

Larmor reads each file into one of its types: an image into a Volume, a gradient table into a GradientTable, a tractogram into a Tractogram and a surface into a Surface. Positions are in world coordinates (RAS+, in mm) whatever the file’s own convention, and the end of the file name chooses the format. File formats lists every format, and Conventions explains the coordinates.

The data on this page

The examples on this page run on a small synthetic image, written to the current folder by the code below: dwi.nii.gz with its gradient files dwi.bval and dwi.bvec, and a mask, mask.nii.gz. It has 12 x 12 x 2 voxels of 2 mm, two b=0 images and 30 directions at each of b=1000 and 2000 s/mm^2. Two bundles cross in its middle, one along x and one along y; the other voxels hold isotropic diffusion. To use your own data, use your files instead.

import numpy as np

from larmor.io import save_volume, write_bvals_bvecs
from larmor.sd import AnisotropicRF
from larmor.types import GradientTable, Volume

rng = np.random.default_rng(0)
directions = rng.normal(size=(60, 3))
directions /= np.linalg.norm(directions, axis=1, keepdims=True)
bvals = np.r_[0, 0, [1000] * 30, [2000] * 30]
bvecs = np.r_[np.zeros((2, 3)), directions]
gtab = GradientTable(bvals, bvecs)

fibre = AnisotropicRF([1.7e-3, 0.3e-3, 0.3e-3], mean_kurtosis=0.5)
along_x, along_y = fibre.signal(bvals, gtab.bvecs, np.eye(3)[:2]).T
isotropic = np.exp(-bvals * 0.8e-3)

signal = np.broadcast_to(isotropic, (12, 12, 2, len(bvals))).copy()
signal[:, 4:8] = along_x
signal[4:8, :] = along_y
signal[4:8, 4:8] = (along_x + along_y) / 2
data = 1000 * signal + rng.normal(0, 10, signal.shape)

affine = np.diag([2.0, 2.0, 2.0, 1.0])
save_volume(Volume(data.astype(np.float32), affine), "dwi.nii.gz")
write_bvals_bvecs(gtab, "dwi.bval", "dwi.bvec", ref_path="dwi.nii.gz")
save_volume(Volume(np.ones((12, 12, 2), np.uint8), affine), "mask.nii.gz")

Images#

from larmor.io import load_header, load_volume, save_volume

mask = load_volume("mask.nii.gz")
print(mask.data.shape, mask.voxel_sizes, mask.orientation)
save_volume(mask, "mask.mif")

grid, shape = load_header("dwi.nii.gz")
print(shape)

lar info mask.nii.gz
lar convert mask.nii.gz mask.mif

A Volume holds the stored values (data) and the affine from voxel indices to world coordinates (affine). as_float() gives the values with the file’s scaling applied. load_header() reads only the grid and the shape, without the data. NIfTI (.nii, .nii.gz) and MRtrix (.mif, .mif.gz, .mih) images are read and written. lar convert --to-ras also reorients an image to its closest RAS orientation, which only permutes and flips the voxel axes.

A diffusion image with its gradients#

load_dwi() reads a diffusion image together with its gradient table. A NIfTI image needs its bvals and bvecs files; a MIF image can hold the table in its header and then needs none.

from larmor.io import load_dwi, save_volume

dwi, gtab = load_dwi("dwi.nii.gz", "dwi.bval", "dwi.bvec")
save_volume(dwi, "dwi.mif", gradients=gtab)
dwi, gtab = load_dwi("dwi.mif")
print(dwi.data.shape, gtab.shells)

lar convert dwi.nii.gz dwi.mif --bvals dwi.bval --bvecs dwi.bvec
lar info dwi.mif

Each lar command that takes a diffusion image has the options --bvals and --bvecs, which a MIF image with a table does not need. Gradient tables and shells covers the tables themselves.

Tractograms#

A tractogram holds all the points of its streamlines in one array, in world coordinates, and the offsets where each streamline starts. Three straight streamlines, saved as a TCK file:

import numpy as np

from larmor.io import save_tractogram
from larmor.types import Tractogram

line = np.linspace([2.0, 0.0, 2.0], [22.0, 0.0, 2.0], 11)
points = np.concatenate([line + [0.0, y, 0.0] for y in (8.0, 10.0, 12.0)])
save_tractogram(Tractogram(points, np.array([0, 11, 22, 33])), "tracks.tck")

TCK and VTK files store world coordinates, while a TRK file stores its points in the voxel frame of an image, so writing a TRK needs that image’s grid.

from larmor.io import load_grid, load_tractogram, save_tractogram

tracts = load_tractogram("tracks.tck")
print(tracts.num_streamlines, tracts.num_points)
save_tractogram(tracts, "tracks.trk", reference=load_grid("dwi.nii.gz"))
save_tractogram(tracts, "tracks.vtk")

lar info tracks.tck
lar convert tracks.tck tracks.trk --reference dwi.nii.gz
lar convert tracks.tck tracks.vtk

For large files, load_tractogram(path, max_streamlines=n, seed=0) keeps a random subset of n streamlines, and iter_tractogram() reads an uncompressed TCK or TRK file batch by batch.

Surfaces#

load_surface() and save_surface() read and write GIfTI and FreeSurfer surfaces; see Surfaces.

See also#