lapy.TriaMesh¶

class lapy.TriaMesh(v, t, fsinfo=None)¶

Class representing a triangle mesh.

This is an efficient implementation of a triangle mesh data structure with core functionality using sparse matrices internally (Scipy).

Parameters:

varray_like: List of lists of 3 float coordinates.
tarray_like: List of lists of 3 int of indices (>= 0) into v array Ordering is important: All triangles should be oriented in the same way (counter-clockwise, when looking from above).
fsinfodict | None: FreeSurfer Surface Header Info.

Attributes

v	(array_like) List of lists of 3 float coordinates.
t	(array_like) List of lists of 3 int of indices (>= 0) into `v` array.
adj_sym	(csc_matrix) Symmetric adjacency matrix as csc sparse matrix.
adj_dir	(csc_matrix) Directed adjacency matrix as csc sparse matrix.
fsinfo	(dict \| None) FreeSurfer Surface Header Info.

Methods

`area`()	Compute the total surface area of triangle mesh.
`avg_edge_length`()	Compute the average edge length of the mesh.
`boundary_loops`()	Compute a tuple of boundary loops.
`centroid`()	Compute centroid of triangle mesh as a weighted average of triangle centers.
`construct_adj_dir_tidx`()	Construct directed adjacency matrix (edge graph) of triangle mesh.
`curvature`([smoothit])	Compute various curvature values at vertices.
`curvature_tria`([smoothit])	Compute min and max curvature and directions (orthogonal and in tria plane).
`edges`([with_boundary])	Compute vertices and adjacent triangle ids for each edge.
`euler`()	Compute the Euler Characteristic.
`has_free_vertices`()	Check if the vertex list has more vertices than what is used in tria.
`is_closed`()	Check if triangle mesh is closed (no boundary edges).
`is_manifold`()	Check if triangle mesh is manifold (no edges with >2 triangles).
`is_oriented`()	Check if triangle mesh is oriented.
`level_length`(vfunc, level)	Compute the length of level sets.
`level_path`(vfunc, level[, get_tria_idx, ...])	Extract levelset as a path.
`map_tfunc_to_vfunc`(tfunc[, weighted])	Map tria function to vertices by attributing 1/3 to each vertex of triangle.
`map_vfunc_to_tfunc`(vfunc)	Map vertex function to triangles by attributing 1/3 to each.
`normal_offset_`(d)	Move vertices along normal by distance `d`.
`normalize_`()	Normalize TriaMesh to unit surface area and centroid at the origin.
`orient_`()	Re-orient triangles of manifold mesh to be consistent.
`read_fssurf`(filename)	Load triangle mesh from FreeSurfer surface geometry file.
`read_off`(filename)	Load triangle mesh from OFF txt file.
`read_vtk`(filename)	Load triangle mesh from VTK txt file.
`refine_`([it])	Refine the triangle mesh by placing new vertex on each edge midpoint.
`rm_free_vertices_`()	Remove unused (free) vertices.
`smooth_`([n])	Smooth the mesh iteratively in-place.
`smooth_vfunc`(vfunc[, n])	Smooth the mesh or a vertex function iteratively.
`tria_areas`()	Compute the area of triangles using Heron's formula.
`tria_normals`()	Compute triangle normals.
`tria_qualities`()	Compute triangle quality for each triangle in mesh where.
`vertex_areas`()	Compute the area associated to each vertex (1/3 of one-ring trias).
`vertex_degrees`()	Compute the vertex degrees (number of edges at each vertex).
`vertex_normals`()	Compute vertex normals.
`volume`()	Compute the volume of closed triangle mesh, summing tetrahedra at origin.
`write_fssurf`(filename)	Save as Freesurfer Surface Geometry file (wrap Nibabel).
`write_vtk`(filename)	Save as VTK file.

Notes

The class has static class methods to read triangle meshes from FreeSurfer, OFF, and VTK file formats.

area()¶

Compute the total surface area of triangle mesh.

Returns:

float: Total surface area.

avg_edge_length()¶

Compute the average edge length of the mesh.

Returns:

float: Average edge length.

boundary_loops()¶

Compute a tuple of boundary loops.

Meshes can have 0 or more boundary loops, which are cycles in the directed adjacency graph of the boundary edges. Works on trias only. Could fail if loops are connected via a single vertex (like a figure 8). That case needs debugging.

Returns:

loopslist of list: List of lists with boundary loops.

centroid()¶

Compute centroid of triangle mesh as a weighted average of triangle centers.

The weight is determined by the triangle area. (This could be done much faster if a FEM lumped mass matrix M is already available where this would be M*v, because it is equivalent with averaging vertices weighted by vertex area)

Returns:

centroidfloat: The centroid of the mesh.
totalareafloat: The total area of the mesh.

