[browser] Bezier(3-degree Spline) Interpolation 2D

bezier-interpolation.js

// make equation coefficients and vectors of x and y for bezier interpolation
export function bezierControlsEquation(knots) {
  // [Theory of Beizer Interpolation]
  // bezier curve:
  // - p(t) = (1-t)^3*P0 + 3(1-t)^2*t*C0 + 3(1-t)t^2*C1 + t^3*P1
  //     - p(0) = P0
  //     - p(1) = P1
  // - p'(t) = -3(t^2-2t+1)*P0 + 3(3t^2-4t+1)*C0 + 3(-3t^2+2t)*C1 + 3t^2*P1
  //     - p'(0) = -3*P0 + 3*C0
  //     - p'(1) = -3*C1 + 3*P1
  // - p''(t) = -3(2t-2)*P0 + 3(6t-4)*C0 + 3 (-6t+2)*C1 + 6t*P1
  //     - p''(0) = 6*P0 - 12*C0 + 6*C1
  //     - p''(1) = 6*C0 - 12*C1 + 6*P1
  //
  // relation of knots and control-points:
  // - k[0]--c[0]--c[1]--k[1]-- ... --k[n-2]--c[2n-4]--c[2n-3]--k[n-1]
  //      -  ... --c[2i-2]--c[2i-1]--k[i]--c[2i]--c[2i+1]-- ...
  //
  // point relation of lines:
  // - p[i  ](t): k[i  ]--c[2i  ]--c[2i+1]--k[i+1]
  // - p[i+1](t): k[i+1]--c[2i+2]--c[2i+3]--k[i+2]
  //
  // conditions of bezier (as 3-degree spline) interpolation (i = 0 .. n-2):
  // - A. p''[0](0) = 0
  // - B. p''[i](1) = p''[i+1](0)
  // - C. p'[i](1) = p'[i+1](0)
  // - D. p''[n-1](1) = 0
  //
  // normalize A:
  // - 6*k[0] - 12*c[0] + 6*c[1] = 0
  // - 2*c[0] - 1*c[1] = k[0]
  //
  // normalize B:
  // - 6*c[2i] - 12*c[2i+1] + 6*k[i+1] = 6*k[i+1] - 12*c[2i+2] + 6*c[2i+3]
  // - 1*c[2i] - 2*c[2i+1] + 2*c[2i+2] - 1*c[2i+3] = 0
  //
  // normalize C:
  // - -3*c[2i+1] + 3*k[i+1] = -3*k[i+1] + 3*c[2i+2]
  // - 1*c[2i+1] + 1*c[2i+2] = 2*k[i+1]
  //
  // normalize D:
  // - 6*c[2n-4] - 12+c[2n-3] + 6*k[n-1] = 0
  // - -1*c[2n-4] + 2*c[2n-3] = k[n-1]
  //
  // normalized conditions of beziere interpolation:
  // - A. 2*c[0]  - 1*c[1]                            = k[0]
  // - B. 1*c[2i] - 2*c[2i+1] + 2*c[2i+2] - 1*c[2i+3] = 0
  // - C.           1*c[2i+1] + 1*c[2i+2]             = 2*k[i+1]
  // - D.                      -1*c[2n-4] + 2*c[2n-3] = k[n-1]
  //
  const n = knots.length;
  const cn = 2 * n - 2;
  const cA = [[2, -1].concat(Array(cn - 2).fill(0))];
  const cD = [Array(cn - 2).fill(0).concat([-1, 2])];
  const cBC = [...Array(n - 2)].flatMap((_, i) => {
    const cB = Array(cn).fill(0), cC = Array(cn).fill(0);
    [cB[2 * i], cB[2 * i + 1], cB[2 * i + 2], cB[2 * i + 3]] = [1, -2, 2, -1];
    [cC[2 * i + 1], cC[2 * i + 2]] = [1, 1];
    return [cB, cC];
  });
  const coeffs = cA.concat(cBC, cD);

  const vList = [...Array(knots[0].length)].map((_, i) => {
    const vA = [knots[0][i]], vD = [knots[n - 1][i]];
    const vBC = [...Array(n - 2)].flatMap((_, j) => [0, 2 * knots[j + 1][i]]);
    return vA.concat(vBC, vD);
  });
  return [coeffs, vList];
}

