Voronoi-Delaunay mesh

from Red Blob Games
DRAFT
15 Jun 2018
1. Delaunator manual https://mapbox.github.io/delaunator/[1]

People know how to work with square grids; I have a page about hexagonal grids. Right now I have a map generation page with voronoi/delaunay. But it's voronoi + how to use voronoi for maps + specific algorithms for rotmg maps. Could the non-rotmg topics be a topic as detailed as the hexagon page?

• how do you generate the points? grid, jitter, poisson disc, lloyd relaxation, sobol sequences, etc. (see https://trello.com/c/5tFjwMD2/209-2018-jul-point-selection-for-voronoi-maps[3])
• triangle quality, inserting circumcenters for narrow triangles
• what are the core concepts? regions, sides, corners, neighbors, and mappings between them all (refer to grid article)
• how can we represent this? dual mesh is one way but polygon mesh is another, or explicit graph
• what are the operations? how are they implemented? (refer to delaunator documentation)
• how do we generate boundary points and/or ghost elements?
• what are the disadvantages compared to tesselations? no more relative positions, no more vector operations (add, subtract, scale, rotate)
• circumcenter vs centroid vs incenter etc
• subdivision (multiple levels of detail, variable density)
• region to wedges
• cellular automata? (e.g. to put beaches between water and land)
• how to tell if it's voronoi or something similar

{ use two column layout? }

Square grids and hex grids (“tesselations[4]”) are common in games. We can generalize from these to Voronoi tesselations[5]. When the diagram's red points are equally spaced, we get a grid. When we use an irregular spacing, we can get a cool organic look. Voronoi can produce both. On this page I'll describe how to work with Voronoi and related structures.

{ canvas might be better than svg for this page; need to make an inventory of diagrams – do I need hover etc. more often, or do I want larger draw counts? it would also be nice to produce images for search engines }

Given these red points, we can construct a Voronoi diagram. {khan academy animation here?}

{ should be possible to smoothly vary from voronoi to centroid in the same diagram, by using x to decide on the mix }

we lose

• tiles are all the same shape, useful for art
• coordinates
• infinite space, because we no longer have coordinates
• fast distances, because we no longer have coordinates
• directions, because polygons have varying shapes

However, we retain neighbors and neighbor-related algorithms like pathfinding, and we gain a cool organic look. We can vary the density and shapes of the cells.

I'll start with Voronoi, but will also describe a related system that has even nicer properties for game maps.

{ might be nice for interactive diagram allow moving points around }

{ show diagrams for square, hex, triangle, hybrid, curved, blue noise, and random }

1  (Gallery)#

Put a bunch of examples here

{ this might look cool, especially , but it's not the primary focus }

http://mkmra2.blogspot.com/2016/10/[8] "smooth voronoi"

2  Dual mesh#

There are two meshes here.

1. Voronoi: there are the red polygons with blue corners, and white sides
2. Delaunay: there are the blue triangles with red corners, and black sides

This is a dual mesh, so the central points in one mesh also act as the vertices for the other mesh. The red points are labels for the regions, and vertices for the triangles. The blue points are labels for the triangles and vertices for the regions. This is weird and it will take some time to get used to.

{ this diagram should not have the centroid/circumcenter slider, but instead a voronoi/delaunay slider }

What should I call the red and blue points? The names “corners” and “vertices” are ambiguous. The word “polygons” also includes triangles. The word “tiles” suggests polygons with the same shape. I've decided to call them regions and triangles.

red points
These are the inputs to Delaunay triangulation, and will get integers 0 ≤ `r` < `numRegions`.
delaunay triangles
These are the outputs of Delaunay triangulation, and will get integers 0 ≤ `t` < `numTriangles`.
black lines
These are also the outputs of the Delaunay triangulation, and get integers 0 ≤ `s` < `numSides`. There are three sides per triangle, so `numSides` = 3 * `numTriangles`. They're represented as half-edges, so adjacent triangles have separate sides instead of sharing a side ID.
blue points
There's one corner of a Voronoi polygons for each of the triangles, so they'll get the same IDs.
voronoi regions
There's one region for each of the red points, so they'll get the same IDs.
white lines
There's one side of a Voronoi regions for each side of the triangles (it may not be obvious), so they'll get the same IDs. These too are half-edges, so adjacent regions will have separate sides instead of sharing a side ID.

The Voronoi regions do not require any additional IDs beyond what the Delaunay triangulation uses.

To summarize, there are three types of ids: regions, triangles, and sides. Region IDs and triangle IDs are unrelated.

NameVoronoiDelaunayId
faces/cells/tiles/polygons0 ≤ `r` < `numRegions`
vertices/triangles0 ≤ `t` < `numTriangles`
edges/borders/sides0 ≤ `s` < `numSides`

{ this section could be condensed - seems a little redundant - but it's a key concept so repetition may be useful }

2.1 Parts

There's no coordinate system for an irregular grid, so instead we can either assign IDs to each part or make objects for each part. I prefer to assign numeric IDs. This lets me store additional data in arrays. For example if I want to store the elevation with each region, I can have a separate array `elevation[r]` indexed by the region number.

{ are these diagrams still useful? could they be combined into one that shows the labels on mouseover? or can I use the left/right trick from above to show all three in the same diagram? }

2.1.1 Regions

"tiles", "regions", "polygons", "cells", "faces", ("vertices" in dual)

2.1.2 Corners

"vertices", "corners", ("triangles" in dual)

2.1.3 Sides

"edges", "sides", "borders"

2.2 Delaunator library#

I'll describe the Delaunator library's data structures here. Other libraries will have different representations. We need ways to get from triangles to sides, sides to triangles, triangles to regions, regions to triangles, regions to sides, and sides to regions. Delaunator provides the building blocks, and then I'll build the rest in the next section.

2.2.1 Triangle and Side IDs

Delaunator uses side (edge) IDs, and puts them in groups of three for each triangle. This means `triangleId` = `floor(sideId / 3)`, and `sideID` = `3*triangleId + {0,1,2}`. It also means `numSides` = `3 * numTriangles`.

This lets us iterate the sides of a triangle. Given a side ID, we can find the triangle it's part of. Given a triangle ID, we can find the three sides of it. Given a side, we can find the next and previous sides in the triangle.

```function s_prev_s(s) { return (s % 3 == 0) ? s+2 : s-1; }
function s_next_s(s) { return (s % 3 == 2) ? s-2 : s+1; }
```

2.2.2`triangles` array

Delaunator gives us a `triangles[sideId]` array, returning the region (vertex) id for the beginning of that side (half edge).

For example `triangles[s1]` is `r1`. (The numbers being the same is a coincidence.) This array lets us get from sides to regions. Since we can already get from triangles to sides, we also can use this to get from triangles to regions.

2.2.3`starts` array

For some algorithms I need the inverse of the `triangles` array. Delaunator doesn't provide this but it provides everything I need to build it, so I build `starts[regionId]` to tell me any side (half edge) that starts at a given region (vertex).

```let starts = new Int32Array(numRegions);
for (let s = 0; s < triangles.length; s++) {
starts[triangles[s]] = s;
}
```

For example `starts[r2]]` can be `s14`. It could've been `s0` or `s5` or `s9` or `s13`. It doesn't matter which (unless you don't use ghosts, which I'll explain later). This gets us from regions to sides and triangles.

1. TODO non-ghost version

There's a problem on the outer hull of the triangulation. We're going to use the `starts` array to iterate through the halfedges, but some halfedges may be -1. The solution is rewind each element by iterating backwards until we hit a -1.

{ need to write code for this and test it thoroughly }

```function s_rewind_s(start_s) {
// TODO: TEST THIS!! I don't think it's right
let s = null,
next_s = start_s;
do {
s = next_s;
let opposite_s = halfedges[s];
next_s = s_prev_s(opposite_s);
} while (s != start_s && next_s != -1);
return s;
}
```

or maybe

```function s_rewind_s(start_s) {
let s = start_s;
while (true) {
let back_s = halfedges[s_prev_s(s)];
if (back_s === -1 || back_s === start_s) {
break;
}
s = back_s;
}
return s;
}
```

and then use

```let starts = new Int32Array(numRegions);
for (let s = 0; s < triangles.length; s++) {
starts[triangles[s]] = s_rewind_s(s);
}
```

2.2.4`halfedges` array

Delaunator provides a `halfedges[sideId]` array, returning another `sideId` for the half edge's “opposite pair”.

For example `halfedges[s2]` is `s5`, and `halfedges[s5]` is `s2`. This is how we'll get from one triangle to adjacent triangles.

3  Traversal algorithms#

Delaunator provides the building blocks. Here I'll describe how I use its data structures for traversing the dual mesh. 9 of these, and my names for them

3.1 TODO Halfedge properties#

Halfedge id to start, end r

Halfedge id to inner, outer t

3.2 Delaunay triangle circulation#

There are two main patterns to look at. Inside a triangle, we can circulate around the sides (halfedges):

This function gives us the three sides of a triangle, which have consecutive IDs starting at 3 times the triangle ID:

```function t_circulate_s(t) {
let out_s = [];
for (let i = 0; i < 3; i++) {
let s = 3*t + i;
out_s.push(s);
}
return out_s;
}
```

Given a side we can find the region point by looking in the `triangles` array. This lets us get the three region points that make up a triangle:

```function t_circulate_r(t) {
return t_circulate_s(t).map(s => triangles[s]);
}
```

To go between triangles, we can jump between a side (halfedge) and its opposite:

This lets us find the adjacent triangles around a given triangle:

```function t_circulate_t(t) {
return t_circulate_s(t).map(s => floor(halfedges[s]/3));
}
```

{ note that this doesn't work without ghosts }

If we start with a triangle ID we can use this loop to find sides, triangles, or regions.

3.3 Voronoi region circulation#

Another useful way to combine these two operations is to circulate one step, jump to the opposite, circulate one step, jump to the opposite, repeat:

This is tricky, so let's list the steps to get the sides leading out of the center region point:

1. Find any halfedge leading out from the starting point. One such halfedge is `s14`. (use the starts array for this)
2. Find `s14`'s opposite, which is `s2`. (use the halfedges array)
3. Circulate around `s2`'s triangle, finding `s0`. (use the s_next_s helper function)
4. Find `s0`'s opposite, which is `s5`.
5. Circulate around `s5`'s triangle, finding `s3`.
6. Find `s3`'s opposite, which is `s8`.
7. Circulate around `s8`'s triangle, finding `s6`.
8. Find `s6`'s opposite, which is `s10`.
9. Circulate around `s10`'s triangle, finding `s11`.
10. Find `s11`'s opposite, which is `s13`.
11. Circulate around `s13`'s triangle, finding `s14`.
12. We've come back to `s14`, so we're done. The outgoing halfedges were `[s14, s0, s3, s6, s11]`. The incoming halfedges were `[s2, s5, s8, s10, s13]`.

Here's what the code looks like for finding outgoing sides from a region point:

```function r_circulate_s(r) {
let out_s = [];
let start_s = starts[r];
let s = start_s;
do {
out_s.push(s);
let opposite_s = halfedges[s];
s = s_next_s(opposite_s);
} while (s != start_s);
return out_s;
}
```

The same loop works for finding the triangles touching a region point, or finding the vertices to draw the voronoi polygon:

```function r_circulate_t(r) {
return r_circulate_s(r).map(s => floor(s/3));
}
```

or finding the neighboring regions adjacent to a given region:

```function r_circulate_r(r) {
return r_circulate_s(r).map(s => triangles[s]);
}
```

That's it! Those are the two basic loops we need to find our way around the dual mesh. In both cases we focus on the halfedges and then construct the region or triangle data from them.

{ animation might be useful here }

{ to test: if we don't have ghosts, then the loop might work nicely if we stop at halfedges[start_s] >= 0? s_prev_s(start_s) : -1 }

3.4 TODO Drawing

Show sample code for drawing the delaunay triangles

Show sample code for drawing the voronoi regions

{ diagrams corresponding to sample code }

3.5 TODO Distances

Have to explore the graph with breadth first search; no coordinates means no fast formula.

{ diagram where you mouseover any region and it shows distances to all other regions }

4  Generation#

4.1 TODO point selection

{ this seems like an independent problem }

4.2 TODO delaunay

I often want to work with the centroid version, not the circumcenters. It's simpler for me to start with Delaunay and then calculate the centroids than for me to use a Voronoi library.

Circumcenters and orthocenters are not guaranteed to be inside the triangle. Centroids and incenters are.

{ should this explain what to look for in a delaunay lib? need the half edge and triangle data }

look for a library that gives you the IDs instead of x,y!

does it matter if all the triangles are oriented the same way? this might be a cute animation, showing not only the circulation order but have that sprite spin when you mouse over a triangle

Some libraries use iterators or arrays so that you can get an integer id.

Some libraries don't give you all the information you need, but you might be able to reconstruct the dual graph.

There are other formats too like http://www.gradientspace.com/tutorials/dmesh3[18]

4.3 TODO triangle centers

4.3.1 circumcenter (voronoi)

{formula for calculating; link to wikipedia}

4.3.2 centroid (barycentric mesh)

{formula for calculating; link to wikipedia}

{ also link to all the other candidates for triangle centers }

4.4 TODO Ghost and solid

Problem: halfedge could be -1. This complicates some of the traversals. {but maybe this is after traversal}

Solution: add "ghost" region on the outside of the mesh. Then add sides and triangles between the ghost region and the solid regions.

Solution: add checks to all the code.

5  Appendix: Relationship between parts#

squares:

• a face has 4 face neighbors, 4 sides, 4 corners
• a corner has 4 corner neighbors, 4 sides, 4 faces
• a side has 2 faces, 2 corners

hexagons:

• a face has 6 face neighbors, 6 sides, 6 corners
• a corner has 3 corner neighbors, 3 sides, 3 faces
• a side has 2 faces, 2 corners

triangles:

• a face has 3 face neighbors, 3 sides, 3 corners
• a corner has 6 corner neighbors, 6 sides, 6 faces
• a side has 2 faces, 2 corners

voronoi:

• a face has N face neighbors, N sides, N corners
• a corner has 3 corner neighbors, 3 sides, 3 faces
• a side has 2 faces, 2 corners

So let's look at our diagram again. I said that square grid corners always have 4 neighbors. But I also said voronoi corners always have 3 neighbors. And I also said square grids have the voronoi property.

All three of these can't be true.

Overlap, maybe centroid slider to show this.

6  Appendix: Centroids and Incenters#

For my map projects I want to put things at the blue points, like roads or castles. When two blue points have the same coordinate, or even are close, things can go wrong.

Let's separate the blue points. We lose the voronoi property and can no longer generate squares. It keeps the same mesh structure and it simplifies other things. More on this later.

I should also try incenter because it seems to have some good properties too (e.g. it's the angle bisector for the wedge of a polygon in a triangle's direction, and the biggest circle inside the triangle is drawn there) and I don't actually know which will be best

Centroid is useful in simulations - if you have f(r) for all the vertices r of a triangle, then the best estimate of f(centroid) is the weighted average of the three f(r) values. And if we build half wedges for polygons then we are bisecting the white edge

Not sure if this will go here, or if it's up near the top.

The square grid is a clear example of the problem. Square grid voronoi produces are square grid. But the triangles tell you connectivity between adjacent polygons. And the triangles tell us that each square has six neighbors, not four!! If we look closely we'll see there are circumcenters that overlap. Two triangles have the same circumcenter. This is messy.

{ diagram with a tiny square grid showing the centroid points joining together when they become a circumcenter }

7  Appendix: Drawing with bezier curves#

Related: https://sighack.com/post/chaikin-curves[19]

8  Appendix: Pixel to polygon#

spatial hash, point-in-poly test

{ diagrams showing what the spatial hash matches, and how we need the point-in-poly to be sure }

{ diagram }

10  Appendix: Subdivision#

Subdivide the triangle, not the regions. Then construct new regions. Not sure if I should have this section.

12  Appendix: Voronoi property#

Voronoi finds all regions closer to a given point than to any other point.

Animation idea from Khan Academy[23] / Peter Collingridge.

13  Appendix: How to tell if it's Voronoi#

{ diagram }

1. Take two adjacent (red) points
2. Draw the black edge between them
3. The white edge should be along the perpendicular bisector

{ should this be a quiz? }

14  NOTES#

• more than structured grids
• less than planar graphs

What to call this?

Twitter suggestions: chickenwire, cellular grid, baked mudflat, a pretty butterfly, foam, bubbles, cellgrid, polygrid, polymesh, convex irregular polygon tiling, tiling, net, voronoi grid, cell structure ngram viewer suggest unstructured grids beats irregular grids, but polygon mesh is also popular.

I decided on "voronoi" because it's what gamedevs will use (especially for procedural maps, which is my main application). https://en.wikipedia.org/wiki/Tessellation#Voronoi_tilings[28]

"Dirichlet tessellations"

"Thiessen polygons"