construct_adj_dir_tidx()¶

Construct directed adjacency matrix (edge graph) of triangle mesh.

The directed adjacency matrix contains the triangle indices (only for non-manifold meshes). Operates only on triangles.

Returns:

csc_matrix: Similar to adj_dir, but stores the tria idx+1 instead of one in the matrix (allows lookup of vertex to tria).

curvature(smoothit=3)¶

Compute various curvature values at vertices.

Note

For the algorithm see e.g. Pierre Alliez, David Cohen-Steiner, Olivier Devillers, Bruno Levy, and Mathieu Desbrun. Anisotropic Polygonal Remeshing. ACM Transactions on Graphics, 2003.

Parameters:

smoothitint: Smoothing iterations on vertex functions.

Returns:

u_minarray of shape (vnum, 3): Minimal curvature directions.
u_maxarray of shape (vnum, 3): Maximal curvature directions.
c_minarray: Minimal curvature.
c_maxarray: Maximal curvature.
c_meanarray: Mean curvature (c_min + c_max) / 2.0m.
c_gaussarray: Gauss curvature c_min * c_maxm.
normalsarray of shape (vnum, 3): Normals.

curvature_tria(smoothit=3)¶

Compute min and max curvature and directions (orthogonal and in tria plane).

First we compute these values on vertices and then smooth there. Finally they get mapped to the trias (averaging) and projected onto the triangle plane, and orthogonalized.

Parameters:

smoothitint: Number of smoothing iterations for curvature computation on vertices.

Returns:

u_minarray: Min curvature direction on triangles.
u_maxarray: Max curvature direction on triangles.
c_minarray: Min curvature on triangles.
c_maxarray: Max curvature on triangles.

edges(with_boundary=False)¶

Compute vertices and adjacent triangle ids for each edge.

Parameters:

with_boundarybool: Also work on boundary half edges, default ignore.

Returns:

vidsarray: Column array with starting and end vertex for each unique inner edge.
tidsarray: 2 column array with triangle containing the half edge from vids[0,:] to vids [1,:] in first column and the neighboring triangle in the second column.
bdrvidsarray: If with_boundary is true: 2 column array with each boundary half-edge.
bdrtidsarray: If with_boundary is true: 1 column array with the associated triangle to each boundary edge.

euler()¶

Compute the Euler Characteristic.

