30:T4f9,md`## Drawing triangles To *fill* the interior the triangles, we can use the loop from above to get the points for each triangle: \`\`\` for (let triangle_id = 0; triangle_id < delaunator.triangles.length / 3; triangle_id++) { let triangle_points = []; for (int i = 0; i < 3; i++) { let halfedge_id = 3*triangle_id + i; let point_id = delaunator.triangles[halfedge_id]; triangle_points.push(points[point_id]); } // fill the polygon with the three points from triangle_points[] } \`\`\` ## Drawing edges To draw the edges of the triangles instead of filling them, we could use that same loop to draw the half-edges, but it will be simpler to skip the triangles altogether and draw the half-edges directly. \`\`\` for (let halfedge_id = 0; halfedge_id < delaunator.triangles.length; halfedge_id++) { if (halfedge_id < delaunator.halfedges[halfedge_id]) { let p = points[halfedge_id]; let q = points[next_halfedge(halfedge_id)]; // draw edge p to q } } \`\`\` Either way, we will end up drawing edges twice, since each edge is made up of two half-edges. To avoid this, check whether the half-edge id is greater than its opposite (this also handles the case where there is no opposite half-edge). `31:[] 2e:{"title":"Notes on Delaunator","cells":[{"id":0,"value":"md`# Notes on Delaunator`\n","pinned":false,"hidden":false,"mode":"ojs"},{"id":4,"value":"md`[Delaunator](https://github.com/mapbox/delaunator) is a fast library for Delaunay triangulation. The output is not the most intuitive so I'm going to describe the concepts here. This was a **DRAFT**. I used ObservableHQ as a way to quickly prototype but wanted to switch to having version control etc., for some of the same reasons as discussed [here](https://news.ycombinator.com/item?id=17839188). \n\n**The new version is [on the Delaunator site](https://mapbox.github.io/delaunator/)**.`","pinned":false,"hidden":false,"mode":"ojs"},{"id":41,"value":"html`

2D points are the input to Delaunator

`","pinned":false,"hidden":false,"mode":"ojs"},{"id":37,"value":"md`## Triangles\n\nAfter constructing a \\`delaunator = Delaunator.from(points)\\` object, it will have a \\`triangles\\` array, where \\`points[triangles[3*t+i]]\\` will return a point from triangle *t*. To construct the three points for triangle *t*,\n\n\\`\\`\\`\nlet triangle = [];\nfor (int i = 0; i < 3; i++) {\n let edge_id = 3*t + i;\n let point_id = delaunator.triangles[edge_id];\n triangle.push(points[point_id]);\n}\n\\`\\`\\`\n`","pinned":false,"hidden":false,"mode":"ojs"},{"id":56,"value":"html`

Delaunay triangulation

`","pinned":false,"hidden":false,"mode":"ojs"},{"id":543,"value":"md`Note that if you used a different constructor such as \\`new Delaunator()\\` your \\`points\\` array will not be the same as in these examples.`","pinned":false,"hidden":false,"mode":"ojs"},{"id":76,"value":"md`## Half-edges\n\nA triangle edge may be shared with another triangle. Instead of thinking about each edge A↔︎B, we will use two *half-edges* A→B and B→A. Having two half-edges is the key to everything this library provides.\n`","pinned":false,"hidden":false,"mode":"ojs"},{"id":102,"value":"html`

Triangle half-edges

`","pinned":false,"hidden":false,"mode":"ojs"},{"id":480,"value":"md`Triangle ids and half-edge ids are related. \n- The half-edges of triangle *t* are \\`3*t\\`, \\`3*t + 1\\`, and \\`3*t + 2\\`. \n- The triangle of half-edge id *e* is \\`floor(e/3)\\`.\n\nIt will be useful to have some helper functions to go from one half-edge to the next and previous half-edges in the same triangle:`","pinned":false,"hidden":false,"mode":"ojs"},{"id":110,"value":"function next_halfedge(e) { return (e % 3 == 2) ? e-2 : e+1; }","pinned":true,"hidden":false,"mode":"ojs"},{"id":112,"value":"function prev_halfedge(e) { return (e % 3 == 0) ? e+2 : e-1; }","pinned":true,"hidden":false,"mode":"ojs"},{"id":147,"value":"md`Half-edges *e* are the indices into both of delaunator's outputs:\n\n* \\`delaunator.triangles[e]\\` returns the point id where the half-edge starts\n* \\`delaunator.halfedges[e]\\` returns the opposite half-edge in the adjacent triangle, or -1 if there is no adjacent triangle\n`","pinned":false,"hidden":false,"mode":"ojs"},{"id":182,"value":"md`## Triangles adjacent to a triangle\n\nTo find the triangles adjacent to *t*, we can run the same loop as earlier, but instead of indexing into \\`points\\` we'll use the half-edge:\n\n\\`\\`\\`\nlet adjacent_triangles = [];\nfor (int i = 0; i < 3; i++) {\n let halfedge_id = 3*t + i;\n let opposite_halfedge_id = delaunator.halfedges[halfedge_id];\n if (opposite_halfedge_id >= 0) {\n let adjacent_triangle_id = Math.floor(opposite_halfedge_id / 3);\n adjacent_triangles.push(adjacent_triangle_id);\n }\n}\n\\`\\`\\`\n`","pinned":false,"hidden":false,"mode":"ojs"},{"id":534,"value":"html`

Triangles adjacent to a triangle

\n\n

TODO: animation for this, and use a different diagram with just a few triangles

`","pinned":false,"hidden":false,"mode":"ojs"},{"id":155,"value":"$30","pinned":false,"hidden":false,"mode":"ojs"},{"id":177,"value":"md`## Triangles adjacent to a point\n\nTo build Voronoi cells, we need to find the triangles adjacent to a point. This is a bit trickier, but the half-edge structures give us what we need. Let's *start* with a half-edge that leads out of the point. We can alternate two steps to loop around:\n\n1. Use \\`next_halfedge(e)\\` to go to the next half-edge in the *current* triangle\n2. Use \\`halfedges[e]\\` to go to the half-edge in the *adjacent* triangle\n`","pinned":false,"hidden":false,"mode":"ojs"},{"id":243,"value":"html`

Circulate around point

`","pinned":false,"hidden":false,"mode":"ojs"},{"id":310,"value":"filter2 = {\n let step = Math.floor(now/1500);\n return (e) => e == animate2[step % animate2.length];\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":323,"value":"// the outgoing *and* incoming edges as we circulate around a point\nfunction animation_order(delaunator) {\n let animate_e = [];\n let start_e = 3;\n let e = start_e;\n do {\n animate_e.push(e);\n e = delaunator.halfedges[e];\n animate_e.push(e);\n e = next_halfedge(e);\n } while (e !== start_e);\n return animate_e;\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":507,"value":"md`**TODO** : show the code, and ideally use it to calculate the animation\n\nOnce we have the half-edges we can get:\n\n1. The triangles adjacent to the point, using \\`Math.floor(e/3)\\`\n2. The points adjacent to the point, using \\`triangles[opposite_e]\\`\n3. The Voronoi region, using the circumcenter of the triangle\n`","pinned":false,"hidden":false,"mode":"ojs"},{"id":317,"value":"// the animation for the \"circulate around point\" diagram\nanimate2 = animation_order(delaunator2)","pinned":false,"hidden":false,"mode":"ojs"},{"id":202,"value":"md`## TODO Reverse index\n\nTypically we won't actually have a half-edge starting from the point. Delaunator only tells us how to go from half-edges to points. We need to go from points to half-edges. Let's build an index for the \\`triangles\\` array. \n\n`","pinned":false,"hidden":false,"mode":"ojs"},{"id":370,"value":"function starts(delaunator) {\n let index = []; // point id to half-edge id\n for (let e = 0; e < delaunator.triangles.length; e++) {\n index[delaunator.triangles[e]] = e;\n }\n return index;\n}","pinned":true,"hidden":false,"mode":"ojs"},{"id":497,"value":"md`## Convex hull\n\nThere's a problem though with the looping code above. Points on the convex hull won't be completely surrounded by triangles. This will make our loop fail. There are three approaches to this:\n\n1. Ignore it. Make sure never to circulate around points on the convex hull.\n2. Change the code. Handle the convex hull points by adding checks for -1.\n - Change the \\`starts\\` array calculation to point to the \"leftmost\" half-edge for points on the convex hull.\n - Change the loop to stop either at \\`start_e\\` or at -1.\n3. Change the data. Remove the convex hull by wrapping the mesh around the \"back\". There will no longer be any -1 halfedges.\n - Add \"ghost\" half-edges that pair up with the ones that point to -1.\n - Add a single ghost point at \"infinity\" that represents the \"back side\" of the triangulation.\n - Add ghost triangles to connect these ghost half-edges to the ghost point.\n\n**TODO** how much do I want to get into the solutions?\n`","pinned":false,"hidden":false,"mode":"ojs"},{"id":189,"value":"md`## Triangle circumcenters\n\nA Voronoi diagram is built by using the *circumcenters* of each of the Delaunay triangles. The circumcenter is often but not always inside the triangle. For some of my projects, I use centroids, which are always inside the triangle and produce rounder looking cells. See my [comparison page](https://beta.observablehq.com/@redblobgames/centroid-and-voronoi-polygons). The code to traverse the Delaunay triangulation to build polygons will be the same whether using circumcenters or centroids.\n`","pinned":false,"hidden":false,"mode":"ojs"},{"id":193,"value":"function circumcenter(a, b, c) {\n // https://en.wikipedia.org/wiki/Circumscribed_circle#Circumcenter_coordinates\n let ad = a[0]*a[0] + a[1]*a[1],\n bd = b[0]*b[0] + b[1]*b[1],\n cd = c[0]*c[0] + c[1]*c[1];\n let D = 2 * (a[0] * (b[1] - c[1]) + b[0] * (c[1] - a[1]) + c[0] * (a[1] - b[1]));\n return [\n 1/D * (ad * (b[1] - c[1]) + bd * (c[1] - a[1]) + cd * (a[1] - b[1])),\n 1/D * (ad * (c[0] - b[0]) + bd * (a[0] - c[0]) + cd * (b[0] - a[0])),\n ]; \n}\n","pinned":false,"hidden":false,"mode":"ojs"},{"id":198,"value":"html`

Circumcenters of triangles

`","pinned":false,"hidden":false,"mode":"ojs"},{"id":532,"value":"md`## Filling Voronoi cells\n\n**TODO** show that it's the same loop above that finds the half-edges that we can use to link up the circumcenters into the voronoi polygons. Again, we have a problem at the convex hull, and we can either ignore those, or we can clip them (out of scope for this document)\n`","pinned":false,"hidden":false,"mode":"ojs"},{"id":222,"value":"html`

Voronoi cells

`","pinned":false,"hidden":false,"mode":"ojs"},{"id":540,"value":"md`**TODO** show code for traversal for drawing`","pinned":false,"hidden":false,"mode":"ojs"},{"id":549,"value":"md`## Drawing Voronoi cells\n\nIf we only need to draw the edges, we can skip constructing the polygons and draw the edges directly. *Each Delaunay triangle half-edge corresponds to one Voronoi polygon half-edge*. We can iterate over the half-edges and construct the line segments:\n\n\\`\\`\\`\nfor (let halfedge_id = 0; halfedge_id < delaunator.triangles.length; halfedge_id++) {\n if (halfedge_id < delaunator.halfedges[halfedge_id]) {\n let p = t_vertex(points, delaunator, Math.round(halfedge_id/3));\n let q = t_vertex(points, delaunator, Math.round(delaunator.halfedges[halfedge_id]/3));\n // draw edge p to q\n }\n}\n\\`\\`\\`\n**TODO** turn this into a regular function and use it to produce the drawing, maybe with a slider to draw one at a time\n`\n","pinned":false,"hidden":false,"mode":"ojs"},{"id":211,"value":"function t_vertex(points, delaunator, t) {\n let a = points[delaunator.triangles[3*t]],\n b = points[delaunator.triangles[3*t + 1]],\n c = points[delaunator.triangles[3*t + 2]];\n return circumcenter(a, b, c);\n}","pinned":true,"hidden":false,"mode":"ojs"},{"id":94,"value":"md`## Appendix\n\nThe rest of this document has the libraries, code, and data needed to produce the diagrams above.`","pinned":false,"hidden":false,"mode":"ojs"},{"id":59,"value":"md`### Data`","pinned":false,"hidden":false,"mode":"ojs"},{"id":12,"value":"points1 = new Poisson([1000, 1000], 70, undefined, undefined, Math.random).fill()","pinned":false,"hidden":false,"mode":"ojs"},{"id":6,"value":"delaunator1 = Delaunator.from(points1)","pinned":false,"hidden":false,"mode":"ojs"},{"id":294,"value":"points2 = [[320, 170], [220, 270], [400, 270], [530, 50], [100, 80], [300, 30]]","pinned":false,"hidden":false,"mode":"ojs"},{"id":297,"value":"delaunator2 = Delaunator.from(points2)","pinned":false,"hidden":false,"mode":"ojs"},{"id":43,"value":"md`### Visualization`","pinned":false,"hidden":false,"mode":"ojs"},{"id":48,"value":"html`CSS`","pinned":false,"hidden":false,"mode":"ojs"},{"id":46,"value":"function redPointsSvg(points) {\n let results = [``];\n for (let p of points) {\n results.push(``);\n }\n results.push(``);\n return results.join(\"\");\n}\n","pinned":false,"hidden":false,"mode":"ojs"},{"id":195,"value":"function bluePointsSvg(points, delaunator) {\n let results = [``];\n for (let t = 0; t < delaunator.triangles.length / 3; t++) {\n let point = t_vertex(points, delaunator, t);\n results.push(``);\n }\n results.push(``);\n return results.join(\"\");\n}\n ","pinned":false,"hidden":false,"mode":"ojs"},{"id":50,"value":"function triangleSvg(points, delaunator) {\n let results = [``];\n for (let e = 0; e < delaunator.triangles.length; e++) {\n if (e > delaunator.halfedges[e]) {\n let a = points[delaunator.triangles[e]], \n b = points[delaunator.triangles[next_halfedge(e)]];\n results.push(``);\n }\n }\n results.push(``);\n return results.join(\"\");\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":105,"value":"function halfedgeSvg(points, delaunator, filter=()=>true) {\n const dn = 4, dt = 15;\n let results = [``]; \n for (let e1 = 0; e1 < delaunator.halfedges.length; e1++) {\n if (!filter(e1)) { continue; }\n let e2 = next_halfedge(e1);\n let a = points[delaunator.triangles[e1]], \n b = points[delaunator.triangles[e2]],\n n = normal(a, b),\n t = tangent(a, b);\n results.push(``);\n }\n results.push(``);\n return results.join(\"\");\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":207,"value":"function polygonSvg(points, delaunator) {\n let results = [``];\n for (let e = 0; e < delaunator.halfedges.length; e++) {\n if (e < delaunator.halfedges[e]) {\n let a = t_vertex(points, delaunator, Math.floor(e / 3)),\n b = t_vertex(points, delaunator, Math.floor(delaunator.halfedges[e] / 3));\n results.push(``);\n }\n }\n results.push(``);\n return results.join(\"\");\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":132,"value":"html`(svg arrowhead)\n`","pinned":false,"hidden":false,"mode":"ojs"},{"id":92,"value":"md`### Utils`","pinned":false,"hidden":false,"mode":"ojs"},{"id":120,"value":"function tangent(a, b) {\n let dx = b[0] - a[0], dy = b[1] - a[1];\n let d = Math.sqrt(dx*dx + dy*dy);\n return [dx/d, dy/d];\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":114,"value":"function normal(a, b) {\n let dx = b[0] - a[0], dy = b[1] - a[1];\n let d = Math.sqrt(dx*dx + dy*dy);\n return [dy/d, -dx/d];\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":8,"value":"md`### Libraries\n\nPoints for the diagrams chosen using the [poisson-disk-sampling](https://github.com/kchapelier/poisson-disk-sampling) library. The Delaunator library builds the triangles.`","pinned":false,"hidden":false,"mode":"ojs"},{"id":10,"value":"Poisson = require('https://bundle.run/poisson-disk-sampling@1.0.2')","pinned":false,"hidden":false,"mode":"ojs"},{"id":30,"value":"Delaunator = require(\"delaunator@2.0.0\")","pinned":false,"hidden":false,"mode":"ojs"}],"files":"$W31","stdlib":"1","titling":"auto"}