Section volumes in the 3D view

/** Outward normal and point of each plane bounding the cell. */
export function cellPlanes(cell) {
  const { min, max } = cell.box, c = centreOf(cell.box), out = [];
  for (const f of FACES) {
    const k = 'xyz'.indexOf(f[0]), hi = f[1] === '1';
    const n0 = [0, 0, 0]; n0[k] = hi ? 1 : -1;
    const p = [...c]; p[k] = hi ? max[k] : min[k];
    const [ta, tb] = cell.tilt?.[f] || [0, 0];
    let n = n0;
    if (ta || tb) {
      const [u, v] = FACE_AXES[f];
      n = [...n0]; n[u] += Math.tan(rad(ta)); n[v] += Math.tan(rad(tb));
      n = unit(n);
    }
    out.push({ n, p, face: f });
  }
  (cell.planes || []).forEach((q, i) => {
    const n = unit([Math.cos(rad(q.el)) * Math.cos(rad(q.az)), Math.cos(rad(q.el)) * Math.sin(rad(q.az)), Math.sin(rad(q.el))]);
    out.push({ n, p: add(c, mul(n, q.off)), extra: i });
  });
  return out;
}

/**
 * The convex polyhedron the planes bound: its corners, each plane's face as a loop of corner indices (counter-clockwise
 * seen from outside; empty when the plane does not touch the cell) and the edges (pairs, once each). Empty cell → no verts.
 */
export function polyhedron(planes) {
  const verts = [], on = [];
  const scale = planes.reduce((m, q) => Math.max(m, Math.abs(q.p[0]), Math.abs(q.p[1]), Math.abs(q.p[2])), 1);
  const tol = scale * 1e-7 + EPS;
  for (let i = 0; i < planes.length; i++) for (let j = i + 1; j < planes.length; j++) for (let k = j + 1; k < planes.length; k++) {
    const x = meet(planes[i], planes[j], planes[k]);
    if (!x || !inside(planes, x, tol * 10)) continue;
    let at = verts.findIndex((v) => Math.hypot(v[0] - x[0], v[1] - x[1], v[2] - x[2]) <= tol * 100);
    if (at < 0) { at = verts.length; verts.push(x); on.push(new Set()); }
    on[at].add(i).add(j).add(k);
  }
  const loops = planes.map((q, i) => {
    const ids = [];
    verts.forEach((v, vi) => { if (on[vi].has(i) || Math.abs(dot(q.n, sub(v, q.p))) <= tol * 100) ids.push(vi); });
    if (ids.length < 3) return [];
    const c = mul(ids.reduce((s, vi) => add(s, verts[vi]), [0, 0, 0]), 1 / ids.length);
    const u = unit(sub(verts[ids[0]], c)), w = cross(q.n, u);
    return ids.map((vi) => { const d = sub(verts[vi], c); return [vi, Math.atan2(dot(d, w), dot(d, u))]; })
      .sort((a, b) => a[1] - b[1]).map((x) => x[0]);
  });
  const seen = new Set(), edges = [];
  for (const L of loops) for (let i = 0; i < L.length; i++) {
    const a = L[i], b = L[(i + 1) % L.length], key = a < b ? a + ',' + b : b + ',' + a;
    if (!seen.has(key)) { seen.add(key); edges.push([a, b]); }
  }
  return { verts, loops, edges };
}

/** A plane through point p with outward normal n, as an extra plane of a cell whose box is `box`. */
export function planeFromPoint(n, p, box) {
  const u = unit(n), c = centreOf(box);
  return { az: ((deg(Math.atan2(u[1], u[0])) % 360) + 360) % 360, el: deg(Math.asin(clamp(u[2], -1, 1))), off: dot(u, sub(p, c)) };
}