The Euler characteristic is the number of vertices minus the number of edges plus the number of triangles (= #V - #E + #T). For example, it is 2 for the sphere and 0 for the torus. This operates only on triangles array.

Returns:

int: Euler characteristic.

has_free_vertices()¶

Check if the vertex list has more vertices than what is used in tria.

Returns:

bool: Whether vertex list has more vertices or not.

is_closed()¶

Check if triangle mesh is closed (no boundary edges).

Operates only on triangles

Returns:

bool: True if no boundary edges in adj matrix.

is_manifold()¶

Check if triangle mesh is manifold (no edges with >2 triangles).

Operates only on triangles

Returns:

bool: True if no edges with > 2 triangles.

is_oriented()¶

Check if triangle mesh is oriented.

True if all triangles are oriented counter-clockwise, when looking from above. Operates only on triangles.

Returns:

bool: True if max(adj_directed)=1.

level_length(vfunc, level)¶

Compute the length of level sets.

Parameters:

vfuncarray: Float vector of values at vertices (here only scalar function 1D).
levelfloat | array: Level set value or array of level values.

Returns:

lengthfloat | array: Length of level set (or array of lengths).

level_path(vfunc, level, get_tria_idx=False, n_points=None)¶

Extract levelset as a path.

Note: Only works for level sets that represent a single path with exactly one start and one endpoint.

Parameters:

vfuncarray: Float vector of values at vertices (here only scalar function 1D).
levelfloat: Level set value.
get_tria_idxbool, default: False: Also return a list of triangle indices for each edge, default False.
n_pointsint: Resample level set into n equidistant points. Cannot be combined with get_tria_idx=True.

Returns:

path: array: Array of shape (n,3) containing coordinates of vertices on level path.
lengthfloat: Length of level set.
tria_idxarray: Array of triangle index for each edge segment on the path (length n-1 if path is length n). Will only be returned if get_tria_idx is True.

map_tfunc_to_vfunc(tfunc, weighted=False)¶

Map tria function to vertices by attributing 1/3 to each vertex of triangle.

Uses vertices and trias.

Parameters:

tfuncarray: Float vector or matrix (#t x N) of values at vertices.
weightedbool, default=False: False, weigh only by 1/3, e.g. to compute vertex areas from tria areas True, weigh by triangle area / 3, e.g. to integrate a function defined on the trias, for example integrating the “one” function will also yield the vertex areas.

Returns:

vfuncarray: Function on vertices vector or matrix (#v x N).

map_vfunc_to_tfunc(vfunc)¶

Map vertex function to triangles by attributing 1/3 to each.

Uses number of vertices and trias

Parameters:

vfuncarray: Float vector or matrix (#t x N) of values at vertices.

Returns:

tfuncarray: Function on trias vector or matrix (#t x N).

normal_offset_(d)¶

Move vertices along normal by distance d.

normal_offset(d) moves vertices along normal by distance d

Parameters:

dint | array: Move distance.

normalize_()¶

Normalize TriaMesh to unit surface area and centroid at the origin.

Modifies the vertices.

orient_()¶

Re-orient triangles of manifold mesh to be consistent.

Re-orients triangles of manifold mesh to be consistent, so that vertices are listed counter-clockwise, when looking from above (outside).

Algorithm:

Construct list for each half-edge with its triangle and edge direction
Drop boundary half-edges and find half-edge pairs
Construct sparse matrix with triangle neighbors, with entry 1 for opposite half edges and -1 for parallel half-edges (normal flip across this edge)
Flood mesh from first tria using triangle neighbor matrix and keeping track of sign
When flooded, negative sign for a triangle indicates it needs to be flipped
If global volume is negative, flip everything (first tria was wrong)

Returns:

flippedint: Number of trias flipped.

classmethod read_fssurf(filename)¶

Load triangle mesh from FreeSurfer surface geometry file.

Parameters:

filenamestr: Filename to load.

Returns:

TriaMesh: Loaded triangle mesh.

classmethod read_off(filename)¶

Load triangle mesh from OFF txt file.

Parameters:

filenamestr: Filename to load.

Returns:

TriaMesh: Loaded triangle mesh.

classmethod read_vtk(filename)¶

Load triangle mesh from VTK txt file.

Parameters:

filenamestr: Filename to load.

Returns:

TriaMesh: Loaded triangle mesh.

refine_(it=1)¶

Refine the triangle mesh by placing new vertex on each edge midpoint.

Thus creates 4 similar triangles from one parent triangle. Modifies mesh in place.

Parameters:

itint: Number of iterations.

rm_free_vertices_()¶

Remove unused (free) vertices.

Free vertices are vertices that are not used in any triangle. They can produce problems when constructing, e.g., Laplace matrices.

Will update v and t in mesh.

Returns:

vkeeparray: Indices (from original list) of kept vertices.
vdelarray: Indices of deleted (unused) vertices.

smooth_(n=1)¶

Smooth the mesh iteratively in-place.

Parameters:

nint: Smoothing iterations.

smooth_vfunc(vfunc, n=1)¶

Smooth the mesh or a vertex function iteratively.

Parameters:

vfuncarray: Float vector of values at vertices, if empty, use vertex coordinates.
nint, default=1: Number of iterations for smoothing.

Returns:

vfuncarray: Smoothed surface vertex function.

tria_areas()¶

Compute the area of triangles using Heron’s formula.

Heron’s formula computes the area of a triangle by using the three edge lengths.

Returns:

areasarray: Array with areas of each triangle.

tria_normals()¶

Compute triangle normals.

Ordering of triangles is important: counterclockwise when looking.

Returns:

narray of shape (n_triangles, 3): Triangle normals.

tria_qualities()¶

Compute triangle quality for each triangle in mesh where.

q = 4 sqrt(3) A / (e1^2 + e2^2 + e3^2 ) where A is the triangle area and ei the edge length of the three edges. Constants are chosen so that q=1 for the equilateral triangle.

Note

This measure is used by FEMLAB and can also be found in: R.E. Bank, PLTMG …, Frontiers in Appl. Math. (7), 1990.

Returns:

array: Array with triangle qualities.

vertex_areas()¶

Compute the area associated to each vertex (1/3 of one-ring trias).

Returns:

vareasarray: Array of vertex areas.

vertex_degrees()¶

Compute the vertex degrees (number of edges at each vertex).

Returns:

vdegarray: Array of vertex degrees.

vertex_normals()¶

Compute vertex normals.

get_vertex_normals(v,t) computes vertex normals: Triangle normals around each vertex are averaged, weighted by the angle that they contribute. Ordering is important: counterclockwise when looking at the triangle from above.

Returns:

narray of shape (n_triangles, 3): Vertex normals.

volume()¶

Compute the volume of closed triangle mesh, summing tetrahedra at origin.

Returns:

volfloat: Total enclosed volume.

write_fssurf(filename)¶

Save as Freesurfer Surface Geometry file (wrap Nibabel).

Parameters:

filenamestr: Filename to save to.

write_vtk(filename)¶

Save as VTK file.

Parameters:

filenamestr: Filename to save to.