File name
Commit message
Commit date
2023-04-18
File name
Commit message
Commit date
File name
Commit message
Commit date
File name
Commit message
Commit date
import Path from "./path.js";
import Polygon from "./polygon.js";
export default class Voronoi {
constructor(delaunay, [xmin, ymin, xmax, ymax] = [0, 0, 960, 500]) {
if (!((xmax = +xmax) >= (xmin = +xmin)) || !((ymax = +ymax) >= (ymin = +ymin))) throw new Error("invalid bounds");
this.delaunay = delaunay;
this._circumcenters = new Float64Array(delaunay.points.length * 2);
this.vectors = new Float64Array(delaunay.points.length * 2);
this.xmax = xmax, this.xmin = xmin;
this.ymax = ymax, this.ymin = ymin;
this._init();
}
update() {
this.delaunay.update();
this._init();
return this;
}
_init() {
const {delaunay: {points, hull, triangles}, vectors} = this;
// Compute circumcenters.
const circumcenters = this.circumcenters = this._circumcenters.subarray(0, triangles.length / 3 * 2);
for (let i = 0, j = 0, n = triangles.length, x, y; i < n; i += 3, j += 2) {
const t1 = triangles[i] * 2;
const t2 = triangles[i + 1] * 2;
const t3 = triangles[i + 2] * 2;
const x1 = points[t1];
const y1 = points[t1 + 1];
const x2 = points[t2];
const y2 = points[t2 + 1];
const x3 = points[t3];
const y3 = points[t3 + 1];
const dx = x2 - x1;
const dy = y2 - y1;
const ex = x3 - x1;
const ey = y3 - y1;
const ab = (dx * ey - dy * ex) * 2;
if (Math.abs(ab) < 1e-9) {
// degenerate case (collinear diagram)
// almost equal points (degenerate triangle)
// the circumcenter is at the infinity, in a
// direction that is:
// 1. orthogonal to the halfedge.
let a = 1e9;
// 2. points away from the center; since the list of triangles starts
// in the center, the first point of the first triangle
// will be our reference
const r = triangles[0] * 2;
a *= Math.sign((points[r] - x1) * ey - (points[r + 1] - y1) * ex);
x = (x1 + x3) / 2 - a * ey;
y = (y1 + y3) / 2 + a * ex;
} else {
const d = 1 / ab;
const bl = dx * dx + dy * dy;
const cl = ex * ex + ey * ey;
x = x1 + (ey * bl - dy * cl) * d;
y = y1 + (dx * cl - ex * bl) * d;
}
circumcenters[j] = x;
circumcenters[j + 1] = y;
}
// Compute exterior cell rays.
let h = hull[hull.length - 1];
let p0, p1 = h * 4;
let x0, x1 = points[2 * h];
let y0, y1 = points[2 * h + 1];
vectors.fill(0);
for (let i = 0; i < hull.length; ++i) {
h = hull[i];
p0 = p1, x0 = x1, y0 = y1;
p1 = h * 4, x1 = points[2 * h], y1 = points[2 * h + 1];
vectors[p0 + 2] = vectors[p1] = y0 - y1;
vectors[p0 + 3] = vectors[p1 + 1] = x1 - x0;
}
}
render(context) {
const buffer = context == null ? context = new Path : undefined;
const {delaunay: {halfedges, inedges, hull}, circumcenters, vectors} = this;
if (hull.length <= 1) return null;
for (let i = 0, n = halfedges.length; i < n; ++i) {
const j = halfedges[i];
if (j < i) continue;
const ti = Math.floor(i / 3) * 2;
const tj = Math.floor(j / 3) * 2;
const xi = circumcenters[ti];
const yi = circumcenters[ti + 1];
const xj = circumcenters[tj];
const yj = circumcenters[tj + 1];
this._renderSegment(xi, yi, xj, yj, context);
}
let h0, h1 = hull[hull.length - 1];
for (let i = 0; i < hull.length; ++i) {
h0 = h1, h1 = hull[i];
const t = Math.floor(inedges[h1] / 3) * 2;
const x = circumcenters[t];
const y = circumcenters[t + 1];
const v = h0 * 4;
const p = this._project(x, y, vectors[v + 2], vectors[v + 3]);
if (p) this._renderSegment(x, y, p[0], p[1], context);
}
return buffer && buffer.value();
}
renderBounds(context) {
const buffer = context == null ? context = new Path : undefined;
context.rect(this.xmin, this.ymin, this.xmax - this.xmin, this.ymax - this.ymin);
return buffer && buffer.value();
}
renderCell(i, context) {
const buffer = context == null ? context = new Path : undefined;
const points = this._clip(i);
if (points === null || !points.length) return;
context.moveTo(points[0], points[1]);
let n = points.length;
while (points[0] === points[n-2] && points[1] === points[n-1] && n > 1) n -= 2;
for (let i = 2; i < n; i += 2) {
if (points[i] !== points[i-2] || points[i+1] !== points[i-1])
context.lineTo(points[i], points[i + 1]);
}
context.closePath();
return buffer && buffer.value();
}
*cellPolygons() {
const {delaunay: {points}} = this;
for (let i = 0, n = points.length / 2; i < n; ++i) {
const cell = this.cellPolygon(i);
if (cell) cell.index = i, yield cell;
}
}
cellPolygon(i) {
const polygon = new Polygon;
this.renderCell(i, polygon);
return polygon.value();
}
_renderSegment(x0, y0, x1, y1, context) {
let S;
const c0 = this._regioncode(x0, y0);
const c1 = this._regioncode(x1, y1);
if (c0 === 0 && c1 === 0) {
context.moveTo(x0, y0);
context.lineTo(x1, y1);
} else if (S = this._clipSegment(x0, y0, x1, y1, c0, c1)) {
context.moveTo(S[0], S[1]);
context.lineTo(S[2], S[3]);
}
}
contains(i, x, y) {
if ((x = +x, x !== x) || (y = +y, y !== y)) return false;
return this.delaunay._step(i, x, y) === i;
}
*neighbors(i) {
const ci = this._clip(i);
if (ci) for (const j of this.delaunay.neighbors(i)) {
const cj = this._clip(j);
// find the common edge
if (cj) loop: for (let ai = 0, li = ci.length; ai < li; ai += 2) {
for (let aj = 0, lj = cj.length; aj < lj; aj += 2) {
if (ci[ai] == cj[aj]
&& ci[ai + 1] == cj[aj + 1]
&& ci[(ai + 2) % li] == cj[(aj + lj - 2) % lj]
&& ci[(ai + 3) % li] == cj[(aj + lj - 1) % lj]
) {
yield j;
break loop;
}
}
}
}
}
_cell(i) {
const {circumcenters, delaunay: {inedges, halfedges, triangles}} = this;
const e0 = inedges[i];
if (e0 === -1) return null; // coincident point
const points = [];
let e = e0;
do {
const t = Math.floor(e / 3);
points.push(circumcenters[t * 2], circumcenters[t * 2 + 1]);
e = e % 3 === 2 ? e - 2 : e + 1;
if (triangles[e] !== i) break; // bad triangulation
e = halfedges[e];
} while (e !== e0 && e !== -1);
return points;
}
_clip(i) {
// degenerate case (1 valid point: return the box)
if (i === 0 && this.delaunay.hull.length === 1) {
return [this.xmax, this.ymin, this.xmax, this.ymax, this.xmin, this.ymax, this.xmin, this.ymin];
}
const points = this._cell(i);
if (points === null) return null;
const {vectors: V} = this;
const v = i * 4;
return V[v] || V[v + 1]
? this._clipInfinite(i, points, V[v], V[v + 1], V[v + 2], V[v + 3])
: this._clipFinite(i, points);
}
_clipFinite(i, points) {
const n = points.length;
let P = null;
let x0, y0, x1 = points[n - 2], y1 = points[n - 1];
let c0, c1 = this._regioncode(x1, y1);
let e0, e1 = 0;
for (let j = 0; j < n; j += 2) {
x0 = x1, y0 = y1, x1 = points[j], y1 = points[j + 1];
c0 = c1, c1 = this._regioncode(x1, y1);
if (c0 === 0 && c1 === 0) {
e0 = e1, e1 = 0;
if (P) P.push(x1, y1);
else P = [x1, y1];
} else {
let S, sx0, sy0, sx1, sy1;
if (c0 === 0) {
if ((S = this._clipSegment(x0, y0, x1, y1, c0, c1)) === null) continue;
[sx0, sy0, sx1, sy1] = S;
} else {
if ((S = this._clipSegment(x1, y1, x0, y0, c1, c0)) === null) continue;
[sx1, sy1, sx0, sy0] = S;
e0 = e1, e1 = this._edgecode(sx0, sy0);
if (e0 && e1) this._edge(i, e0, e1, P, P.length);
if (P) P.push(sx0, sy0);
else P = [sx0, sy0];
}
e0 = e1, e1 = this._edgecode(sx1, sy1);
if (e0 && e1) this._edge(i, e0, e1, P, P.length);
if (P) P.push(sx1, sy1);
else P = [sx1, sy1];
}
}
if (P) {
e0 = e1, e1 = this._edgecode(P[0], P[1]);
if (e0 && e1) this._edge(i, e0, e1, P, P.length);
} else if (this.contains(i, (this.xmin + this.xmax) / 2, (this.ymin + this.ymax) / 2)) {
return [this.xmax, this.ymin, this.xmax, this.ymax, this.xmin, this.ymax, this.xmin, this.ymin];
}
return P;
}
_clipSegment(x0, y0, x1, y1, c0, c1) {
while (true) {
if (c0 === 0 && c1 === 0) return [x0, y0, x1, y1];
if (c0 & c1) return null;
let x, y, c = c0 || c1;
if (c & 0b1000) x = x0 + (x1 - x0) * (this.ymax - y0) / (y1 - y0), y = this.ymax;
else if (c & 0b0100) x = x0 + (x1 - x0) * (this.ymin - y0) / (y1 - y0), y = this.ymin;
else if (c & 0b0010) y = y0 + (y1 - y0) * (this.xmax - x0) / (x1 - x0), x = this.xmax;
else y = y0 + (y1 - y0) * (this.xmin - x0) / (x1 - x0), x = this.xmin;
if (c0) x0 = x, y0 = y, c0 = this._regioncode(x0, y0);
else x1 = x, y1 = y, c1 = this._regioncode(x1, y1);
}
}
_clipInfinite(i, points, vx0, vy0, vxn, vyn) {
let P = Array.from(points), p;
if (p = this._project(P[0], P[1], vx0, vy0)) P.unshift(p[0], p[1]);
if (p = this._project(P[P.length - 2], P[P.length - 1], vxn, vyn)) P.push(p[0], p[1]);
if (P = this._clipFinite(i, P)) {
for (let j = 0, n = P.length, c0, c1 = this._edgecode(P[n - 2], P[n - 1]); j < n; j += 2) {
c0 = c1, c1 = this._edgecode(P[j], P[j + 1]);
if (c0 && c1) j = this._edge(i, c0, c1, P, j), n = P.length;
}
} else if (this.contains(i, (this.xmin + this.xmax) / 2, (this.ymin + this.ymax) / 2)) {
P = [this.xmin, this.ymin, this.xmax, this.ymin, this.xmax, this.ymax, this.xmin, this.ymax];
}
return P;
}
_edge(i, e0, e1, P, j) {
while (e0 !== e1) {
let x, y;
switch (e0) {
case 0b0101: e0 = 0b0100; continue; // top-left
case 0b0100: e0 = 0b0110, x = this.xmax, y = this.ymin; break; // top
case 0b0110: e0 = 0b0010; continue; // top-right
case 0b0010: e0 = 0b1010, x = this.xmax, y = this.ymax; break; // right
case 0b1010: e0 = 0b1000; continue; // bottom-right
case 0b1000: e0 = 0b1001, x = this.xmin, y = this.ymax; break; // bottom
case 0b1001: e0 = 0b0001; continue; // bottom-left
case 0b0001: e0 = 0b0101, x = this.xmin, y = this.ymin; break; // left
}
// Note: this implicitly checks for out of bounds: if P[j] or P[j+1] are
// undefined, the conditional statement will be executed.
if ((P[j] !== x || P[j + 1] !== y) && this.contains(i, x, y)) {
P.splice(j, 0, x, y), j += 2;
}
}
if (P.length > 4) {
for (let i = 0; i < P.length; i+= 2) {
const j = (i + 2) % P.length, k = (i + 4) % P.length;
if (P[i] === P[j] && P[j] === P[k]
|| P[i + 1] === P[j + 1] && P[j + 1] === P[k + 1])
P.splice(j, 2), i -= 2;
}
}
return j;
}
_project(x0, y0, vx, vy) {
let t = Infinity, c, x, y;
if (vy < 0) { // top
if (y0 <= this.ymin) return null;
if ((c = (this.ymin - y0) / vy) < t) y = this.ymin, x = x0 + (t = c) * vx;
} else if (vy > 0) { // bottom
if (y0 >= this.ymax) return null;
if ((c = (this.ymax - y0) / vy) < t) y = this.ymax, x = x0 + (t = c) * vx;
}
if (vx > 0) { // right
if (x0 >= this.xmax) return null;
if ((c = (this.xmax - x0) / vx) < t) x = this.xmax, y = y0 + (t = c) * vy;
} else if (vx < 0) { // left
if (x0 <= this.xmin) return null;
if ((c = (this.xmin - x0) / vx) < t) x = this.xmin, y = y0 + (t = c) * vy;
}
return [x, y];
}
_edgecode(x, y) {
return (x === this.xmin ? 0b0001
: x === this.xmax ? 0b0010 : 0b0000)
| (y === this.ymin ? 0b0100
: y === this.ymax ? 0b1000 : 0b0000);
}
_regioncode(x, y) {
return (x < this.xmin ? 0b0001
: x > this.xmax ? 0b0010 : 0b0000)
| (y < this.ymin ? 0b0100
: y > this.ymax ? 0b1000 : 0b0000);
}
}