Introduction

When a hospital MRI scanner finishes its work, the result is not a picture — it is a grid of numbers. Each tiny voxel (a 3D pixel) stores a density value: how much the tissue at that point resisted the magnetic field. Bone is dense; air is not; soft tissue lives somewhere between.

To make the scan useful, doctors need a surface: the boundary between one material and another. Draw a threshold value and ask — which voxels are above it and which are below? The surface you want is the shape that separates them.

That problem looks impossibly hard at first. A real scan has hundreds of millions of voxels; the surface can be any shape imaginable. But in 1987, William Lorensen and Harvey Cline discovered a beautiful shortcut: no matter how complicated the surface is, each individual cube of eight voxels can only be in one of 256 possible states. Pre-compute what triangles to draw for each state, store the answer in a lookup table, and then march through the grid one cube at a time — that is the Marching Cubes algorithm.

It is fast, it is simple, and it powered the very first real-time 3D reconstructions of medical scans. Today it still runs — in a slightly refined form — under nearly every 3D volume viewer in medicine, science, and game development.

Build a Surface Live

The demo below shows a 2D version of the algorithm (called marching squares) on a small scalar field. Each cell in the grid stores a value; the color shows how large it is. Drag the threshold slider and watch the algorithm classify each cell, look up the right edge-crossing pattern, and stitch the contour together.

<div class="controls">
  <label for="thresh">Threshold: <span id="thresh-val">0.5</span></label>
  <input type="range" id="thresh" min="0" max="1" step="0.01" value="0.5">
  <button id="animate" type="button">▶ Animate step by step</button>
  <button id="reset-btn" type="button" class="ghost">Reset</button>
</div>
<canvas id="canvas" width="420" height="280"></canvas>
<div id="info" class="info">Drag the threshold slider to extract a contour.</div>

