cedalion.dataclasses.geometry

Dataclasses for representing geometric objects.

Functions

affine_transform_from_numpy(transform, ...)	Wrap a 4×4 numpy array as a cedalion AffineTransform.
robust_geodesic_distance(surface, seed_vertices)	Compute geodesic distances from seed vertices robustly across multiple submeshes.

Classes

PointType(value)	Categorical label for points in a LabeledPoints array.
PycortexSurface(mesh, crs, units)	Provides the geodesic functionality of Pycortex.
SimpleMesh(pts, polys)	Lightweight mesh container for vertex positions and face indices.
Surface(mesh, crs, units, vertex_coords, ...)	Abstract base class for 3-D triangulated surfaces.
TrimeshSurface(mesh, crs, units, ...)	A surface represented by a trimesh object.
VTKSurface(mesh, crs, units, vertex_coords, ...)	A surface backed by a VTK vtkPolyData mesh.
Voxels(voxels, crs, units)	3D voxels represented by a np.array.

class cedalion.dataclasses.geometry.PointType(value)

Bases: Enum

Categorical label for points in a LabeledPoints array.

Used as the "type" coordinate to distinguish optode types, anatomical landmarks, and EEG electrodes within the same point cloud.

UNKNOWN = 0

SOURCE = 1

DETECTOR = 2

LANDMARK = 3

ELECTRODE = 4

class cedalion.dataclasses.geometry.Surface(

mesh: ~typing.Any,

crs: str,

units: ~pint.registry.Unit,

vertex_coords: dict[str,

~numpy._typing._array_like._Buffer | ~numpy._typing._array_like._SupportsArray[~numpy.dtype[~typing.Any]] | ~numpy._typing._nested_sequence._NestedSequence[~numpy._typing._array_like._SupportsArray[~numpy.dtype[~typing.Any]]] | bool | int | float | complex | str | bytes | ~numpy._typing._nested_sequence._NestedSequence[bool | int | float | complex | str | bytes]] = <factory>,

)

Bases: ABC

Abstract base class for 3-D triangulated surfaces.

Concrete subclasses (TrimeshSurface, VTKSurface, PycortexSurface) wrap different mesh backends while sharing a common interface for coordinate access, KD-tree queries, and affine transforms.

mesh

Backend mesh object (type depends on subclass).

Type:

Any

crs

Name of the coordinate reference system (e.g. "ras").

Type:

str

units

Physical units of the vertex coordinates.

Type:

pint.registry.Unit

vertex_coords

Optional dict of extra per-vertex coordinate arrays (e.g. parcel labels) to attach as xarray coordinates.

Type:

dict[str, numpy._typing._array_like._Buffer | numpy._typing._array_like._SupportsArray[numpy.dtype[Any]] | numpy._typing._nested_sequence._NestedSequence[numpy._typing._array_like._SupportsArray[numpy.dtype[Any]]] | bool | int | float | complex | str | bytes | numpy._typing._nested_sequence._NestedSequence[bool | int | float | complex | str | bytes]]

mesh: Any

crs: str

units: Unit

vertex_coords: dict[str, _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes]]

abstract property vertices: Annotated[DataArray, DataArraySchema(dims='label', coords='label', 'label', 'type')]

Vertices of the mesh as a LabeledPoints array.

abstract property nvertices: int

Number of vertices in the mesh.

abstract property nfaces: int

Number of triangular faces in the mesh.

abstractmethod copy()

abstractmethod apply_transform(transform: DataArray)

Return a new surface with all vertices transformed by transform.

Parameters:

transform – 4×4 affine AffineTransform.

Returns:

A new surface instance of the same concrete type.

property kdtree

KD-tree built lazily from vertex coordinates for fast spatial queries.

snap(

points: Annotated[DataArray, DataArraySchema(dims='label', coords='label', 'label', 'type')],

)

Project points onto their nearest vertex on this surface.

Parameters:

points – LabeledPoints array in the same CRS and units as this surface.

Returns:

LabeledPoints array with the same labels as points but coordinates replaced by those of the closest surface vertex.

Raises:

• CRSMismatchError – If points are not in the same CRS as the surface.

• ValueError – If points have different units from the surface.

class cedalion.dataclasses.geometry.Voxels(voxels: ndarray, crs: str, units: Unit)

Bases: object

3D voxels represented by a np.array.

voxels

