Hazard

Island Generator November 9, 2016

I love maps, so when I saw a procedurally generated island on Job’s screen, I naturally rolled up next to him to find out more. He showed me a couple example islands, and I grilled him on his algorithm so I could draw my own.

Job’s technique doesn’t create any topography, only coastal outlines. Let’s start with a regular polygon. Here’s code to generate the vertices of a regular n-gon of a given size:


(def tau (* 2 js/Math.PI))

(defn n-gon [radius sides]
  (assert (< 2 sides) "Must have at least 3 sides")
  (for [angle (range 0 (- tau 0.000001) (/ tau sides))]
    [(* radius (js/Math.cos angle))
     (* radius (js/Math.sin angle))]))

(n-gon 50 3)

Notice that range stops just shy of 2π. Floating point math sometimes produces a final angle indistinguishably close to 2π, adding one extra vertex to our n-gon. To see for yourself, change (- tau 0.000001) to tau and set the number of sides to 6. You’ll see an extra point in the output very close to the starting point.

Let’s draw our points.


(draw (n-gon 100 3))

Let’s turn our regular n-gon into an irregular n-gon by sliding each vertex toward or away from the center of the shape.


(defn irregular [pts]
  (for [[x y] pts
        :let [d (+ 1 (rand))]]
    [(* d x) (* d y)]))

(def irregular-ngon (irregular (n-gon 100 3)))

(draw irregular-ngon)

One thing I learned working with Job is that he reveres vector math. Rightly so. Pushing and pulling the points of the polygon along the line they form with the center is as easy as treating each vertex in the polygon as a vector, then scaling it. In this case, I chose to scale it between 0.5 and 1.5 times its original length. No trig necessary.

To turn our irregular n-gon into something more island-like, we’ll subdivide each edge by offsetting its midpoint some distance perpendicular to the edge. Below is a function that, given two points, returns an offset midpoint. The midpoint is offset by a random distance in either direction, no more than half the length of the original line.


(defn offset-midpoint [[x1 y1] [x2 y2]]
  (let [mx (/ (+ x1 x2) 2)
        my (/ (+ y1 y2) 2)
        vx (- (- y1 y2))
        vy (- x1 x2)
        d (- (rand) 0.5)]
    [(+ mx (* d vx))
     (+ my (* d vy))]))

offset-midpoint takes two points and calculates the coordinates of the midpoint as the average of the source points’ coordinates. vx and vy describe a vector perpendicular to the line between the source points. It’s a two-step calculation folded into one pair of expressions: 1) get the vector from the first point to the second, then 2) rotate the vector 90° by swapping the coordinates of the vector and negating one of them.

Let’s subdivide our edges to add sides to our polygon. I’ll give you the code first and explain it afterward.


(defn loop-around [pts]
  (concat pts [(first pts)]))

(defn subdivide-with [f pts]
  (->> pts
    (partition 2 1)
    (mapcat (fn [[p1 p2]]
              [p1 (f p1 p2)]))))

(def subdivided-ngon
  (subdivide-with
    offset-midpoint
    (loop-around irregular-ngon)))

(vec subdivided-ngon)

The utility function loop-around adds the first point in the sequence to the end. So far we haven’t had to worry about connecting the end of our polygon to the beginning because the drawing function and HTML5 canvas API completed the loop for us. Because we convert the sequence of points to pairs of points using partition, we need to circle back to the beginning ourselves. If we don’t, we won’t get the final pair of points representing the last edge of the polygon. (We’ll see later why I chose to pass the completed loop as an argument to subdivide-with rather than calling loop-around within subdivide-with.)

The work inside subdivide-with is done by mapcat. The function passed to mapcat takes our partitioned pair of points (edge endpoints) and returns two points: the edge’s first endpoint and the result of passing the endpoints to some other function, in this case offset-midpoint. Notice that the edge’s second endpoint isn’t returned. It will be included when the next pair of points is processed, and we don’t want to include it twice. This strategy leaves off the final point (the one we put on the end with loop-around), but the draw function still completes the polygon for us.


(draw subdivided-ngon)

Let’s subdivide again.


(draw (subdivide-with offset-midpoint (loop-around subdivided-ngon)))

It’s starting to look more like a proper island! Let’s go a lot deeper.


(draw (->> subdivided-ngon
        loop-around (subdivide-with offset-midpoint)
        loop-around (subdivide-with offset-midpoint)
        loop-around (subdivide-with offset-midpoint)
        loop-around (subdivide-with offset-midpoint)
        loop-around (subdivide-with offset-midpoint)
        loop-around (subdivide-with offset-midpoint)
        loop-around (subdivide-with offset-midpoint)
        loop-around (subdivide-with offset-midpoint)
        loop-around (subdivide-with offset-midpoint)))

We’ve subdivided the edges 10 times, and it looks great!

As you generate islands, you might notice that the polygon’s edges sometimes cross other edges. This is an advantage as it gives much more interesting islands. By setting the fill rule to “evenodd”, overlapping parts of the island turn into lakes and rivers!