export function bezierLoopControlsEquation(knots) {
  // Bezier Interpolation as Loop
  // Loop knots:
  // - k[n] = k[0]
  // - c[2n] = c[0]
  // - p[n-1](t) = (1-t)^3k[n-1] + 3(1-t)^2tc[2n-2] + 3(1-t)t^2c[2n-1] +t^3k[0]
  //
  // normalized conditions of beziere interpolation:
  // - B. 1*c[2i] - 2*c[2i+1] + 2*c[2i+2] - 1*c[2i+3] = 0
  // - C.           1*c[2i+1] + 1*c[2i+2]             = 2*k[i+1]
  //
  // Edge case of the last i=n-1:
  // - B0. 1*c[2n-2] - 2*c[2n-1] + 2*c[0] - 1*c[1] = 0
  // - C0.             1*c[2i-1] + 1*c[0]          = 2*k[0]
  //
  const n = knots.length;
  const cn = 2 * n;
  const cBC = [...Array(n - 1)].flatMap((_, i) => {
    const cB = Array(cn).fill(0), cC = Array(cn).fill(0);
    [cB[2 * i], cB[2 * i + 1], cB[2 * i + 2], cB[2 * i + 3]] = [1, -2, 2, -1];
    [cC[2 * i + 1], cC[2 * i + 2]] = [1, 1];
    return [cB, cC];
  });
  const cB = Array(cn).fill(0), cC = Array(cn).fill(0);
  [cB[2 * n - 2], cB[2 * n - 1], cB[0], cB[1]] = [1, -2, 2, -1];
  [cC[2 * n - 1], cC[0]] = [1, 1];
  const coeffs = cBC.concat([cB, cC]);

  const vList = [...Array(knots[0].length)].map((_, i) => {
    const vBC = [...Array(n - 1)].flatMap((_, j) => [0, 2 * knots[j + 1][i]]);
    return vBC.concat([0, 2 * knots[0][i]]);
  });
  return [coeffs, vList];
}


export function solve(coeffs, vList) {
  // gaussian elimination
  console.assert(vList.every(v => v.length === coeffs.length));
  console.assert(coeffs.every(c => c.length === coeffs.length));
  //console.assert(coeffs.every((c, i) => c[i] !== 0));
  const cs = coeffs.map(l => l.slice()), vs = vList.map(v => v.slice());
  for (let i = 0; i < cs.length; i++) {
    for (let j = i + 1; j < cs.length; j++) {
      if (cs[i][i] === 0) continue;
      const f = cs[j][i] / cs[i][i];
      cs[j][i] = 0;
      for (let k = i + 1; k < cs[j].length; k++) cs[j][k] -= f * cs[i][k];
      for (const v of vs) v[j] -= f * v[i];
    }
  }
  for (let i = cs.length - 1; i >= 0; i--) {
    const d = cs[i][i];
    cs[i][i] = 1;
    for (let j = i - 1; j >= 0; j--) cs[i][j] /= d;
    for (const v of vs) v[i] /= d;
    for (let j = i - 1; j >= 0; j--) {
      const f = cs[j][i];
      cs[j][i] = 0;
      for (let k = i - 1; k >= 0; k--) cs[j][k] -= f * cs[i][k];
      for (const v of vs) v[j] -= f * v[i];
    }
  }
  return vs;
}

index.html

<!doctype html>
<html>
  <head>
    <meta charset="utf-8">
    <link rel="icon" href="data:image/x-icon;">
    <script type="module" src="./main.js"></script>
  </head>
  <body>
    <h1>Bezier Interpolation</h1>
    <ul><li>Click to put several knots for bezier interpolatation</li></ul>
  </body>
</html>

main.js

import {bezierControlsEquation, bezierLoopControlsEquation, solve}
from "./bezier-interpolation.js";