The voxels.

Type:

np.ndarray

crs

The coordinate reference system of the voxels.

Type:

str

units

The units of the voxels.

Type:

pint.Unit

voxels: ndarray

crs: str

units: Unit

property vertices: Annotated[DataArray, DataArraySchema(dims='label', coords='label', 'label', 'type')]

property nvertices: int

apply_transform(

transform: DataArray,

) → Voxels

property kdtree

class cedalion.dataclasses.geometry.TrimeshSurface(

mesh: ~trimesh.base.Trimesh,

crs: str,

units: ~pint.registry.Unit,

vertex_coords: dict[str,

~numpy._typing._array_like._Buffer | ~numpy._typing._array_like._SupportsArray[~numpy.dtype[~typing.Any]] | ~numpy._typing._nested_sequence._NestedSequence[~numpy._typing._array_like._SupportsArray[~numpy.dtype[~typing.Any]]] | bool | int | float | complex | str | bytes | ~numpy._typing._nested_sequence._NestedSequence[bool | int | float | complex | str | bytes]] = <factory>,

)

Bases: Surface

A surface represented by a trimesh object.

mesh

The trimesh object representing the surface.

Type:

trimesh.Trimesh

crs

The coordinate reference system of the surface.

Type:

str

units

The units of the surface.

Type:

pint.Unit

mesh: Trimesh

property vertices: Annotated[DataArray, DataArraySchema(dims='label', coords='label', 'label', 'type')]

Vertices of the mesh as a LabeledPoints array.

property nvertices: int

Number of vertices in the mesh.

property nfaces: int

Number of triangular faces in the mesh.

copy() → TrimeshSurface

apply_transform(

transform: DataArray,

) → TrimeshSurface

Apply an affine transformation to this surface.

Parameters:

transform (cdt.AffineTransform) – The affine transformation to apply.

Returns:

The transformed surface.

Return type:

TrimeshSurface

decimate(

face_count: int,

) → TrimeshSurface

Use quadric decimation to reduce the number of vertices.

Parameters:

face_count – the number of faces of the decimated mesh

Returns:

The surface with a decimated mesh

smooth(lamb: float) → TrimeshSurface

Apply a Taubin smoothing filter to the mesh.

Parameters:

lamb – Taubin lambda parameter controlling the amount of smoothing (typically between 0 and 1).

Returns:

A new TrimeshSurface with a smoothed mesh.

get_vertex_normals(

points: Annotated[DataArray, DataArraySchema(dims='label', coords='label', 'label', 'type')],

normalized=True,

)

Return vertex normals at surface vertices nearest to points.

Parameters:

• points – LabeledPoints query positions in the same CRS and units as this surface.

• normalized – If True (default), return unit-length normals.

Returns:

DataArray of shape (n_points, 3) containing the normal vectors at the nearest surface vertices.

Raises:

ValueError – If normalization is requested but any normal has zero length.

fix_vertex_normals()

Flip any vertex normals that point inward (towards the mesh centroid).

Returns:

A new TrimeshSurface with consistently outward-pointing vertex normals.

classmethod from_vtksurface(vtk_surface: VTKSurface)

crs: str

units: Unit

vertex_coords: dict[str, ArrayLike]

class cedalion.dataclasses.geometry.VTKSurface(

mesh: ~vtkmodules.vtkCommonDataModel.vtkPolyData,

crs: str,

units: ~pint.registry.Unit,

vertex_coords: dict[str,

~numpy._typing._array_like._Buffer | ~numpy._typing._array_like._SupportsArray[~numpy.dtype[~typing.Any]] | ~numpy._typing._nested_sequence._NestedSequence[~numpy._typing._array_like._SupportsArray[~numpy.dtype[~typing.Any]]] | bool | int | float | complex | str | bytes | ~numpy._typing._nested_sequence._NestedSequence[bool | int | float | complex | str | bytes]] = <factory>,

)

Bases: Surface

A surface backed by a VTK vtkPolyData mesh.

mesh

The underlying vtk.vtkPolyData object.

Type:

vtkmodules.vtkCommonDataModel.vtkPolyData

crs

Name of the coordinate reference system.

Type:

str

units

Physical units of the vertex coordinates.

Type:

pint.registry.Unit

mesh: vtkPolyData

property vertices: Annotated[DataArray, DataArraySchema(dims='label', coords='label', 'label', 'type')]