* { box-sizing: border-box; margin: 0; padding: 0; }
body { font-family: system-ui, sans-serif; color: #222; background: #fff; }
.controls { display: flex; flex-wrap: wrap; align-items: center; gap: .5rem .8rem; margin-bottom: .6rem; }
label { font-size: .85rem; font-weight: 600; white-space: nowrap; }
input[type=range] { width: 140px; accent-color: #1d3557; }
button { font: 600 13px system-ui; padding: .38rem .8rem; border: 1px solid #1d3557;
         background: #1d3557; color: #fff; border-radius: 7px; cursor: pointer; }
button.ghost { background: #fff; color: #1d3557; }
canvas { display: block; border: 1px solid #d0d8e0; border-radius: 8px; max-width: 100%; }
.info { font-size: .82rem; color: #556; margin-top: .5rem; min-height: 1.3em; }

const fmt = (s, v) => s.replace(/%(\w+)%/g, (_, k) => String(v[k] ?? ''));
// Marching Squares demo — 2D analogue of Marching Cubes.
// A 7×5 scalar field; each cell has 4 corners. We classify each corner as
// "inside" (value >= threshold) or "outside" (value < threshold), form a 4-bit index (0..15), look up which edges are crossed, interpolate the crossing point on each edge, and draw contour line segments.
// 
//
// Edge numbering per cell (corners: 0=TL 1=TR 2=BR 3=BL):
//   edge 0: top    (corners 0-1)
//   edge 1: right  (corners 1-2)
//   edge 2: bottom (corners 2-3)
//   edge 3: left   (corners 3-0)
//
// LOOKUP TABLE: for each of 16 cases, list pairs of edges to connect.
// (symmetric ambiguous cases 6 and 9 use the simpler saddle-point resolution)
const EDGES = [
  [],           // 0000 — all outside
  [[2,3]],      // 0001 — BL inside
  [[1,2]],      // 0010 — BR inside
  [[1,3]],      // 0011 — BL+BR inside
  [[0,1]],      // 0100 — TR inside
  [[0,3],[1,2]],// 0101 — TR+BL (ambiguous — two segments)
  [[0,2]],      // 0110 — TR+BR inside
  [[0,3]],      // 0111 — TR+BR+BL inside
  [[0,3]],      // 1000 — TL inside
  [[0,2]],      // 1001 — TL+BL inside
  [[0,1],[2,3]],// 1010 — TL+BR (ambiguous)
  [[0,1]],      // 1011 — TL+BR+BL inside
  [[1,3]],      // 1100 — TL+TR inside
  [[1,2]],      // 1101 — TL+TR+BL inside
  [[2,3]],      // 1110 — TL+TR+BR inside
  []            // 1111 — all inside
];

// A 8×6 grid of scalar values (rows×cols of corners = (ROWS+1)×(COLS+1))
const COLS = 7, ROWS = 5;
// Generate a smooth-ish field using a sum of Gaussian bumps
const FIELD = [];
const bumps = [
  {cx:2.5,cy:1.5,a:0.9,sx:2,sy:1.5},
  {cx:5,  cy:3.5,a:0.8,sx:1.8,sy:1.8},
  {cx:1,  cy:3,  a:0.7,sx:1.5,sy:1.5},
  {cx:6,  cy:1,  a:0.6,sx:1.2,sy:1.2}
];
for (let r = 0; r <= ROWS; r++) {
  FIELD[r] = [];
  for (let c = 0; c <= COLS; c++) {
    let v = 0;
    for (const b of bumps) {
      const dx = (c - b.cx) / b.sx, dy = (r - b.cy) / b.sy;
      v += b.a * Math.exp(-(dx*dx + dy*dy));
    }
    FIELD[r][c] = Math.min(1, Math.max(0, v));
  }
}

const canvas = document.getElementById('canvas');
const ctx = canvas.getContext('2d');
const W = canvas.width, H = canvas.height;
const PAD = 18;
const CW = (W - 2*PAD) / COLS;   // cell width
const CH = (H - 2*PAD) / ROWS;   // cell height

function cornerX(col) { return PAD + col * CW; }
function cornerY(row) { return PAD + row * CH; }

// Interpolate crossing point on an edge given threshold t
function interp(r0,c0,r1,c1,t) {
  const v0 = FIELD[r0][c0], v1 = FIELD[r1][c1];
  const frac = (v0 === v1) ? 0.5 : (t - v0) / (v1 - v0);
  return {
    x: cornerX(c0) + frac * (cornerX(c1) - cornerX(c0)),
    y: cornerY(r0) + frac * (cornerY(r1) - cornerY(r0))
  };
}

// Edge index -> pair of corner (row,col) offsets within a cell (row, col)
// Corners: 0=TL(r,c) 1=TR(r,c+1) 2=BR(r+1,c+1) 3=BL(r+1,c)
const CORNER = [[0,0],[0,1],[1,1],[1,0]];
const EDGE_CORNERS = [[0,1],[1,2],[2,3],[3,0]];

function edgePoint(cellRow, cellCol, edgeIdx, t) {
  const [ca, cb] = EDGE_CORNERS[edgeIdx];
  const [ra, ca2] = CORNER[ca], [rb, cb2] = CORNER[cb];
  return interp(cellRow+ra, cellCol+ca2, cellRow+rb, cellCol+cb2, t);
}

let threshold = 0.5;
let animStep = -1;  // -1 = show all
let animTimer = null;

function valueColor(v) {
  // Blue (low) to orange (high)
  const r = Math.round(30  + v * 200);
  const g = Math.round(80  + v * 100);
  const b = Math.round(180 - v * 150);
  return `rgb(${r},${g},${b})`;
}

function draw(step) {
  ctx.clearRect(0, 0, W, H);

  // Draw cells with background color = average corner value
  for (let r = 0; r < ROWS; r++) {
    for (let c = 0; c < COLS; c++) {
      const avg = (FIELD[r][c]+FIELD[r][c+1]+FIELD[r+1][c]+FIELD[r+1][c+1])/4;
      ctx.fillStyle = valueColor(avg);
      ctx.globalAlpha = 0.55;
      ctx.fillRect(cornerX(c)+1, cornerY(r)+1, CW-2, CH-2);
      ctx.globalAlpha = 1;
    }
  }

  // Draw grid lines
  ctx.strokeStyle = '#8899aa';
  ctx.lineWidth = 0.5;
  for (let r = 0; r <= ROWS; r++) {
    ctx.beginPath(); ctx.moveTo(cornerX(0), cornerY(r)); ctx.lineTo(cornerX(COLS), cornerY(r)); ctx.stroke();
  }
  for (let c = 0; c <= COLS; c++) {
    ctx.beginPath(); ctx.moveTo(cornerX(c), cornerY(0)); ctx.lineTo(cornerX(c), cornerY(ROWS)); ctx.stroke();
  }

  // Draw corner dots colored by inside/outside
  for (let r = 0; r <= ROWS; r++) {
    for (let c = 0; c <= COLS; c++) {
      const inside = FIELD[r][c] >= threshold;
      ctx.beginPath();
      ctx.arc(cornerX(c), cornerY(r), 4, 0, 2*Math.PI);
      ctx.fillStyle = inside ? '#e63946' : '#adb8c5';
      ctx.fill();
    }
  }

  // Draw contour segments (all or up to step)
  let cellIdx = 0;
  let drawnCount = 0;
  for (let r = 0; r < ROWS; r++) {
    for (let c = 0; c < COLS; c++) {
      const highlight = (step >= 0 && cellIdx === step);
      const process   = (step < 0 || cellIdx <= step);
      cellIdx++;

      // Compute index
      const tl = FIELD[r][c]   >= threshold ? 8 : 0;
      const tr = FIELD[r][c+1] >= threshold ? 4 : 0;
      const br = FIELD[r+1][c+1] >= threshold ? 2 : 0;
      const bl = FIELD[r+1][c] >= threshold ? 1 : 0;
      const idx = tl|tr|br|bl;

      if (highlight) {
        ctx.fillStyle = 'rgba(255,220,50,0.35)';
        ctx.fillRect(cornerX(c)+1, cornerY(r)+1, CW-2, CH-2);
        // Show index label
        ctx.fillStyle = '#1d3557';
        ctx.font = 'bold 11px monospace';
        ctx.textAlign = 'center';
        ctx.fillText(fmt('case %idx%', { idx: idx }), cornerX(c)+CW/2, cornerY(r)+CH/2);
      }

      if (!process) continue;

      const segs = EDGES[idx];
      for (const [e0,e1] of segs) {
        const p0 = edgePoint(r, c, e0, threshold);
        const p1 = edgePoint(r, c, e1, threshold);
        ctx.beginPath();
        ctx.moveTo(p0.x, p0.y);
        ctx.lineTo(p1.x, p1.y);
        ctx.strokeStyle = '#e63946';
        ctx.lineWidth = 2.5;
        ctx.stroke();
        drawnCount++;
      }
    }
  }

  const totalCells = ROWS * COLS;
  if (step < 0) {
    document.getElementById('info').textContent =
      fmt('All %totalCells% cells processed. %drawnCount% contour segment(s) drawn.', { totalCells: totalCells, drawnCount: drawnCount });
  } else {
    document.getElementById('info').textContent =
      fmt('Step %step%%v% / %totalCells%  — cell highlighted in yellow (case index: %v2%)', { step: step, v: 1, totalCells: totalCells, v2: (() => {
        const r = Math.floor(step / COLS), c = step % COLS;
        const tl = FIELD[r][c]   >= threshold ? 8 : 0;
        const tr = FIELD[r][c+1] >= threshold ? 4 : 0;
        const br = FIELD[r+1][c+1] >= threshold ? 2 : 0;
        const bl = FIELD[r+1][c] >= threshold ? 1 : 0;
        return (tl|tr|br|bl);
      })() });
  }
}

function stopAnim() {
  if (animTimer) { clearInterval(animTimer); animTimer = null; }
  document.getElementById('animate').textContent = '▶ Animate step by step';
}

document.getElementById('thresh').addEventListener('input', function() {
  threshold = parseFloat(this.value);
  document.getElementById('thresh-val').textContent = threshold.toFixed(2);
  stopAnim();
  animStep = -1;
  draw(-1);
});

document.getElementById('animate').addEventListener('click', function() {
  if (animTimer) { stopAnim(); animStep = -1; draw(-1); return; }
  animStep = 0;
  this.textContent = '■ Stop';
  draw(animStep);
  animTimer = setInterval(() => {
    animStep++;
    if (animStep >= ROWS * COLS) { stopAnim(); animStep = -1; draw(-1); return; }
    draw(animStep);
  }, 180);
});

document.getElementById('reset-btn').addEventListener('click', function() {
  stopAnim();
  animStep = -1;
  threshold = 0.5;
  document.getElementById('thresh').value = '0.5';
  document.getElementById('thresh-val').textContent = '0.50';
  draw(-1);
});

draw(-1);
Result
Console

Notice two things. First, once the threshold is set, each cell is decided independently — the algorithm never backtracks or searches. Second, cells near the boundary flip between "inside" and "outside" as you drag the slider, instantly updating the contour. This is the lookup-table trick: a handful of pre-drawn cases, stamped onto the grid one cell at a time. The 3D version works the same way with 256 cases instead of 16.

The Real Complexity

Marching Cubes is one of the rare algorithms whose complexity is exactly as good as it can possibly be.

  • Each cube is independent. The algorithm looks at the eight corner values of one cube, computes a single 8-bit index (0–255), fetches the pre-stored triangle list from a table, and moves on. No search, no backtracking, no communication with other cubes.
  • O(n)O(n) time, O(1)O(1) per cube. If the volume has n voxels, the algorithm takes O(n)O(n) time in total and constant time per cube. You cannot do better: you must at least look at every voxel once.
  • The lookup table has exactly 256 entries — one per distinct combination of eight binary corner states (inside/outside). By symmetry, many are rotations or reflections of each other; Lorensen and Cline reduced them to 15 canonical cases in their original paper.
  • Ambiguous faces are the only complication. When opposite corners of a face are inside and the other two are outside, two different triangle tilings are both valid, and choosing the wrong one can leave a hole in the surface. Several follow-up algorithms (Marching Tetrahedra, Dual Contouring) trade a small increase in complexity for a guaranteed hole-free result.
  • Parallelism is trivial. Because cubes are independent, marching cubes is embarrassingly parallel — the entire volume can be split across GPU cores, which is why modern volume renderers run it in real time on medical-grade datasets.

Contrast this with problems like P vs NP where checking solutions is easy but finding them is expensive. Marching cubes is the opposite extreme: both finding and checking are linear. The hard part was realizing that 256 cases are enough.

Where It Matters

The 256-case lookup table shows up in a surprisingly wide range of fields:

  • Medical imaging: CT and MRI scanners produce voxel grids; marching cubes turns them into the smooth 3D organ models that surgeons use for pre-operative planning and that radiologists annotate on screen.
  • Scientific visualization: weather simulations, computational fluid dynamics, and molecular dynamics all produce scalar fields — density, pressure, temperature — and marching cubes extracts the surface at any threshold in one linear pass.
  • Video-game terrain and destruction: voxel-based games (Minecraft's smoother cousins, teardown simulations) store the world as a density grid and use marching cubes — or its variants — to generate renderable meshes in real time as terrain changes.
  • Metaballs and implicit surfaces: 3D artists define objects as mathematical formulas ("the set of points closer than r to this skeleton"); marching cubes samples the formula on a grid and extracts the visible surface.
  • 3D printing slicing: converting a solid model into printable layers requires intersecting the model with horizontal planes — essentially marching squares repeated per layer.

The algorithm is closely related to ideas in convex hull computation and mesh generation — anywhere a discrete sample must be turned into a geometric boundary.

Conclusion

Lorensen and Cline's key insight was deceptively simple: however wild the surface, each cube of eight voxels carries only one bit of information per corner, giving at most 256 distinct cases. Pre-tabulate the answer for all 256 and the rest is just iteration — march through the grid, stamp the right triangles, done.

That simplicity is the algorithm's superpower. It runs in linear time, parallelises perfectly, and can reconstruct a surface at any density threshold without ever searching or backtracking. Forty years after its publication, the same idea still drives the 3D models on every radiologist's screen.

The next time you see a smooth 3D scan of a brain or a lung, picture 256 tiny cookie-cutters stamping triangles across millions of cubes — one lookup per cube, no exceptions. That is marching cubes, and it is one of the cleanest algorithmic ideas in all of computer graphics.

Share this article

Pick a channel — or use your device's native share sheet.

Comments

Loading comments...

https://www.kipuhub.com/en/article/marching-cubes/Content licensed under CC BY-NC 4.0.