content/blog/2016/06/diamond-square.markdown @ de64f1b0d5c0
Add public to .hgignore.
| author | Steve Losh <steve@stevelosh.com> |
|---|---|
| date | Fri, 07 Oct 2016 12:50:31 +0000 |
| parents | e7bc59b9ebda |
| children | fdf01e99fd51 |
+++ title = "Terrain Generation with Diamond Square" snip = "Improving on Midpoint Displacement." date = 2016-06-27T13:35:00Z draft = false +++ <script data-cfasync="false" src="/media/js/three.min.js"></script> <script data-cfasync="false" src="/media/js/TrackballControls.js"></script> <script data-cfasync="false" src="/media/js/terrain3.js"></script> In the last two posts we looked at implementing the Midpoint Displacement algorithm for procedurally generating terrain. Today we're going to look at a similar algorithm called Diamond Square that fixes some problems with Midpoint Displacement. The full series of posts so far: * [Midpoint Displacement](/blog/2016/02/midpoint-displacement/) * [Recursive Midpoint Displacement](/blog/2016/03/recursive-midpoint-displacement/) * [Diamond Square](/blog/2016/06/diamond-square/) Midpoint Displacement is simple and fast, but it has some flaws. If you look at the end result you'll probably start to notice "seams" in the terrain that appear on perfectly square chunks. This happens because when you're calculating the midpoints of the edges of each square you're only averaging two sources of information. [Diamond Square][] is an algorithm that works a lot like Midpoint Displacement, but it adds an extra step to ensure that (almost) every point uses *four* sources of data. This reduces the visual artifacts a lot without much extra effort. Let's [draw the owl][owl]. [Diamond Square]: https://en.wikipedia.org/wiki/Diamond-square_algorithm [owl]: https://i.imgur.com/RadSf.jpg {{% toc %}} ## Overview The high-level description of Diamond Square goes something like this: 1. Initialize the corners to random values. 2. Set the center of the heightmap to the average of the corners (plus jitter). 3. Set the midpoints of the edges of the heightmap to the average of the four points on the "diamond" around them (plus jitter). 4. Repeat steps 2-4 on successively smaller chunks of the heightmap until you bottom out at 3x3 chunks. ## Initialization Initialization is pretty simple: ```clojure (defn ds-init-corners [heightmap] (let [last (heightmap-last-index heightmap)] (heightmap-set! heightmap 0 0 (rand)) (heightmap-set! heightmap 0 last (rand)) (heightmap-set! heightmap last 0 (rand)) (heightmap-set! heightmap last last (rand)))) ``` Just set the corners to random values, exactly like we did in Midpoint Displacement. <div id="demo-1" class="threejs"></div> ## Square In the "square" step of Diamond Square, we take a square region of the heightmap and set the center of it to the average of the four corners (plus jitter). <pre class="lineart"> Column 0 1 2 3 4 0 1 2 3 4 0 1 2 3 4 ┌─┬─┬─┬─┬─┐ ┌─┬─┬─┬─┬─┐ ┌─┬─┬─┬─┬─┐ 0 │1│ │ │ │8│ 0 │1│ │ │ │8│ 0 │1│ │ │ │8│ ├─┼─┼─┼─┼─┤ ├─╲─┼─┼─╱─┤ ├─╲─┼─┼─╱─┤ 1 │ │ │ │ │ │ 1 │ │╲│ │╱│ │ 1 │ │╲│ │╱│ │ R ├─┼─┼─┼─┼─┤ ├─┼─◢─◣─┼─┤ ├─┼─◢─◣─┼─┤ o 2 │ │ │◉│ │ │ 2 │ │ │◉│ │ │ 2 │ │ │3│ │ │ w ├─┼─┼─┼─┼─┤ ├─┼─◥─◤─┼─┤ ├─┼─◥─◤─┼─┤ 3 │ │ │ │ │ │ 3 │ │╱│ │╲│ │ 3 │ │╱│ │╲│ │ ├─┼─┼─┼─┼─┤ ├─╱─┼─┼─╲─┤ ├─╱─┼─┼─╲─┤ 4 │0│ │ │ │3│ 4 │0│ │ │ │3│ 4 │0│ │ │ │3│ └─┴─┴─┴─┴─┘ └─┴─┴─┴─┴─┘ └─┴─┴─┴─┴─┘ </pre> The Wikipedia article calls this the "diamond" step, but that seems backwards to me. We're operating on a *square* region of the heightmap, so this should be the *square* step. Since we're no longer working with heightmap "slices" like we were in the previous post, we'll need a way to specify the square chunk of the heightmap we're working on. I chose to use a center point (`x` and `y` coordinates) and a radius: ```clojure (defn ds-square [heightmap x y radius spread] (let [new-height (jitter (average4 (heightmap-get heightmap (- x radius) (- y radius)) (heightmap-get heightmap (- x radius) (+ y radius)) (heightmap-get heightmap (+ x radius) (- y radius)) (heightmap-get heightmap (+ x radius) (+ y radius))) spread)] (heightmap-set! heightmap x y new-height))) ``` <div id="demo-2" class="threejs"></div> ## Diamond In the "diamond" step of Diamond Square, we take a diamond-shaped region and set the center of it to the average of the four corners: <pre class="lineart"> Column 0 1 2 3 4 5 6 7 8 0 1 2 3 4 5 6 7 8 ┌─┬─┬─┬─┬─┬─┬─┬─┬─┐ ┌─┬─┬─┬─┬─┬─┬─┬─┬─┐ 0 │1│ │ │ │3│ │ │ │3│ 0 │1│ │ │ │3│ │ │ │3│ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ ├─┼─┼─┼─╱─╲─┼─┼─┼─┤ 1 │ │ │ │ │ │ │ │ │ │ 1 │ │ │ │╱│ │╲│ │ │ │ R ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ ├─┼─┼─╱─┼─┼─╲─┼─┼─┤ o 2 │ │ │3│ │ │ │4│ │ │ 2 │ │ │3│ │◉│ │4│ │ │ w ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ ├─┼─┼─╲─┼─┼─╱─┼─┼─┤ 3 │ │ │ │ │ │ │ │ │ │ 3 │ │ │ │╲│ │╱│ │ │ │ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ ├─┼─┼─┼─╲─╱─┼─┼─┼─┤ 4 │5│ │ │ │5│ │ │ │5│ 4 │5│ │ │ │5│ │ │ │5│ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ 5 │ │ │ │ │ │ │ │ │ │ 5 │ │ │ │ │ │ │ │ │ │ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ 6 │ │ │7│ │ │ │6│ │ │ 6 │ │ │7│ │ │ │6│ │ │ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ 7 │ │ │ │ │ │ │ │ │ │ 7 │ │ │ │ │ │ │ │ │ │ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ 8 │9│ │ │ │7│ │ │ │7│ 8 │9│ │ │ │7│ │ │ │7│ └─┴─┴─┴─┴─┴─┴─┴─┴─┘ └─┴─┴─┴─┴─┴─┴─┴─┴─┘ Column 0 1 2 3 4 5 6 7 8 0 1 2 3 4 5 6 7 8 ┌─┬─┬─┬─┬─┬─┬─┬─┬─┐ ┌─┬─┬─┬─┬─┬─┬─┬─┬─┐ 0 │1│ │ │ │3│ │ │ │3│ 0 │1│ │ │ │3│ │ │ │3│ ├─┼─┼─┼─╱┬╲─┼─┼─┼─┤ ├─┼─┼─┼─╱┬╲─┼─┼─┼─┤ 1 │ │ │ │╱│││╲│ │ │ │ 1 │ │ │ │╱│││╲│ │ │ │ R ├─┼─┼─╱─┼▼┼─╲─┼─┼─┤ ├─┼─┼─╱─┼▼┼─╲─┼─┼─┤ o 2 │ │ │3├─▶◉◀─┤4│ │ │ 2 │ │ │3├─▶4◀─┤4│ │ │ w ├─┼─┼─╲─┼▲┼─╱─┼─┼─┤ ├─┼─┼─╲─┼▲┼─╱─┼─┼─┤ 3 │ │ │ │╲│││╱│ │ │ │ 3 │ │ │ │╲│││╱│ │ │ │ ├─┼─┼─┼─╲┴╱─┼─┼─┼─┤ ├─┼─┼─┼─╲┴╱─┼─┼─┼─┤ 4 │5│ │ │ │5│ │ │ │5│ 4 │5│ │ │ │5│ │ │ │5│ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ 5 │ │ │ │ │ │ │ │ │ │ 5 │ │ │ │ │ │ │ │ │ │ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ 6 │ │ │7│ │ │ │6│ │ │ 6 │ │ │7│ │ │ │6│ │ │ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ 7 │ │ │ │ │ │ │ │ │ │ 7 │ │ │ │ │ │ │ │ │ │ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ 8 │9│ │ │ │7│ │ │ │7│ 8 │9│ │ │ │7│ │ │ │7│ └─┴─┴─┴─┴─┴─┴─┴─┴─┘ └─┴─┴─┴─┴─┴─┴─┴─┴─┘ </pre> Once again we'll use a center point and radius to specify the region inside the heightmap: ```clojure (defn ds-diamond [heightmap x y radius spread] (let [new-height (jitter (safe-average (heightmap-get-safe heightmap (- x radius) y) (heightmap-get-safe heightmap (+ x radius) y) (heightmap-get-safe heightmap x (- y radius)) (heightmap-get-safe heightmap x (+ y radius))) spread)] (heightmap-set! heightmap x y new-height))) ``` <div id="demo-3" class="threejs"></div> ## Edge Cases You might have noticed that we used `heightmap-get-safe` and `safe-average` in the `ds-diamond` code. This is to handle the diamonds on the edge of the heightmap. Those diamonds only have three points to average: <pre class="lineart"> Column 0 1 2 3 4 0 1 2 3 4 ┌─┬─┬─┬─┬─┐ ┌─┬─┬─┬─┬─┐ 0 │1│ │ │ │8│ 0 │1│ │ │ │8│ ├─┼─┼─┼─┼─┤ ├─┼─┼─┼─╱┬╲ 1 │ │ │ │ │ │ 1 │ │ │ │╱│││╲ R ├─┼─┼─┼─┼─┤ ├─┼─┼─╱─┼▼┤ ╲ o 2 │ │ │3│ │ │ 2 │ │ │3├─▶◉◀──? w ├─┼─┼─┼─┼─┤ ├─┼─┼─╲─┼▲┤ ╱ 3 │ │ │ │ │ │ 3 │ │ │ │╲│││╱ ├─┼─┼─┼─┼─┤ ├─┼─┼─┼─╲┴╱ 4 │0│ │ │ │3│ 4 │0│ │ │ │3│ └─┴─┴─┴─┴─┘ └─┴─┴─┴─┴─┘ </pre> So we just ignore the missing point: ```clojure (defn safe-average [a b c d] (let [total 0 count 0] (when a (add! total a) (inc! count)) (when b (add! total b) (inc! count)) (when c (add! total c) (inc! count)) (when d (add! total d) (inc! count)) (/ total count))) (defn heightmap-get-safe [heightmap x y] (let [last (heightmap-last-index heightmap)] (when (and (<= 0 x last) (<= 0 y last)) (heightmap-get heightmap x y)))) ``` Technically this means they have less information to work with than the rest of the points, and you might see some slight artifacts. But in practice it's not too bad, and it's a lot nicer than midpoint displacement where *every* line step uses only *two* points. Another way to handle this is to also wrap those edge coordinates around and pull the point from the other side of the heightmap. This is a bit more work but gives you an extra benefit: the heightmap will be tileable (**update**: like everything where someone doesn't show you running code, it's [not actually that simple](/blog/2016/08/lisp-jam-postmortem/#tiling-diamond-square)). ## Iteration Unfortunately we can't use recursion to iterate for Diamond Square like we did for Midpoint Displacement. We have to interleave the diamond and square steps, performing each round of squares on the *entire* map before moving on to a round of diamonds. If we don't, bad things will happen. Let's look at an example to see why. First we perform the initial square step: <pre class="lineart"> Column 0 1 2 3 4 0 1 2 3 4 ┌─┬─┬─┬─┬─┐ ┌─┬─┬─┬─┬─┐ 0 │1│ │ │ │8│ 0 │1│ │ │ │8│ ├─┼─┼─┼─┼─┤ ├─╲─┼─┼─╱─┤ 1 │ │ │ │ │ │ 1 │ │╲│ │╱│ │ R ├─┼─┼─┼─┼─┤ ├─┼─◢─◣─┼─┤ o 2 │ │ │ │ │ │ 2 │ │ │3│ │ │ w ├─┼─┼─┼─┼─┤ ├─┼─◥─◤─┼─┤ 3 │ │ │ │ │ │ 3 │ │╱│ │╲│ │ ├─┼─┼─┼─┼─┤ ├─╱─┼─┼─╲─┤ 4 │0│ │ │ │3│ 4 │0│ │ │ │3│ └─┴─┴─┴─┴─┘ └─┴─┴─┴─┴─┘ </pre> That's fine. Then we recurse into the first of the four diamonds: <pre class="lineart"> Column 0 1 2 3 4 0 1 2 3 4 ┌─┬─┬─┬─┬─┐ ┌─┬─┬─┬─┬─┐ 0 │1│ │ │ │8│ 0 │1├─▶4◀─┤8│ ├─┼─┼─┼─┼─┤ ├─╲─┼▲┼─╱─┤ 1 │ │ │ │ │ │ 1 │ │╲│││╱│ │ R ├─┼─┼─┼─┼─┤ ├─┼─╲┴╱─┼─┤ o 2 │ │ │3│ │ │ 2 │ │ │3│ │ │ w ├─┼─┼─┼─┼─┤ ├─┼─┼─┼─┼─┤ 3 │ │ │ │ │ │ 3 │ │ │ │ │ │ ├─┼─┼─┼─┼─┤ ├─┼─┼─┼─┼─┤ 4 │0│ │ │ │3│ 4 │0│ │ │ │3│ └─┴─┴─┴─┴─┘ └─┴─┴─┴─┴─┘ </pre> Everything's still good. Then we recurse into the square: <pre class="lineart"> Column 0 1 2 3 4 0 1 2 3 4 ┌─┬─┬─┬─┬─┐ ╔═════╗─┬─┐ 0 │1│ │4│ │8│ 0 ║1│ │4║ │8│ ├─┼─┼─┼─┼─┤ ║─◢─◣─║─┼─┤ 1 │ │ │ │ │ │ 1 ║ │ │ ║ │ │ R ├─┼─┼─┼─┼─┤ ║─◥─◤─║─┼─┤ o 2 │ │ │3│ │ │ 2 ║?│ │3║ │ │ w ├─┼─┼─┼─┼─┤ ╚═════╝─┼─┤ 3 │ │ │ │ │ │ 3 │ │ │ │ │ │ ├─┼─┼─┼─┼─┤ ├─┼─┼─┼─┼─┤ 4 │0│ │ │ │3│ 4 │0│ │ │ │3│ └─┴─┴─┴─┴─┘ └─┴─┴─┴─┴─┘ </pre> But wait, there's a corner missing! We needed to complete that left-hand diamond step to get it, but that's still waiting on the stack for this line of recursion to bottom out. ### A "Strided" Iteration Construct So we need to use normal iteration. We can start by writing the code to do a single round of diamonds or squares on the heightmap. First we'll make a nice looping construct called `do-stride`: ```clojure (defmacro do-stride [varnames start-form end-form stride-form & body] (let [stride (gensym "stride") start (gensym "start") end (gensym "end") build (fn build [vars] (if (empty? vars) `(do ~@body) (let [varname (first vars)] `(loop [~varname ~start] (when (< ~varname ~end) ~(build (rest vars)) (recur (+ ~varname ~stride)))))))] ; Fix the numbers once outside the nested loops, ; and then build the guts. `(let [~start ~start-form ~end ~end-form ~stride ~stride-form] ~(build varnames)))) ``` This scary-looking macro builds the nested looping structure we'll need to perform the iteration for our diamonds and squares. It abstracts away the tedium of writing out nested `for` loops so we can focus on the rest of the algorithm. `do-stride` takes a list of variables to bind, a number to start at, a number to end before, and a body. For each variable it starts at `start`, runs `body`, increments by `stride`, and loops until the value is at or above `end`. For example: * `(do-stride [x] 0 5 1 (console.log x))` will log `0 1 2 3 4`. * `(do-stride [x] 0 5 2 (console.log x))` will log `0 2 4`. * `(do-stride [x y] 0 3 2 (console.log [x y]))` will log `[0, 0] [0, 2] [2, 0] [2, 2]`. ### Square Iteration Now that we've got this, we can iterate our square step over the heightmap for a specific radius: ```clojure (defn ds-squares [heightmap radius spread] (do-stride [x y] radius (heightmap-resolution heightmap) (* 2 radius) (ds-square heightmap x y radius spread))) ``` We use `do-stride` to start at the radius and loop until we fall off the heightmap, stepping by *twice* the radius each time. Because our heightmaps are always square we can use the same values for both dimensions. `do-stride` can take zero or more variables and will just do the right thing, so we don't need to worry about it. The stride is double the radius because as we step between squares we move over the right half of the first *and* the left half of the next: <pre class="lineart"> Column 0 1 2 3 4 5 6 7 8 0 1 2 3 4 5 6 7 8 ╔═════════╗─┬─┬─┬─┐ ╔═════════╗─┬─┬─┬─┐ 0 ║1│ │ │ │3║ │ │ │3│ 0 ║1│ │ │ │3║ │ │ │3│ ║─┼─┼─┼─┼─║─┼─┼─┼─┤ ║─┼─┼─┼─┼─║─┼─┼─┼─┤ 1 ║ │ │ │ │ ║ │ │ │ │ 1 ║ │ │ │ │ ║ │ │ │ │ R ║─┼─┼─┼─┼─║─┼─┼─┼─┤ radius ┼─┼─║─┼─┼─┼─┤ o 2 ║ │ │◉│ │ ║ │◉│ │ │ 2 ║◀━━▶◉◀━━━━━▶◉│ │ │ w ║─┼─┼─┼─┼─║─┼─┼─┼─┤ ║─┼─┼─ stride ┼─┼─┤ 3 ║ │ │ │ │ ║ │ │ │ │ 3 ║ │ │ │ │ ║ │ │ │ │ ║─┼─┼─┼─┼─║─┼─┼─┼─┤ ║─┼─┼─┼─┼─║─┼─┼─┼─┤ 4 ║5│ │ │ │5║ │ │ │5│ 4 ║5│ │ │ │5║ │ │ │5│ ╚═════════╝─┼─┼─┼─┤ ╚═════════╝─┼─┼─┼─┤ 5 │ │ │ │ │ │ │ │ │ │ 5 │ │ │ │ │ │ │ │ │ │ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ 6 │ │ │◉│ │ │ │◉│ │ │ 6 │ │ │◉│ │ │ │◉│ │ │ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ 7 │ │ │ │ │ │ │ │ │ │ 7 │ │ │ │ │ │ │ │ │ │ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ 8 │9│ │ │ │7│ │ │ │7│ 8 │9│ │ │ │7│ │ │ │7│ └─┴─┴─┴─┴─┴─┴─┴─┴─┘ └─┴─┴─┴─┴─┴─┴─┴─┴─┘ </pre> ### Diamond Iteration The diamonds are a little more complicated. First of all, in the `y` dimension we need to start at `0` instead of `radius` (these will be the top and bottom edge cases we looked at). Also, in the `x` direction all the even iterations need to start at `radius`, while the odd iterations should start at `0`. Hopefully a picture will make it a bit more clear: <pre class="lineart"> Column 0 1 2 3 4 5 6 7 8 ┌─┬─┬─┬─┬─┬─┬─┬─┬─┐ 0 │1│ │◉│ │3│ │◉│ │3│ y iteration: 0, start at x = radius ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ 1 │ │ │ │ │ │ │ │ │ │ R ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ o 2 │◉│ │3│ │◉│ │4│ │◉│ y iteration: 1, start at x = 0 w ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ 3 │ │ │ │ │ │ │ │ │ │ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ 4 │5│ │◉│ │5│ │◉│ │5│ y iteration: 2, start at x = radius ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ 5 │ │ │ │ │ │ │ │ │ │ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ 6 │◉│ │7│ │◉│ │6│ │◉│ y iteration: 3, start at x = 0 ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ 7 │ │ │ │ │ │ │ │ │ │ ├─┼─┼─┼─┼─┼─┼─┼─┼─┤ 8 │9│ │◉│ │7│ │◉│ │7│ y iteration: 4, start at x = radius └─┴─┴─┴─┴─┴─┴─┴─┴─┘ </pre> The code for the diamonds is the hairiest bit of the algorithm. Make sure you understand it before moving on. `shift` will be `radius` for even iterations and `0` for odd ones: ```clojure (defn ds-diamonds [heightmap radius spread] (let [size (heightmap-resolution heightmap)] (do-stride [y] 0 size radius (let [shift (if (even? (/ y radius)) radius 0)] (do-stride [x] shift size (* 2 radius) (ds-diamond heightmap x y radius spread)))))) ``` ### Top-Level Iteration Finally we can put it all together. We initialize the corners and starting values, then loop over smaller and smaller radii and just call `ds-squares` and `ds-diamonds` to do the heavy lifting: ```clojure (defn diamond-square [heightmap] (let [initial-spread 0.3 spread-reduction 0.5 center (heightmap-center-index heightmap) size (aget heightmap.shape 0)] (ds-init-corners heightmap) (loop [radius center spread initial-spread] (when (>= radius 1) (ds-squares heightmap radius spread) (ds-diamonds heightmap radius spread) (recur (/ radius 2) (* spread spread-reduction)))) (sanitize heightmap))) ``` ## Result The end result looks like this: <div id="demo-final" class="threejs"></div> It looks a lot like the result from Midpoint Displacement, but without the ugly seams around the square chunks. The code for these blog posts is a bit of a mess because I've been copy/pasting to show the partially-completed demos after each step. I've created a little [single-page demo][ymir] with completed versions of the various algorithms you can play with, and [that code][ymir-code] should be a lot more readable than the hacky code for these posts. [ymir]: http://ymir.stevelosh.com/ [ymir-code]: http://bitbucket.org/sjl/ymir/