2f:T10ae,As explained in our [previous notebook](https://observablehq.com/@allmaps/a-simple-polygon-triangulation-algorithm), we need a triangulation algorithm in our Allmaps render pipeline, and it's hard to find one that suits our needs. In the previous notebook we created one ourself building on [poly2tri.js](https://github.com/r3mi/poly2tri.js). But alas, that library hasn't been updated for 7 years and uses undesirable approaches we can't cope with any longer: it's not published as an ECMAScript Module (ESM) and makes impractical use of the `global` object. This was creating too much hassle for us, so we needed to look at a different approach. Instead of patching the error in the poly2tri.js library and re-publishing it (like others have apparently done: [patch](https://github.com/jarvis394/poly2tri.js/commit/d96a6f2a33b51cc3c6bfe09cc464a473343986e6), [package](https://www.npmjs.com/package/re-poly2tri)), we decided to put together a now polygon triangulation approach based on three other modern libraries which are all published as ESM in TypeScript or with types available. Our new approach is the following: - [**Delaunator**](https://github.com/mapbox/delaunator) is given the points of our input polygon and returns a Delaunay triangulation from them. This is a good first step, but not exactly what we need yet: we are not guaranteed that the edges of out input polygon are still present in the Delaunay triangulation based only on the points of out input polygon, and furthermore if our polygon is concave it's Delaunay triangulation will contain triangles outside of it. (You can see examples of both below). We solve this using two other packages. - [**Constrainautor**](https://github.com/kninnug/Constrainautor) takes in the output of Delaunator plus a set of edges to respect. We pass all the edges of our polygon as edges to respect. It then 'flips' the edges of the Delaunay triangulation where needed to make sure our polygon edges are present in the (now constrained) Delaunay triangulation. See [this example of a flip](https://observablehq.com/@fil/hello-constrainautor) and this [practical application on some interesting polygons](https://observablehq.com/@fil/polygon-triangulation). - [**Robust Point In Polygon**](https://github.com/mikolalysenko/robust-point-in-polygon) can then check for all resulting triangles if their midpoint is inside the polygon. We need one more thing to satisfy our original goal: we wanted an algorithm that would triangulate an input polygon and automatically produce **smaller, well conditioned triangles** in larger polygons (this, by the way, is why [Earcut](https://github.com/mapbox/earcut) was never an option: it produces long, pointy trianges). Previously we achieved this by creating a grid and adding grid points as 'Steiner points' to poly2tri.js. We can use a similar approach here: we again create a grid within the bbox of our polygon and add these grid points to the Delaunator input (if they are inside our polygon) and also again split up ('interpolate') the polygon edges in smaller pieces when they are longer then the grid distance. An extra advantage of this new approach is that it is more accepting of weird polygon shapes, for which Poly2tri.js sometimes threw errors. Possible improvements (which will be made in the new version of our [Allmaps Triangulate](https://github.com/allmaps/allmaps/tree/main/packages/triangulate) package but are outside of the scope of this notebook): speed up the 'triangle in polygon' check by only checking the triangles with at least one point on the edge, and expand the approach to support holes. Our final approach checks off all our boxes: - ✔️ It's based on al algorithm that implements constrained Delaunay triangulation of points set. We create smaller, well conditioned triangles using extra points from a grid. - ✔️ Written in JavaScript (or TypeScript) and on NPM, **published as an ESM module**. - ✔️ Small: total package size (minified + gziped): 3.1kB + 5.1kB + 2kB - ✔️ Robust enough: better than Poly2tri.js! - ✔️ Fast enough: our input polygons are small and few, so we don't need record speeds, but Robust Point In Polygon and especially Delaunator are known for their correctness and impressive speeds!30:[] 2d:{"title":"A modern polygon triangulation algorithm","cells":[{"id":0,"value":"# A modern polygon triangulation algorithm","pinned":false,"hidden":false,"mode":"md"},{"id":873,"value":"$2f","pinned":false,"hidden":false,"mode":"md"},{"id":884,"value":"## Try it out!","pinned":false,"hidden":false,"mode":"md"},{"id":854,"value":"# ✎👇","pinned":false,"hidden":false,"mode":"md"},{"id":683,"value":"viewof polygon = polygonInput({\n title: \"Draw a polygon\",\n description:\n \"Drawing three or more points creates a polygon. It's triangulation is then shown below.\"\n})","pinned":false,"hidden":false,"mode":"ojs"},{"id":13,"value":"{\n let width = 400;\n let height = Math.round((210 / 400) * width);\n\n let ctx = DOM.context2d(width, height);\n\n function render() {\n ctx.clearRect(0, 0, width, width);\n if (useGrid) {\n drawGrid(ctx, grid);\n }\n if (polygon.length > 0) {\n drawPolygon(ctx, polygon);\n }\n if (polygon.length > 2 && approach == \"Poly2tri\") {\n drawPoly2triTriangles(ctx, poly2triTriangles);\n }\n if (polygon.length > 2 && approach == \"Delaunator\") {\n drawDelaunatorTriangles(ctx, delaunatorTriangles);\n }\n if (polygon.length > 2 && approach == \"Constrainautor\") {\n drawDelaunatorTriangles(ctx, constrainautorTriangles);\n }\n if (polygon.length > 2 && approach == \"Constrainautor + inside\") {\n drawDelaunatorTriangles(ctx, insideTriangles);\n }\n }\n\n render();\n\n return ctx.canvas;\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":800,"value":"Play around with these settings to better understand what's happening","pinned":false,"hidden":false,"mode":"md"},{"id":1404,"value":"viewof approach = Inputs.select(\n [\"Poly2tri\", \"Delaunator\", \"Constrainautor\", \"Constrainautor + inside\"],\n { label: \"Select triangulation approach\", value: \"Constrainautor + inside\" }\n)","pinned":false,"hidden":false,"mode":"ojs"},{"id":803,"value":"viewof useGrid = Inputs.toggle({ label: \"Use Grid\", value: true })","pinned":false,"hidden":false,"mode":"ojs"},{"id":5,"value":"viewof gridSize = Inputs.range([5, 80], {\n value: 20,\n step: 1,\n label: \"Grid size\"\n})","pinned":false,"hidden":false,"mode":"ojs"},{"id":805,"value":"viewof usePolygonInterpolation = Inputs.toggle({\n label: \"Use Polygon Interpolation\",\n value: true\n})","pinned":false,"hidden":false,"mode":"ojs"},{"id":818,"value":"---","pinned":false,"hidden":false,"mode":"md"},{"id":666,"value":"### Polygon","pinned":false,"hidden":false,"mode":"md"},{"id":679,"value":"import { polygonInput } from \"@harrystevens/polygon-input\"","pinned":false,"hidden":false,"mode":"ojs"},{"id":61,"value":"drawPolygon = function (ctx, polygon, f) {\n ctx.fillStyle = \"white\";\n for (let i = 0; i + 1 < polygon.length; i++) {\n ctx.beginPath();\n ctx.arc(polygon[i][0], polygon[i][1], 4, 0, 2 * Math.PI, false);\n ctx.moveTo(polygon[i][0], polygon[i][1]);\n ctx.lineTo(polygon[i + 1][0], polygon[i + 1][1]);\n ctx.stroke(); // Render the path\n ctx.fill();\n }\n ctx.beginPath();\n var i = polygon.length - 1;\n ctx.arc(polygon[i][0], polygon[i][1], 4, 0, 2 * Math.PI, false);\n ctx.moveTo(polygon[i][0], polygon[i][1]);\n ctx.lineTo(polygon[0][0], polygon[0][1]);\n ctx.stroke(); // Render the path\n ctx.fill();\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":671,"value":"### Grid","pinned":false,"hidden":false,"mode":"md"},{"id":773,"value":"getBbox = function (polygon) {\n // polygon = [[3,5], [7,10], [20, 1]]\n const result = [Infinity, Infinity, -Infinity, -Infinity];\n polygon.forEach((point) => {\n if (result[0] > point[0]) {\n result[0] = point[0];\n }\n if (result[1] > point[1]) {\n result[1] = point[1];\n }\n if (result[2] < point[0]) {\n result[2] = point[0];\n }\n if (result[3] < point[1]) {\n result[3] = point[1];\n }\n });\n return result;\n // result = [minX, minY, maxX, maxY]\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":9,"value":"grid = {\n let grid = [];\n let bbox = getBbox(polygon);\n for (let x = bbox[0] + gridSize, i = 0; x <= bbox[2]; i++, x += gridSize) {\n for (let y = bbox[1] + gridSize, j = 0; y <= bbox[3]; j++, y += gridSize) {\n grid.push([x, y]);\n }\n }\n return grid;\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":11,"value":"drawGrid = function (ctx, grid, f) {\n for (let i = 0; i < grid.length; i++) {\n ctx.fillStyle = \"grey\";\n ctx.beginPath();\n ctx.arc(grid[i][0], grid[i][1], 2, 0, 2 * Math.PI, false);\n ctx.fill();\n }\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":584,"value":"### Triangulation imports","pinned":false,"hidden":false,"mode":"md"},{"id":83,"value":"poly2tri = import('https://cdn.skypack.dev/poly2tri@1.5.0?min')","pinned":false,"hidden":false,"mode":"ojs"},{"id":932,"value":"robustPointInPolygon = import(\n \"https://cdn.skypack.dev/robust-point-in-polygon@1.0.3\"\n)","pinned":false,"hidden":false,"mode":"ojs"},{"id":1073,"value":"Delaunator = require(\"delaunator@5.0.0\")","pinned":false,"hidden":false,"mode":"ojs"},{"id":1071,"value":"Constrainautor = require(\"@kninnug/constrainautor\")","pinned":false,"hidden":false,"mode":"ojs"},{"id":576,"value":"### Triangulation helper functions","pinned":false,"hidden":false,"mode":"md"},{"id":324,"value":"getSegmentDistance = function (segment) {\n // segment = [[x0, y0], [x1, y1]]\n // return lenght of segment\n return Math.sqrt(\n Math.pow(segment[1][0] - segment[0][0], 2) +\n Math.pow(segment[1][1] - segment[0][1], 2)\n );\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":339,"value":"getSegmentAngle = function (segment) {\n // segment = [[x0, y0], [x1, y1]]\n // return angle of segment (in radians, signed)\n return Math.atan2(\n segment[1][1] - segment[0][1],\n segment[1][0] - segment[0][0]\n );\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":329,"value":"getNextPoint = function (point, distance, angle) {\n // point = [x0, y0], distance = number, angle = number in radians, signed\n // compute next point from point along angle over distance\n return [\n point[0] + Math.cos(angle) * distance,\n point[1] + Math.sin(angle) * distance\n ];\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":298,"value":"interpolateSegment = function (segment, distance) {\n // segment = [[3,5], [7,10]]\n // distance = 1\n\n var currentPoint = segment[0];\n var result = [currentPoint];\n\n while (getSegmentDistance([currentPoint, segment[1]]) > distance) {\n var nextPoint = getNextPoint(\n currentPoint,\n distance,\n getSegmentAngle(segment)\n );\n result.push(nextPoint);\n currentPoint = nextPoint;\n }\n // segment[1] is not pushed\n return result;\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":516,"value":"interpolatePolygon = function (polygon, distance) {\n // polygon = [[3,5], [7,10], [20, 1]]\n // distance = 1\n\n // close polygon\n polygon = [...polygon, polygon[0]];\n\n var result = [];\n for (let i = 0; i < polygon.length - 1; i++) {\n result = result.concat(\n interpolateSegment([polygon[i], polygon[i + 1]], distance)\n );\n }\n\n return result;\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":115,"value":"drawPoly2triTriangles = function (ctx, triangles) {\n triangles.forEach(function (t) {\n ctx.fillStyle = \"hsl(\" + 360 * Math.random() + \", 50%, 50%, 50%)\";\n ctx.beginPath();\n ctx.moveTo(t.getPoint(0).x, t.getPoint(0).y);\n ctx.lineTo(t.getPoint(1).x, t.getPoint(1).y);\n ctx.lineTo(t.getPoint(2).x, t.getPoint(2).y);\n ctx.fill();\n });\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":1089,"value":"drawDelaunatorTriangles = function (ctx, triangles) {\n triangles.forEach(function (t) {\n ctx.fillStyle = \"hsl(\" + 360 * Math.random() + \", 50%, 50%, 50%)\";\n ctx.beginPath();\n ctx.moveTo(t[0][0], t[0][1]);\n ctx.lineTo(t[1][0], t[1][1]);\n ctx.lineTo(t[2][0], t[2][1]);\n ctx.fill();\n });\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":607,"value":"### Poly2tri triangulation","pinned":false,"hidden":false,"mode":"md"},{"id":488,"value":"TRIANG = {\n function TRIANG(polygon, points, distance) {\n // initialise Constrained Delaunay Triangulation with polygon\n let polygonUsed = usePolygonInterpolation\n ? interpolatePolygon(polygon, distance)\n : polygon;\n var swctx = new poly2tri.SweepContext(\n polygonUsed.map((p) => new poly2tri.Point(p[0], p[1]))\n );\n if (useGrid) {\n // add grid points as Steiner points if they are inside the polygon\n for (var i = 0; i < points.length; i++) {\n if (robustPointInPolygon.default(polygon, points[i]) == -1) {\n swctx.addPoint(new poly2tri.Point(points[i][0], points[i][1]));\n }\n }\n }\n // triangulate\n swctx.triangulate();\n let triangles = swctx.getTriangles();\n return triangles;\n }\n return TRIANG;\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":478,"value":"poly2triTriangles = polygon.length > 2 ? TRIANG(polygon, grid, gridSize) : []","pinned":false,"hidden":false,"mode":"ojs"},{"id":1443,"value":"### New, modern triangulation","pinned":false,"hidden":false,"mode":"md"},{"id":1334,"value":"usedPolygon = usePolygonInterpolation\n ? interpolatePolygon(polygon, gridSize)\n : polygon","pinned":true,"hidden":false,"mode":"ojs"},{"id":1373,"value":"gridPointsInPolygon = grid.filter((point) => {\n if (robustPointInPolygon.default(polygon, point) == -1) {\n return true;\n }\n})","pinned":true,"hidden":false,"mode":"ojs"},{"id":1077,"value":"usedPolygonAndGrid = useGrid\n ? [...usedPolygon, ...gridPointsInPolygon]\n : usedPolygon","pinned":true,"hidden":false,"mode":"ojs"},{"id":1079,"value":"delaunatorTriangles = {\n const delaunator = Delaunator.from(usedPolygonAndGrid)\n \n const coordinates = [];\n for (let i = 0; i < delaunator.triangles.length; i += 3) {\n coordinates.push([\n usedPolygonAndGrid[delaunator.triangles[i]],\n usedPolygonAndGrid[delaunator.triangles[i + 1]],\n usedPolygonAndGrid[delaunator.triangles[i + 2]]\n ]);\n }\n return coordinates;\n}","pinned":true,"hidden":false,"mode":"ojs"},{"id":1169,"value":"edgeIndices = {\n const indices = [];\n for (let i = 0; i < usedPolygon.length - 1; i++) {\n indices.push([i, i + 1]);\n }\n indices.push([usedPolygon.length - 1, 0]);\n return indices;\n}","pinned":true,"hidden":false,"mode":"ojs"},{"id":1163,"value":"constrainautorTriangles = {\n if (polygon.length < 3) {\n return [];\n }\n const del = Delaunator.from(usedPolygonAndGrid);\n const con = new Constrainautor(del, edgeIndices);\n\n const coordinates = [];\n for (let i = 0; i < con.del.triangles.length; i += 3) {\n coordinates.push([\n usedPolygonAndGrid[con.del.triangles[i]],\n usedPolygonAndGrid[con.del.triangles[i + 1]],\n usedPolygonAndGrid[con.del.triangles[i + 2]]\n ]);\n }\n return coordinates;\n}","pinned":true,"hidden":false,"mode":"ojs"},{"id":1107,"value":"insideTriangles = {\n const coordinates = [];\n for (let i = 0; i < constrainautorTriangles.length; i++) {\n const triangleMidPoint = [\n (1 / 3) * constrainautorTriangles[i][0][0] +\n (1 / 3) * constrainautorTriangles[i][1][0] +\n (1 / 3) * constrainautorTriangles[i][2][0],\n (1 / 3) * constrainautorTriangles[i][0][1] +\n (1 / 3) * constrainautorTriangles[i][1][1] +\n (1 / 3) * constrainautorTriangles[i][2][1]\n ];\n\n if (robustPointInPolygon.default(polygon, triangleMidPoint) == -1) {\n coordinates.push(constrainautorTriangles[i]);\n }\n }\n return coordinates;\n}","pinned":true,"hidden":false,"mode":"ojs"}],"files":"$W30","stdlib":"1","titling":"auto"}