Marching Cubes

Marching cubes is the workhorse of iso-surface extraction: it visits every grid cell and cuts triangles out of it based on which corners are inside the surface. It is the default choice for most fields.

import torch
import isoext
from isoext import viewer

# Setup: create a grid with a sphere
grid = isoext.UniformGrid([64, 64, 64])
points = grid.get_points()
grid.set_values(points.norm(dim=-1) - 0.7)

Basic Usage

vertices, faces = isoext.marching_cubes(grid, level=0.0, method="vega")
  • grid: A UniformGrid or SparseGrid with values set

  • level: The iso-value to extract (default: 0.0)

  • method: Algorithm variant ("vega", "lewiner", "nagae", or "lorensen")

vertices, faces = isoext.marching_cubes(grid)

print(f"Vertices: {vertices.shape}")  # (N, 3) float32
print(f"Faces: {faces.shape}")  # (M, 3) uint32

viewer.embed(vertices, faces)
Vertices: torch.Size([9168, 3])
Faces: torch.Size([18332, 3])

How It Works

Marching cubes visits every cell independently. Each of a cell’s 8 corners is either inside or outside the surface, giving 256 possible sign patterns, and a lookup table maps each pattern to a set of triangles. The triangle vertices sit on the cell’s edges. If one end of an edge is inside the surface and the other is outside, the surface has to cross that edge somewhere in between, and interpolating the two values estimates where: with level-relative endpoint values \(v_a\) and \(v_b\), the crossing sits at the fraction \(t = v_a / (v_a - v_b)\) of the way along the edge. An endpoint whose value is close to zero is close to the surface, so it pulls the crossing toward itself.

The simplest pattern has one corner inside, the red dot below. The surface crosses the three edges that connect it to outside corners; the gold dots mark those crossings. For this pattern the lookup table produces one triangle, and the three crossings are its vertices:

cell = isoext.UniformGrid([2, 2, 2])
cell.set_values(torch.tensor(
    [-0.5, 0.3, 0.4, 0.5, 0.6, 0.4, 0.3, 0.2], device="cuda"
).reshape(2, 2, 2))

crossings = isoext.get_intersection(cell).get_points()
v, f = isoext.marching_cubes(cell)
print(f"{f.shape[0]} triangle(s)")
viewer.embed(v, f, color="steelblue", flat_shading=True, side="double", height=300,
             grid=cell, draw=lambda s: viewer.add_points(s, crossings, point_size=0.09))
1 triangle(s)

Algorithm Variants

The method argument selects between four variants that share the same interface:

  • vega (default) resolves topological ambiguities by following the field’s interpolant, using the corrected interior test.

  • lewiner follows the interpolant with the older interior test; kept for comparison with other implementations.

  • nagae uses fixed tables; always closed and slightly faster.

  • lorensen is the original 1987 algorithm; ambiguous cells can leave small cracks.

See Marching Cubes Variants for the differences in detail.

for method in ["nagae", "lorensen", "lewiner", "vega"]:
    v, f = isoext.marching_cubes(grid, method=method)
    print(f"{method:9s} {f.shape[0]:,} triangles")
nagae     18,332 triangles
lorensen  18,332 triangles
lewiner   18,332 triangles
vega      18,332 triangles

Iso-Level

The level parameter controls which iso-surface to extract. For signed distance fields, level=0 gives the surface. Other values give offset surfaces:

# Extract at different iso-levels
for level in [-0.2, 0.0, 0.2]:
    v, f = isoext.marching_cubes(grid, level=level)
    print(f"level={level:+.1f}: {v.shape[0]:,} vertices (radius ≈ {0.7 - level:.1f})")
level=-0.2: 4,728 vertices (radius ≈ 0.9)
level=+0.0: 9,168 vertices (radius ≈ 0.7)
level=+0.2: 15,072 vertices (radius ≈ 0.5)

Saving Meshes

Use write_obj to save the mesh to an OBJ file:

vertices, faces = isoext.marching_cubes(grid)
isoext.write_obj("mesh.obj", vertices, faces)
print("Saved to mesh.obj")
Saved to mesh.obj