Vertices of the mesh as a LabeledPoints array.

property nvertices: int

Number of vertices in the mesh.

property nfaces: int

Number of triangular faces in the mesh.

apply_transform(transform: DataArray)

Return a new surface with all vertices transformed by transform.

Parameters:

transform – 4×4 affine AffineTransform.

Returns:

A new surface instance of the same concrete type.

copy()

classmethod from_trimeshsurface(tri_mesh: TrimeshSurface)

decimate(

reduction: float,

**kwargs,

) → VTKSurface

Use VTK’s decimate_pro method to reduce the number of vertices.

Parameters:

• reduction – Reduction factor. A value of 0.9 will leave 10% of the original number of vertices.

• **kwargs – additional keyword arguments are passed to decimate_pro

Returns:

The surface with a decimated mesh

crs: str

units: Unit

vertex_coords: dict[str, ArrayLike]

class cedalion.dataclasses.geometry.SimpleMesh(pts: ndarray, polys: ndarray)

Bases: object

Lightweight mesh container for vertex positions and face indices.

pts

Vertex positions, shape (n_vertices, 3).

Type:

numpy.ndarray

polys

Triangle face indices, shape (n_faces, 3).

Type:

numpy.ndarray

pts: ndarray

polys: ndarray

class cedalion.dataclasses.geometry.PycortexSurface(

mesh: SimpleMesh,

crs: str,

units: Unit,

)

Bases: Surface

Provides the geodesic functionality of Pycortex.

References

Functions in this class are based on the implementation in Pycortex (Gao et al. [GHLG15]). Gao JS, Huth AG, Lescroart MD and Gallant JL (2015) Pycortex: an interactive surface visualizer for fMRI. Front. Neuroinform. 9:23. doi: 10.3389/fninf.2015.00023

mesh: SimpleMesh

property vertices: Annotated[DataArray, DataArraySchema(dims='label', coords='label', 'label', 'type')]

Vertices of the mesh as a LabeledPoints array.

property nvertices: int

Number of vertices in the mesh.

property nfaces: int

Number of triangular faces in the mesh.

copy()

apply_transform(

transform: DataArray,

) → PycortexSurface

Return a new surface with all vertices transformed by transform.

Parameters:

transform – 4×4 affine AffineTransform.

Returns:

A new surface instance of the same concrete type.

decimate(

face_count: int,

) → PycortexSurface

get_vertex_normals(

points: Annotated[DataArray, DataArraySchema(dims='label', coords='label', 'label', 'type')],

normalized=True,

)

property ppts: ndarray

3D matrix of points in each face.

n faces x 3 per face x 3 coords per point.

property connected: csr_matrix

Sparse matrix of vertex-face associations.

property adj: csr_matrix

Sparse vertex adjacency matrix.

property face_normals: ndarray

Normal vector for each face.

property vertex_normals: ndarray

Normal vector for each vertex (average of normals for neighboring faces).

property face_areas: ndarray

Area of each face.

property cotangent_weights: ndarray

Cotangent of angle opposite each vertex in each face.

property laplace_operator

Laplace-Beltrami operator for this surface.

A sparse adjacency matrix with edge weights determined by the cotangents of the angles opposite each edge. Returns a 4-tuple (B, D, W, V) where D is the ‘lumped mass matrix’, W is the weighted adjacency matrix, and V is a diagonal matrix that normalizes the adjacencies. The ‘stiffness matrix’, A, can be computed as V - W.

The full LB operator can be computed as D^{-1} (V - W).

B is the finite element method (FEM) ‘mass matrix’, which replaces D in FEM analyses.

See ‘Discrete Laplace-Beltrami operators for shape analysis and segmentation’ by Reuter et al., 2009 for details.

property avg_edge_length

Average length of all edges in the surface.

surface_gradient(scalars, at_verts=True)

Gradient of a function with values scalars at each vertex on the surface.

If at_verts, returns values at each vertex. Otherwise, returns values at each face.

Parameters:

• scalars – 1D ndarray, shape (total_verts,) a scalar-valued function across the cortex.

• at_verts – bool, optional If True (default), values will be returned for each vertex. Otherwise, values will be returned for each face.

Returns:

2D ndarray, shape (total_verts,3) or (total_polys,3)

Contains the x-, y-, and z-axis gradients of the given scalars at either each vertex (if at_verts is True) or each face.

