Resources Numerical accuracy Calculation policy

Terrain from contours and terrain grids

Builds a ground surface from contour lines, or a flat grid to sculpt.

How the precision is kept

The triangulation is Bowyer–Watson in plan with 64-bit numbers. The in-circle decision compares a squared distance with a circumcircle radius computed in floating point, not with an exact predicate. Heights are carried through unchanged; points are merged on a 1 cm grid in plan before triangulating.

Open item. Subtract the box center from the points before the circumcircles are computed, so the squares of large coordinates do not use up the digits.

Calculations under this policy

3 calculations in 1 source file. The code is read from the released app when this page loads. It is the calculation people run, not a copy kept for this page.

terrain.js ModelCraft 0.8.1 64-bit 3 calculations Reviewed Nov 6, 2026

Tests bin/test-modelcraft.mjs

Delaunay triangulation

64-bit

Triangles over a set of points in plan.

D = 2 (aₓ(b_y − c_y) + bₓ(c_y − a_y) + cₓ(a_y − b_y)); center = (Σ |p|² × Δ) ÷ D; inside when |p − center|² < r² (1 − 10⁻¹²)

The circumcenter uses the squares of the absolute coordinates, so for a site far from the origin the products are large and their differences lose digits. The decision is a floating-point comparison with a relative margin of 1e-12, not an exact predicate. A finished triangle whose cross product is below 1e-14 (absolute) is dropped.

terrain.js
lines 15 to 53
ModelCraft 0.8.1

export function delaunay(pts) {
  const n = pts.length;
  if (n < 3) return [];
  let x0 = Infinity, y0 = Infinity, x1 = -Infinity, y1 = -Infinity;
  for (const [x, y] of pts) { x0 = Math.min(x0, x); y0 = Math.min(y0, y); x1 = Math.max(x1, x); y1 = Math.max(y1, y); }
  const d = Math.max(x1 - x0, y1 - y0) || 1, cx = (x0 + x1) / 2, cy = (y0 + y1) / 2;
  const P = pts.concat([[cx - 20 * d, cy - d], [cx, cy + 20 * d], [cx + 20 * d, cy - d]]);
  const circ = (a, b, c) => {
    const [ax, ay] = P[a], [bx, by] = P[b], [qx, qy] = P[c];
    const D = 2 * (ax * (by - qy) + bx * (qy - ay) + qx * (ay - by));
    if (Math.abs(D) < 1e-18) return { x: 0, y: 0, r2: Infinity };
    const a2 = ax * ax + ay * ay, b2 = bx * bx + by * by, c2 = qx * qx + qy * qy;
    const x = (a2 * (by - qy) + b2 * (qy - ay) + c2 * (ay - by)) / D, y = (a2 * (qx - bx) + b2 * (ax - qx) + c2 * (bx - ax)) / D;
    return { x, y, r2: (ax - x) ** 2 + (ay - y) ** 2 };
  };
  let tris = [{ v: [n, n + 1, n + 2], c: circ(n, n + 1, n + 2) }];
  for (let i = 0; i < n; i++) {
    const [px, py] = P[i];
    const bad = [], keep = [];
    for (const tr of tris) ((px - tr.c.x) ** 2 + (py - tr.c.y) ** 2 < tr.c.r2 * (1 - 1e-12) ? bad : keep).push(tr);
    // the hole's rim: edges of the bad triangles that only one of them has
    const count = new Map();
    for (const tr of bad) for (let k = 0; k < 3; k++) { const a = tr.v[k], b = tr.v[(k + 1) % 3], key = a < b ? a + ',' + b : b + ',' + a; count.set(key, (count.get(key) || 0) + 1); }
    tris = keep;
    for (const tr of bad) for (let k = 0; k < 3; k++) {
      const a = tr.v[k], b = tr.v[(k + 1) % 3];
      if (count.get(a < b ? a + ',' + b : b + ',' + a) === 1) tris.push({ v: [a, b, i], c: circ(a, b, i) });
    }
  }
  const out = [];
  for (const { v } of tris) {
    if (v.some((k) => k >= n)) continue;
    const [a, b, c] = v, [ax, ay] = P[a], [bx, by] = P[b], [qx, qy] = P[c];
    const cr = (bx - ax) * (qy - ay) - (by - ay) * (qx - ax);
    if (Math.abs(cr) < 1e-14) continue;
    out.push(cr > 0 ? [a, b, c] : [a, c, b]);
  }
  return out;
}

