3 Jun 2017

Note: the /x/ series of posts on my site are mostly unpolished mini posts. The polished version of this will come later, maybe in a few months.

The previous blog post was about Voronoi and Delaunay graphs, how I used them for procedurally generating maps in 2010, and how I want to use them differently for my next procedural map generator. I use these unstructured grids instead of regular grids to add variety and interestingness to the maps.

I need a way to represent this graph, including many of the relationships listed here, except adapted to this dual graph:

Back in 2010, I represented each of these as a separately allocated array. I didn’t know any better. The people who work in computational geometry have much better ways to represent these kinds of “dual graphs”. The structure of the graph is the same whether I’m using Voronoi (circumcentral coordinates) or the graphs I ended up with (barycentric / centroid coordinates). I only need to represent one half of the dual graph, since the parts correspond:

There are many possible representations, including quad edges, split edges, winged edges, half edges, and others. I decided to use the “half edge” data structure for this project. Read about them on flipcode.com or martindevans.me.

#Structure

Let’s focus on the triangles. The red points are the vertices. They come from the seed points, so I’m going to call give them an “s-index”, which also identifies the white polygons. Let’s give the triangles a “t-index”, which also identifies the blue points. Edges will be represented as directed (“half”) edges, each with an “e-index”, which also lets us construct the white edges.

How will I represent this graph? For the structure I need the relationships between seeds and edges:

  1. edges: edges → seeds that are the start point of the edge. For example edge e2 is s{{graph.s_begin(2)}}s{{graph.s_end(2)}} so edges[e2] = s{{graph.edges[2]}} (where the edge starts).
  2. opposites: edges → edges that are the same location, opposite direction. For example edge e2 is s{{graph.s_begin(2)}}s{{graph.s_end(2)}} and e{{graph.opposites[2]}} is s{{graph.s_begin(graph.opposites[2])}}s{{graph.s_end(graph.opposites[2])}} so opposites[e2] = e{{graph.opposites[2]}}. Similarly, opposites[e{{graph.opposites[2]}}] = e2. Some edges won’t have opposites.
  3. starts: seeds → edges : an edge that starts with this seed. For example, seed s2 has edges e2 and e3 coming from it, so either of those can go into the array. For interior seed points it doesn’t matter whether which. For exterior seed points we want the “leftmost” edge, so starts[s2] = e2.

I also have two arrays for the geometry of the graph:

  1. s_coords[s] is the x,y for each red point s
  2. t_coords[t] is the x,y for each blue point t

Most of the procedural map generation steps will work with the structure; I will use the geometry when I’m drawing the map.

#Example

There are exactly three edges for each triangle. That means I can work with edges most of the time: t_index(e_index) = int(e_index/3). Here’s the edges mapping for the above diagram:

Here’s the opposites mapping:

And here’s the starts mapping:

Since the s-indices and e-indices are small integers, these mappings can be stored as a dense array of integers. This is memory efficient and also cache friendly. For example, opposites[e2] = e5 will be stored as opposites[2] = 5. If there’s no opposite, store -1.

Bonus: in GL_TRIANGLES indexed mode drawing, edges will be the indices and s_coords will be the vertices.

#Operations

