2e:null 2f:["$","title",null,{"children":"Contouring Unstructured Data Using Meandering Triangles | Corey Guastini | Observable"}] 34:T4b1,function interpolate(cat, v1, v2, v3, lev) { let segment, p1, p2; //make sure larger values are always to the same side of the segment //to ensure proper connection of head to tail switch (cat) { default: segment = null; break; case 0: segment = null; break; case 1: p1 = ipolate(lev, v1, v3); p2 = ipolate(lev, v1, v2); segment = {p1:p1, p2:p2, visited:false}; break; case 2: p1 = ipolate(lev, v2, v1); p2 = ipolate(lev, v2, v3); segment = {p1:p1, p2:p2, visited:false}; break; case 3: p1 = ipolate(lev, v1, v3); p2 = ipolate(lev, v2, v3); segment = {p1:p1, p2:p2, visited:false}; break; case 4: p1 = ipolate(lev, v3, v2); p2 = ipolate(lev, v3, v1); segment = {p1:p1, p2:p2, visited:false}; break; case 5: p1 = ipolate(lev, v3, v2); p2 = ipolate(lev, v1, v2); segment = {p1:p1, p2:p2, visited:false}; break; case 6: p1 = ipolate(lev, v2, v1); p2 = ipolate(lev, v3, v1); segment = {p1:p1, p2:p2, visited:false}; break; case 7: segment = null; break; } return segment; }35:T795,md`Then, for each contour level desired, the code loops through each triangle and places it into one of eight categories. The eight categories are defined based on whether the value at each vertex of the triangle is greater than or equal to the given contour value. [Wikipedia](https://en.wikipedia.org/wiki/Marching_squares) describes the categories well under the "Contouring triangle meshes" heading. If the code determines that a contour runs through a triangle, it performs interpolation along the triangle edges to determine the starting and ending points of the contour segment. Because the interpolation is always along a triangle edge, its formulation is simple, but a more in-depth exploration of [interpolating in a triangle](https://codeplea.com/triangular-interpolation) helped me visualize the process. At the conclusion of checking each triangle for a given contour level, we're left with a collection of contour segments that need to be joined into contours. The joining can happen via brute force (checking each segment to see if it has a starting point equal to another segment's ending point), but that's slow. Instead, the segments are entered into two [d3 quadtrees](https://github.com/d3/d3-quadtree), one generated using segment starting point and the other using segment ending point. Quadtrees are used because they make nearest-neighbor searching much faster than brute force. For a given segment ending point, the starting-point quadtree is searched for the closest point. If the returned closest starting point is equal to the search ending point, then the segment should be added to the contour line. This search is repeated until there are no more starting points equal to an ending point. The same larger process is then repeated except the code now searches the ending-point quadtree for the point closest to the contour starting point. The final result is an object of contours that can be easily rendered.`36:Tfd2,function* drawing() { replay; const EPSILON = Math.pow(2, -52); context.clearRect(0,0,context.canvas.width,context.canvas.height); context2.clearRect(0,0,context2.canvas.width,context2.canvas.height); drawCircles(context); drawTriangles(context); //object to hold all contours for all levels let allContours = {}; //holds all contour segments which will be joined into contour lines let levSegs; //quadtrees let quadstart, quadend; //contour levels to compute and plot let contourLevels = linspace(-1,1,11); for (let i = 0 ; i < contourLevels.length ; ++i) { allContours[contourLevels[i]]=[]; levSegs = createSegments(contourLevels[i]); //if no segments for contour level, continue to next level if (levSegs.length == 0) continue; //create two quadtrees, one generated using segment starting point and the other using segment ending point quadstart = d3.quadtree(levSegs, function(d) { return d.p1.x ;}, function(d) { return d.p1.y; }); quadend = d3.quadtree(levSegs, function(d) { return d.p2.x ;}, function(d) { return d.p2.y; }); //initialize contour line let contour = []; let testSeg; //loop through all segments for (let j = 0; j < levSegs.length ; ++j) { //if segment has already been tested for inclusion into contour, continue if (levSegs[j].visited == true) continue; contour.push(levSegs[j].p1); allContours[contourLevels[i]].push(contour); testSeg = levSegs[j]; testSeg.visited = true; //find segment with starting point closest to this segment ending point let closest = quadstart.find(testSeg.p2.x, testSeg.p2.y); //if found segment sufficiently close (location should be equal, but floating point error means we can't compare using equality) and not already tested for inclusion into contour while (Math.abs(closest.p1.x - testSeg.p2.x) < EPSILON && Math.abs(closest.p1.y - testSeg.p2.y) < EPSILON && !closest.visited) { closest.visited = true; contour.push(closest.p1); testSeg = closest; closest = quadstart.find(testSeg.p2.x, testSeg.p2.y); //replace contour with updated contour allContours[contourLevels[i]][allContours[contourLevels[i]].length-1] = contour; //render updated contour yield renderContour(context2, allContours); } //check if contour needs to be closed if (Math.abs(testSeg.p2.x - contour[0].x) < EPSILON && Math.abs(testSeg.p2.y - contour[0].y) < EPSILON) { contour.push(contour[0]); contour.push(contour[1]); allContours[contourLevels[i]][allContours[contourLevels[i]].length-1] = contour; yield renderContour(context2, allContours); } //reverse contour to test the other direction let reversedContour = contour.reverse(); testSeg = levSegs[j]; closest = quadend.find(testSeg.p1.x, testSeg.p1.y); while (Math.abs(closest.p2.x - testSeg.p1.x) < EPSILON && Math.abs(closest.p2.y - testSeg.p1.y) < EPSILON && !closest.visited) { closest.visited = true; reversedContour.push(closest.p1); testSeg = closest; closest = quadend.find(testSeg.p1.x, testSeg.p1.y); allContours[contourLevels[i]][allContours[contourLevels[i]].length-1] = reversedContour; yield renderContour(context2, allContours); } //check if contour needs to be closed if (Math.abs(testSeg.p2.x - reversedContour[0].x) < EPSILON && Math.abs(testSeg.p2.y - reversedContour[0].y) < EPSILON) { reversedContour.push(reversedContour[0]); reversedContour.push(reversedContour[1]); allContours[contourLevels[i]][allContours[contourLevels[i]].length-1] = reversedContour; yield renderContour(context2, allContours); } //reset contour contour = []; } } //renderContour(context2,allContours); yield allContours; }37:[] 32:{"id":"419459ccee76eb3c","slug":"contouring-unstructured-data-using-meandering-triangles","title":"Contouring Unstructured Data Using Meandering Triangles","thumbnail":"8156f5287a3d4c978e9ad0f770940b7aa946db750cb24ecaa695fa5cfec3c91b","access_level":"public","create_time":"$D2019-12-20T16:33:45.245Z","update_time":"$D2022-02-10T14:20:47.724Z","view_time":null,"first_public_time":"$D2019-12-24T00:50:49.903Z","has_unpublished_edits":null,"fork":false,"stars":10,"starred":false,"owner":{"type":"team","id":"da36068e9c88f131","login":"cguastini","name":"Corey Guastini","avatar_url":"https://avatars-next.observableusercontent.com/7c6d91c2a664c9f78ab78b1e0bf2b539441aaa1d9430ae24e9886ce825f3c1b0","tier":"starter_2024"},"version":801,"publish_version":801,"role":null,"authored":false,"authors":[{"type":"individual","id":"56bb3b86bffb2836","login":"cguastini","name":"Corey Guastini","avatar_url":"https://avatars-next.observableusercontent.com/7c6d91c2a664c9f78ab78b1e0bf2b539441aaa1d9430ae24e9886ce825f3c1b0","tier":null}],"members":[{"type":"individual","id":"56bb3b86bffb2836","login":"cguastini","name":"Corey Guastini","avatar_url":"https://avatars-next.observableusercontent.com/7c6d91c2a664c9f78ab78b1e0bf2b539441aaa1d9430ae24e9886ce825f3c1b0","tier":null,"role":"editor"}],"comments":[],"custom_thumbnail":true,"trashed":false,"lapsed":null,"subscribed":null,"forks":2,"collection":null,"type":"notebook","body":{"title":"Contouring Unstructured Data Using Meandering Triangles","cells":[{"id":0,"value":"md`# Contouring Unstructured Data Using Meandering Triangles`","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":301,"value":"viewof numpoints = slider({\n min: 500,\n max: 10000,\n step: 500,\n value: 5000,\n format: \",\"\n})","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":24,"value":"md`Field sampled using ${numpoints} random points`","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":715,"value":"viewof replay = html`