Let’s do some refactoring. Calling loop-around and subdivide-with over and over has several disadvantages. While the code is awkward, we could clean it up with a loop/recur. The bigger problem is that we’re walking a growing sequence over and over. There’s also a conceptual mismatch: we really don’t care how many times we subdivide, we just want the edges to shrink far enough that we don’t see straight, unnatural-looking lines.

To fix all three of these problems, let’s redefine subdivide-with to recursively subdivide edges until they’re short enough.


(defn dist [[x1 y1] [x2 y2]]
  (js/Math.sqrt
    (+ (js/Math.pow (- x1 x2) 2)
       (js/Math.pow (- y1 y2) 2))))

(defn subdivide-with [f pred pts]
  (->> pts
    (partition 2 1)
    (mapcat (fn [[p1 p2]]
              (if (pred p1 p2)
                (subdivide-with f pred [p1 (f p1 p2) p2])
                [p1])))))

I added a dist function to calculate the distance between two points. We’ll use this to check if we should subdivide the edge between two points or whether they’re close enough and we should leave them as is.

subdivide-with now takes three arguments instead of two. The new second argument is a predicate to check if two points should be subdivided. Instead of having the function passed to mapcat unconditionally return the first endpoint and an offset midpoint, the function first asks if the two points should be subdivided, then returns either a subdivided sequence or the first endpoint. The sequence of points given in the recursive call to subdivide-with is where we do the actual subdividing; it introduces the offset midpoint.

Now you can see why I chose to put loop-around outside subdivide-with: I’m using subdivide-with recursively, and I don’t want each subdivided segment to loop back to the first point.

Note how abstract subdivide-with is. If I were to change the representation of points from two-item vectors to a record, nothing within subdivide-with would have to change. In fact, it doesn’t even show any sign that we’re dealing with points. While I don’t have any alternate uses for subdivide-with on hand, it’s helpful not having to worry about both points and recursive subdivision within a single function.

Let’s see how our refactoring looks:


(draw (subdivide-with
        offset-midpoint
        #(< 2 (dist %1 %2))
        (loop-around irregular-ngon)))

Looks good! Generate some islands and see how you like it. There are a couple parameters you can play with, such as the distance midpoints are offset and whether they tend to be pushed into the island or out from it. You can tweak both of these behaviors in the definition of offset-midpoint.

As you play, you might notice that all the islands have about the same coastal roughness. In fact, each island has about the same roughness all the way around. It looks a little bit like a fractal, probably because it is one. Shawn, who has experience in procedural generation, pointed this out and directed my attention to the coast of Ireland:

The eastern and western shores of the Emerald Isle are very different. The east is quite smooth compared to the west, which has lots of islands, bays, and finger-like projections.

The last improvement I want to make is to give our coastlines this kind of real-world regional variation.

To do that, we’ll need to add some metadata to each point. Let’s wrap each of our irregular n-gon’s vertices in a map, including in it two random values: one for whether midpoints are more likely to be offset inward or outward, and another for how far midpoints can be offset.


(defn with-roughness [pt]
  {:position pt
   :balance (rand)
   :max-offset (+ 0.05 (rand))})

(def roughened-ngon (map with-roughness irregular-ngon))

(vec roughened-ngon)

We’ll need a new midpoint function that understands this structure. Each midpoint will average the :balance and :max-offset values of its parents and use them in the calculation of its offset.


(defn offset-midpoint-rough
  [{[x1 y1] :position b1 :balance mo1 :max-offset}
   {[x2 y2] :position b2 :balance mo2 :max-offset}]
  (let [mx (/ (+ x1 x2) 2)
        my (/ (+ y1 y2) 2)
        vx (- (- y1 y2))
        vy (- x1 x2)
        balance (/ (+ b1 b2) 2)
        max-offset (/ (+ mo1 mo2) 2)
        d (* max-offset (- (rand) balance))]
   {:position [(+ mx (* d vx))
               (+ my (* d vy))]
    :balance balance
    :max-offset max-offset}))

Let’s see what our islands look like now.


(draw
  (map :position
       (subdivide-with
               offset-midpoint-rough
               #(< 2 (dist (:position %1) (:position %2)))
               (loop-around roughened-ngon))))

You can see how the roughness of the coastline is different on different parts of the island.

Let’s try it with a random number of sides on the original polygon.


(defn island []
  (map :position
       (subdivide-with
               offset-midpoint-rough
               #(< 2 (dist (:position %1) (:position %2)))
               (loop-around
                 (map with-roughness
                      (irregular (n-gon 100 (+ 3 (rand-int 17)))))))))

(draw (island))

Excellent!

This algorithm generates a wide range of island shapes, from simple, rounded islands with smooth coastlines to craggy-looking islands full of deep coves and rocky shoals, to long, branching archipelagoes. I especially like the coastal rivers and lagoons that occasionally appear.

My original island generator is here, the only meaningful differences being that islands are zoomed to fit the window and you can bookmark individual islands for coming back to or for sharing.

The original code is on GitHub. You’re welcome to suggest improvements.

UPDATE: John Earnest has implemented both algorithms in K using iKe.