Try mousing over or tapping any red or blue point:

  • Blue point (t) to one triangle (e1 e2 e3) it’s contained in:
    3*t+i for 0 ≤ i < 3
  • Blue point (t) to three polygons (s1 s2 s3) it is a corner of:
    edges[3*t+i] for 0 ≤ i < 3
  • Blue point (t) to three blue points it connects to:
    t_index(opposites[3*t+i]) for 0 ≤ i < 3 (*)
  • Red point (s) to outgoing edges (e1 e2 …):
    e for e ∈ e_loop(s)
  • Red point (s) to incoming edges (e1 e2 …):
    PREV(e) for e ∈ e_loop(s)
  • Red point (s) to red points (s1 s2 …) it connects to (**):
    s_end(e) for e ∈ e_loop(s)
  • Red point (s) to many triangles (t1 t2 …) it’s a corner of:
    t_index(e) for e ∈ e_loop(s)
  • Half edge (e) to white edge’s two blue endpoints: (no demo)
    t_index(e), t_index(opposites[e]) (*)
  • Half edge (e) to black edge’s two red endpoints: (no demo)
    edges[e], edges[NEXT(e)]
    function PREV(e) { return (e % 3 == 0) ? e+2 : e-1; }
    function NEXT(e) { return (e % 3 == 2) ? e-2 : e+1; }
    function e_loop(s) {
        const e0 = starts[s];
        let e = e0, e1 = -1;
        do {
            yield e;
            e1 = opposites[e];
            e = NEXT(e1);
        } while (e1 >= 0 && e != e0);
    }

(*) There are special cases we have to handle at the exterior edges because they don’t have opposites.

(**) This algorithm as described doesn’t actually work for exterior red points. See the next section for a way to fix it.

Although there are a lot of operations listed above, they can be put into two basic patterns:

  1. Circulate around a blue point using the 3 indices for a triangle’s black edges, using the NEXT (counterclockwise) or PREV (clockwise) functions.
  2. Circulate around a red point using the e_loop function. This will visit (clockwise) each black edge connected to the red point.

I’m still looking for good names for these functions; I should look at libraries like jCAE and CGAL and OpenMesh. I don’t yet know whether I want to think about blue and red points, with the triangles and polygons secondary, or if I want to think about triangles and (red) vertices, with the blue points and polygons secondary.

#Special cases

In the previous section there are some functions marked with (*) where we might walk outside the exterior and get a null (-1). Exterior seeds also complicate things for the starts[] array, since we have to treat interior and exterior seeds differently. The red point to red points it connects to function is missing one red point if you start with an exterior seed. Ugh! Most half edges have an opposite but some don’t, so you have to be careful, and sometimes you can’t even get the information you want.

With grids there’s a trick to avoid some of these situations: add an extra row and column to the outer boundaries of the map. When you try to access array[x-1][y] you access a special boundary grid location that is filled in with default values.

There’s a similar trick for polygon meshs which I found in several papers: make new triangles around the outer boundary of the mesh. It goes by different names: “outside” vertex, “infinity” vertex, “ghost” vertex. If there are N half edges that are missing an opposite:

  1. Create N ghost triangles. These will be on the outer edges of the map.
  2. Create 3*N ghost edges. Of each three, one will fill in the missing half edge. The other two are used to build the triangles.
  3. Create 1 ghost seed. This seed will be outside the map area, and is the third vertex for the added triangles.

Think of this putting an extra seed on the back of your map, and then putting extra triangles around the boundary of your map that wrap around to the back.

Ok, so now we’ve replaced all the checks for -1 by checks for things that go to the back of the map. Isn’t it still annoying? Yes, but it’s less annoying.

Since all the extra edges are added to the end of the existing edges, it’s easy to know whether it’s a ghost edge: e > numSolidEdges. Similarly, the ghost seed is added to the end of the seed array, so it’s easy to identify.

Filtering out all the ghosts seems annoying. But so does dealing with nulls / -1s. Do I need to filter out ghost structures at every step? No, I think most of the map generation code won’t care which structures are solid and which are ghost. I only have to filter the ghosts at the very end, when drawing the map. I think that’s a win.

#More reading

I’m a computational graphics newbie. I’m learning about this stuff while working on the map generator project. I don’t know all the different mesh formats and which works best for my application, but this seems like it’ll be pretty good. Here’s some reading:

I usually use D3-voronoi for my projects but this time I’m trying the Delaunator library to construct the half edge structure. Not only is it very fast, it returns information in exactly the format I wanted for the edges and opposites arrays. Nice. I’ll have to extend these arrays if I want to use the “infinity” point to reduce special case code.