Terrain from contour edges

64-bit

Collects the contour points, merges close ones and triangulates them.

key = round(100 x), round(100 y); thinning step = n ÷ 8000

Points are merged by rounding x and y to a 0.01 grid; the first point seen in a cell is kept with its own unrounded coordinates and height. Two points closer than 0.01 on either side of a cell edge are not merged. Above 8000 points every (n ÷ 8000)-th point is kept.

terrain.js
lines 61 to 78
ModelCraft 0.8.1

export function terrainFromContours(model, def, eids, { name = 'Terrain' } = {}) {   // i18n-ignore: a fallback name
  const seen = new Map();
  for (const eid of eids) {
    const e = def.edges.get(eid);
    if (!e) continue;
    for (const v of [e.a, e.b]) {
      const p = def.verts.get(v), k = Math.round(p[0] * 100) + ',' + Math.round(p[1] * 100);
      if (!seen.has(k)) seen.set(k, p);
    }
  }
  let pts = [...seen.values()];
  if (pts.length > 8000) { const step = pts.length / 8000; pts = Array.from({ length: 8000 }, (_, i) => pts[Math.floor(i * step)]); }
  const tris = delaunay(pts.map((p) => [p[0], p[1]]));
  if (!tris.length) return null;
  const d = meshDef(model, tris.map((tr) => tr.map((i) => pts[i].slice())), model.uniqueName(name), 'group');
  softInside(d);
  return d;
}

Terrain grid

64-bit

A flat grid of triangles of a given width, depth and spacing.

nₓ = round(width ÷ spacing), 1…200; x = x₀ + width × i ÷ nₓ

Each grid point is computed from its index by one multiplication and one division, so points do not drift along a row. The cell size is width ÷ nₓ, which equals the requested spacing only when the width is a whole number of spacings.

terrain.js
lines 80 to 93
ModelCraft 0.8.1

export function terrainGrid(model, { width = 20, depth = 20, spacing = 1, origin = [0, 0, 0], name = 'Terrain' } = {}) {   // i18n-ignore: a fallback name
  const nx = Math.max(1, Math.min(200, Math.round(width / spacing))), ny = Math.max(1, Math.min(200, Math.round(depth / spacing)));
  const P = (i, j) => [origin[0] + (width * i) / nx, origin[1] + (depth * j) / ny, origin[2]];
  const polys = [];
  for (let i = 0; i < nx; i++) for (let j = 0; j < ny; j++) {
    // alternate the diagonals, so slopes sculpted later have no grain
    const a = P(i, j), b = P(i + 1, j), c = P(i + 1, j + 1), e = P(i, j + 1);
    if ((i + j) % 2) { polys.push([a, b, c], [a, c, e]); } else { polys.push([a, b, e], [b, c, e]); }
  }
  const d = meshDef(model, polys, model.uniqueName(name), 'group');
  // meshDef merges nothing; flat triangles stay separate faces (each one moves on its own when sculpted)
  softInside(d);
  return d;
}

In detail

Method: Bowyer–Watson Delaunay triangulation over the contour points in plan, inside a super-triangle 20 times the extent; grid points by division of the width and depth.
Number types: 64-bit numbers.
Tolerances: In-circle: squared distance below r² × (1 − 1e-12), relative. Circumcircle denominator below 1e-18 and triangle cross product below 1e-14 are treated as degenerate, both absolute. Points merge on a grid of 0.01 in plan.
Limits: At most 8000 points; above that the points are thinned evenly. The circumcenter is computed from the points' absolute coordinates.

Approach: 64-bit, direct

Every value is a 64-bit number from input to result. Formulas are closed form or a direct solve, and every comparison uses a stated tolerance.

Every policy that takes this approach