{
  const canvas = document.createElement("canvas");
  canvas.width = canvas.height = 512;
  canvas.style = "border: solid; display: block;";
  const pop = document.createElement("button");
  pop.textContent = "remove the last knot";
  const clear = document.createElement("button");
  clear.textContent = "clear";
  const showControls = document.createElement("input");
  showControls.type = "checkbox";
  showControls.checked = true;
  const showControlsLabel = document.createElement("label");
  showControlsLabel.append(showControls, "show control points");;
  const loop = document.createElement("input");
  loop.type = "checkbox";
  loop.checked = true;
  const loopLabel = document.createElement("label");
  loopLabel.append(loop, "loop");;
  document.body.append(canvas, loopLabel, showControlsLabel, pop, clear);

  const c2d = canvas.getContext("2d");
  const knots = [];
  const redraw = () => paint(c2d, knots, showControls.checked, loop.checked);
  canvas.addEventListener("mouseup", ev => {
    const rect = ev.currentTarget.getBoundingClientRect();
    const x = ev.clientX - rect.left, y = ev.clientY - rect.top;
    knots.push([x, y]);
    redraw();
  });
  pop.addEventListener("click", ev => {
    knots.pop();
    redraw();
  });
  clear.addEventListener("click", ev => {
    knots.splice(0, knots.length);
    redraw();
  });
  showControls.addEventListener("click", redraw);
  loop.addEventListener("click", redraw);
}

function paint(c2d, knots, showControls = true, loop = false) {
  c2d.clearRect(0, 0, c2d.canvas.width, c2d.canvas.height);
  paintKnots(c2d, knots);
  if (knots.length < 2) return;
  const [coeffs, vList] = loop ? bezierLoopControlsEquation(knots) :
        bezierControlsEquation(knots);
  //console.log(coeffs, vList);
  const [cx, cy] = solve(coeffs, vList);
  //console.log(cx, cy);
  const k = loop ? knots.concat([knots[0]]) : knots;
  if (showControls) paintControls(c2d, k, cx, cy);
  paintCurves(c2d, k, cx, cy);
}

function paintControls(c2d, knots, cx, cy) {
  c2d.strokeStyle = "blue";
  c2d.beginPath();
  for (let i = 0; i < knots.length - 1; i++) {
    const c1 = 2 * i, c2 = 2 * i + 1;
    const x1 = knots[i][0], y1 = knots[i][1];
    const x2 = knots[i + 1][0], y2 = knots[i + 1][1];
    c2d.moveTo(x1, y1);
    c2d.lineTo(cx[c1], cy[c1]);
    c2d.moveTo(cx[c2], cy[c2]);
    c2d.lineTo(x2, y2);
  }
  c2d.moveTo(cx[cx.length - 1], cy[cy.length - 1]);
  c2d.lineTo(knots[knots.length - 1][0], knots[knots.length - 1][1]);
  c2d.strokeStyle = "blue";
  c2d.stroke();

  for (let i = 0; i < cx.length; i++) {
    c2d.beginPath();
    c2d.arc(cx[i], cy[i], 3, 0, 2 * Math.PI, true);
    c2d.closePath();
    c2d.fillStyle = "blue";
    c2d.fill();
  }
}

function paintCurves(c2d, knots, cx, cy) {
  c2d.beginPath();
  c2d.moveTo(knots[0][0], knots[0][1]);
  for (let i = 1; i < knots.length; i++) {
    const c1 = 2 * i - 2, c2 = 2 * i - 1, x = knots[i][0], y = knots[i][1];
    c2d.bezierCurveTo(cx[c1], cy[c1], cx[c2], cy[c2], x, y);
  }
  c2d.strokeStyle = "black";
  c2d.stroke();
}

function paintKnots(c2d, knots) {
  for (const [x, y] of knots) {
    c2d.beginPath();
    c2d.arc(x, y, 2, 0, 2 * Math.PI, true);
    c2d.closePath();
    c2d.fillStyle = "black";
    c2d.fill();
  }
}

demo: https://gist.githack.com/bellbind/af27993852c8c5e33e450bb9f5c75f8c/raw/index.html