Return type:

gradu

geodesic_distance(verts, m=1.0, fem=False)

Calcualte the minimum mesh geodesic distance (in mm).

The geodesic distance is calculated from each vertex in surface to any vertex in the collection verts.

Geodesic distance is estimated using heat-based method (see ‘Geodesics in Heat’, Crane et al, 2012). Diffusion of heat along the mesh is simulated and then used to infer geodesic distance. The duration of the simulation is controlled by the parameter m. Larger values of m will smooth & regularize the distance computation. Smaller values of m will roughen and will usually increase error in the distance computation. The default value of 1.0 is probably pretty good.

This function caches some data (sparse LU factorizations of the laplace-beltrami operator and the weighted adjacency matrix), so it will be much faster on subsequent runs.

The time taken by this function is independent of the number of vertices in verts.

Parameters:

• verts – 1D array-like of ints Set of vertices to compute distance from. This function returns the shortest distance to any of these vertices from every vertex in the surface.

• m – float, optional Reverse Euler step length. The optimal value is likely between 0.5 and 1.5. Default is 1.0, which should be fine for most cases.

• fem – bool, optional Whether to use Finite Element Method lumped mass matrix. Wasn’t used in Crane 2012 paper. Doesn’t seem to help any.

Returns:

1D ndarray, shape (total_verts,) Geodesic distance (in mm) from each vertex in the surface to the closest vertex in verts.

geodesic_path(a, b, max_len=1000, d=None, **kwargs)

Finds the shortest path between two points a and b.

This shortest path is based on geodesic distances across the surface. The path starts at point a and selects the neighbor of a in the graph that is closest to b. This is done iteratively with the last vertex in the path until the last point in the path is b.

Other Parameters in kwargs are passed to the geodesic_distance function to alter how geodesic distances are actually measured

Parameters:

• a – int Vertex that is the start of the path

• b – int Vertex that is the end of the path

• d – array array of geodesic distances, will be computed if not provided

• max_len – int, optional, default=1000 Maximum path length before the function quits. Sometimes it can get stuck in loops, causing infinite paths.

• m – float, optional Reverse Euler step length. The optimal value is likely between 0.5 and 1.5. Default is 1.0, which should be fine for most cases.

• fem – bool, optional Whether to use Finite Element Method lumped mass matrix. Wasn’t used in Crane 2012 paper. Doesn’t seem to help any.

• kwargs – other arugments are passed to self.geodesic_distance

Returns:

list

List of the vertices in the path from a to b

Return type:

path

property graph

NetworkX undirected graph representing this Surface.

property iter_surfedges

classmethod from_trimeshsurface(

tri_mesh: TrimeshSurface,

)

classmethod from_vtksurface(vtk_surface: VTKSurface)

crs: str

units: Unit

vertex_coords: dict[str, ArrayLike]

cedalion.dataclasses.geometry.affine_transform_from_numpy(

transform: ndarray,

from_crs: str,

to_crs: str,

from_units: str,

to_units: str,

) → DataArray

Wrap a 4×4 numpy array as a cedalion AffineTransform.

The resulting DataArray has pint units of to_units / from_units and dimension names [to_crs, from_crs], matching the convention used throughout cedalion for affine transforms.

Parameters:

• transform – 4×4 affine transformation matrix as a numpy array.

• from_crs – Name of the source coordinate reference system.

• to_crs – Name of the target coordinate reference system.

• from_units – Unit string for the source space (e.g. "mm").

• to_units – Unit string for the target space (e.g. "m").

Returns:

A 4×4 AffineTransform DataArray with dims [to_crs, from_crs] and units to_units / from_units.

cedalion.dataclasses.geometry.robust_geodesic_distance(

surface: TrimeshSurface,

seed_vertices: list[int] | int,

m: float = 10.0,

) → ndarray

Compute geodesic distances from seed vertices robustly across multiple submeshes.

Parameters:

• surface (TrimeshSurface) – The surface.

• seed_vertices (list[int] | int) – Seed vertices (single int or list of ints).

• m (float) – Geodesic distance parameter.

Returns:

Distance array of shape (num_vertices,), inf for vertices

not reachable (not on the same submesh as any seed).

Return type:

np.ndarray

Initial Contributors:

• Thomas Fischer | t.fischer.1@campus.tu-berlin.de | 2025