# Terrain Generation with Midpoint Displacement

Posted on February 19th, 2016.

I'm taking the Game Engine Architecture class at Reykjavík University. My group just finished our midterm project where we played around with procedural terrain generation in Unity. It was a lot of fun and I really enjoyed creating terrain out of numbers, so I thought I'd write up a little introduction to a few of the algorithms.

In this post we'll be looking at Midpoint Displacement.

The full series of posts so far:

## Terrain Generation

Procedural generation is a pretty big field. The Wikipedia article gives a pretty brief overview. The Procedural Content Generation Wiki is an entire wiki devoted to procedural generation. There's also a small but active subreddit for interested folks.

Today we're just going to look at a single algorithm for generating realistic-looking terrain called "Midpoint Displacement". It's relatively simple (compared to other algorithms) but it produces terrain that actually looks kind of cool.

## Resources, Code, and Examples

There are a lot of resources for learning how terrain generation algorithms work. For midpoint displacement I used the Wikipedia article, the PCG Wiki article, and academic papers like this one.

Unfortunately, while there are a *lot* of places that describe algorithms like
this, there are relatively few that explain it *thoroughly*. Academic papers in
particular tend to have a really bad case of
draw-the-rest-of-the-fucking-owl-syndrome. They tend to slap an equation
or two on the page and move on without showing any code.

As anyone who's done much programming knows, there's a big difference between
saying something like "set the midpoints of the diamonds to the average of the
corners" and making a computer actually *do* that. How do we iterate over the
array (row-major or column-major)? Does it matter (yes: if you care about cache
coherency)? What about the edges where there are fewer corners? How random
(and how fast) is our random number generator (and how much do we care)? How do
we write all this code so we can actually read it six months later?

So I'm going to try to do my part to improve the situation by not just
describing the algorithm, but showing *actual code* that runs and generates
a terrain.

The full code is here, but we'll be going through the important parts below.

### Wisp

I'm using a language called Wisp. It's a small, Clojure-inspired dialect of Lisp that compiles down to vanilla Javascript.

I went with Javascript because it means I can actually put demos inline in this post. Sometimes it's so much easier to understand something when you can actually see it running.

I've tried to write the code so that it's readable even if you don't know Lisp.
The point of me showing you this code is not to give you something to copy and
paste. The point is that until you *actually write and run* the code, you don't
know all the edge cases and problems that can happen. So by making examples
that actually run I feel a bit more confident that I'm explaining things enough.

I haven't focused on performance at all. If you want a fast implementation of this, you'll probably want to use a language with real arrays (i.e. contiguous chunks of memory holding a bunch of unboxed floats or integers).

### Three.js

Three.js is a Javascript library for 3D rendering. I'm using it
here to get something on the screen you can play with. I'm not going to cover
the code to put the terrains on the screen because it's really an entirely
separate problem to *generating* the terrain, and it's mostly uninteresting
boilerplate for these little demos.

I haven't been able to get these demos working on my iPhone. Sorry about that, but I'm not much of a frontend developer. You'll have to use a computer to see them. Pull requests are welcome.

## Heightmaps

Let's get started. The main data structure we're going to be using is a heightmap.

### Overview

Heightmaps are essentially just big two-dimensional arrays of numbers in a given range that represent the height of the ground at a given point. For example, a heightmap with a small "mountain" in the middle might look like:

```
[[0, 0, 2, 0, 0],
[0, 2, 5, 1, 0],
[2, 5, 9, 4, 1],
[0, 1, 4, 1, 0],
[0, 0, 1, 0, 0]]
```

In practice, heightmaps are often represented as a single-dimensional array instead to ensure that all the data is together in one big hunk of memory:

```
[0, 0, 2, 0, 0,
0, 2, 5, 1, 0,
2, 5, 9, 4, 1,
0, 1, 4, 1, 0,
0, 0, 1, 0, 0]
```

The type and range of the numbers in the map varies by program and programmer. Some programs use nonzero integers, some use all the integers, some use floating point, etc.

The heightmaps in this post will contain floats from `0.0`

to `1.0`

. I like
using floats between zero and one because it's easy to roughly visualize them in
your head (e.g. `0.1`

is ten percent of the maximum height).

One common way to visualize heightmaps is by turning them into greyscale images
with one pixel per number, with black pixels for lowest value in the map's range
and white for the highest. The Wikipedia article has an example of
this. This is handy because you can also write some code to *read* images, and
then you can use any image editor to draw your own heightmaps by hand.

We're going to visualize the heightmaps in 3D. It's a lot more fun, and they start to actually look like real terrains.

### Code

We're going to start by writing a few functions that will hide away the internal details of how our heightmaps are implemented, so we can talk about them in nice high-level terms for the rest of the program. We'll also need a couple of helper functions along the way.

We'll start with something to create a heightmap:

`(`*defn* make-heightmap [exponent]
(*let* [resolution (+ 1 (Math.pow 2 exponent))]
(l (+ "Creating " resolution " by " resolution " heightmap..."))
(def heightmap
(new Array (* resolution resolution)))
(set! heightmap.resolution resolution)
(set! heightmap.exponent exponent)
(set! heightmap.last (- resolution 1))
(zero-heightmap heightmap)))

`make-heightmap`

takes an exponent and creates a heightmap. For reasons that
will become clear later, all of our heightmaps must be square, and they'll need
to have \(2^n + 1\) rows and columns. This means we can make heightmaps of
sizes like 3x3, 5x5, 9x9, 17x17, 33x33, etc. So our function just takes the
\(n\) to use and creates a correctly-sized array.

(If you've ever poked around with Unity's Terrain objects you might recognize these "powers of two plus one".)

We create a new array (1-dimensional) and store a few pieces of extra data on it
for later use. `last`

is the index of the last row/column, so for a 3x3 array
it would be `2`

.

`l`

is just a simple little logging function:

`(`*defn* l [v]
(console.log v))

We'll need to zero out the heightmap as well:

`(`*defn* zero-heightmap [heightmap]
(do-times i heightmap.length
(set! (aget heightmap i) 0.0))
heightmap)

`zero-heightmap`

just iterates from `0`

to `heightmap.length`

and sets each
element to `0.0`

. Wisp doesn't have syntax to loop from 0 to some number, but
it's a Lisp, so we can make some:

`(`*defmacro* when [condition & body]
`(*if* ~condition
(do ~@body)))
(*defmacro* do-times [varname limit & body]
(*let* [end (gensym)]
`(*let* [~end ~limit]
(*loop* [~varname 0]
(when (< ~varname ~end)
~@body
(recur (+ 1 ~varname)))))))

If you're not a Lisp person, don't sweat this too much. Just trust me that
`(do-times i 10 (console.log i))`

will print the numbers `0`

to `9`

.

Now that we can create a heightmap, we'll probably want to be able to read its elements:

`(`*defmacro* heightmap-get [hm x y]
`(aget ~hm (+ (* ~y (.-resolution ~hm)) ~x)))
(*defn* heightmap-get-safe [hm x y]
(when (and (<= 0 x hm.last)
(<= 0 y hm.last))
(heightmap-get hm x y)))

`heightmap-get`

handles the nasty business of translating our human-thinkable
`x`

and `y`

coordinates into an index into the single-dimensional array. It's
not particularly fast (especially if the compiler won't inline it), but we're
not worried about speed for this post.

`heightmap-get-safe`

is a version of `heightmap-get`

that will do bounds
checking for us and return `nil`

if we ask for something out of range.
Javascript will return `undefined`

if you ask for a non-existent index, but
because we're storing the array as one big line of data, using a `y`

value
that's too big will actually "wrap around" to the next row, so it's safer to
just check it explicitly.

Of course we'll also want to be able to set values in our heightmap:

`(`*defmacro* heightmap-set! [hm x y val]
`(set! (heightmap-get ~hm ~x ~y) ~val))

And finally, we'll want to be able to "normalize" our array. During the course
of running our algorithm, it's possible for the values to become less than `0.0`

or greater than `1.0`

. We want our heightmaps to contain only values between
zero and one, so rather than change the algorithm we can just loop through the
array at the end and "squeeze" the range down to 0 and 1 (or "stretch" it if
it ends up smaller):

`(`*defn* normalize [hm]
(*let* [max (- Infinity)
min Infinity]
(do-times i hm.length
(*let* [el (aget hm i)]
(when (< max el) (set! max el))
(when (> min el) (set! min el))))
(*let* [span (- max min)]
(do-times i hm.length
(set! (aget hm i)
(/ (- (aget hm i) min)
span))))))

## Random Noise

Now that we've got a nice little interface to heightmaps we can move on to actually making some terrain! Before we dive into midpoint displacement, let's get a black triangle on the screen. We'll start by just setting every point in the heightmap to a completely random value.

`(`*defn* rand []
(Math.random))
(*defn* rand-around-zero [spread]
(- (* spread (rand) 2) spread))
(*defn* random-noise [heightmap]
(do-times i heightmap.length
(set! (aget heightmap i) (rand))))

Let's see it! Hit the "refresh" button to create and render the terrain. You can keep hitting it to get fresh terrains.

Not really a very fun landscape to explore, but at least we've got something running. Onward to the real algorithm.

## Midpoint Displacement

Midpoint displacement is a simple, fast way to generate decent-looking terrain. It's not perfect (you'll notice patterns once you start viewing the mesh from different angles) but it's a good place to start.

The general process goes like this:

- Initialize the four corners of the heightmap to random values.
- Set the midpoints of each edge to the average of the two corners it's between, plus or minus a random amount.
- Set the center of the square to the average of those edge midpoints you just set, again with a random jitter.
- Recurse on the four squares inside this one, reducing the jitter.

That's the standard description. Now let's draw the fucking owl.

### Initialize the Corners

The first step is pretty easy: we just need to take our heightmap:

Column 0 1 2 3 4 ┌───┬───┬───┬───┬───┐ 0 │ │ │ │ │ │ ├───┼───┼───┼───┼───┤ 1 │ │ │ │ │ │ R ├───┼───┼───┼───┼───┤ o 2 │ │ │ │ │ │ w ├───┼───┼───┼───┼───┤ 3 │ │ │ │ │ │ ├───┼───┼───┼───┼───┤ 4 │ │ │ │ │ │ └───┴───┴───┴───┴───┘

Find the corners:

Column 0 1 2 3 4 ┌───┬───┬───┬───┬───┐ 0 │ ◉ │ │ │ │ ◉ │ ├───┼───┼───┼───┼───┤ 1 │ │ │ │ │ │ R ├───┼───┼───┼───┼───┤ o 2 │ │ │ │ │ │ w ├───┼───┼───┼───┼───┤ 3 │ │ │ │ │ │ ├───┼───┼───┼───┼───┤ 4 │ ◉ │ │ │ │ ◉ │ └───┴───┴───┴───┴───┘

And shove random values in them:

Column 0 1 2 3 4 ┌───┬───┬───┬───┬───┐ 0 │ 2 │ │ │ │ 8 │ ├───┼───┼───┼───┼───┤ 1 │ │ │ │ │ │ R ├───┼───┼───┼───┼───┤ o 2 │ │ │ │ │ │ w ├───┼───┼───┼───┼───┤ 3 │ │ │ │ │ │ ├───┼───┼───┼───┼───┤ 4 │ 0 │ │ │ │ 3 │ └───┴───┴───┴───┴───┘

(I'm using 0 to 10 instead of 0 to 1 because it's easier to read in the ASCII-art diagrams.)

The code is about what you'd expect:

`(`*defn* mpd-init-corners [heightmap]
(heightmap-set! heightmap 0 0 (rand))
(heightmap-set! heightmap 0 heightmap.last (rand))
(heightmap-set! heightmap heightmap.last 0 (rand))
(heightmap-set! heightmap heightmap.last heightmap.last (rand)))
(*defn* midpoint-displacement-d1 [heightmap]
(mpd-init-corners heightmap))

That was easy.

### Set the Edges

For the next step, we need to take our square heightmap with filled-in corners and use the corners to fill in the middle element of each edge:

Column 0 1 2 3 4 ┌───┬───┬───┬───┬───┐ 0 │ 2 │ │ ◉ │ │ 8 │ ├───┼───┼───┼───┼───┤ 1 │ │ │ │ │ │ R ├───┼───┼───┼───┼───┤ o 2 │ ◉ │ │ │ │ ◉ │ w ├───┼───┼───┼───┼───┤ 3 │ │ │ │ │ │ ├───┼───┼───┼───┼───┤ 4 │ 0 │ │ ◉ │ │ 3 │ └───┴───┴───┴───┴───┘ Column 0 1 2 3 4 ┌───┬───┬───┬───┬───┐ 0 │ 2 ├───▶ ◉ ◀───┤ 8 │ ├─┬─┼───┼───┼───┼─┬─┤ 1 │ │ │ │ │ │ │ │ R ├─▼─┼───┼───┼───┼─▼─┤ o 2 │ ◉ │ │ │ │ ◉ │ w ├─▲─┼───┼───┼───┼─▲─┤ 3 │ │ │ │ │ │ │ │ ├─┴─┼───┼───┼───┼─┴─┤ 4 │ 0 ├───▶ ◉ ◀───┤ 3 │ └───┴───┴───┴───┴───┘ Column 0 1 2 3 4 ┌───┬───┬───┬───┬───┐ 0 │ 2 ├───▶ 5 ◀───┤ 8 │ ├─┬─┼───┼───┼───┼─┬─┤ 1 │ │ │ │ │ │ │ │ R ├─▼─┼───┼───┼───┼─▼─┤ o 2 │ 1 │ │ │ │5.5│ w ├─▲─┼───┼───┼───┼─▲─┤ 3 │ │ │ │ │ │ │ │ ├─┴─┼───┼───┼───┼─┴─┤ 4 │ 0 ├───▶1.5◀───┤ 3 │ └───┴───┴───┴───┴───┘

This is why the size of the heightmap has to be \(2^n + 1\) — the \(+ 1\)
ensures that the sides are odd lengths, which makes sure they *have* a middle
element for us to set.

Let's make a couple of helper functions:

`(`*defn* jitter [value spread]
(+ value (rand-around-zero spread)))
(*defn* midpoint [a b]
(/ (+ a b) 2))
(*defn* average2 [a b]
(/ (+ a b) 2))
(*defn* average4 [a b c d]
(/ (+ a b c d) 4))

`jitter`

is just a nice way of saying "take this point and move it a random
amount up or down".

The other three functions are pretty self-explanatory. Also, two of them are exactly the same function! There's a reason for this.

Lately I've started digging into SICP again and watching the lectures. Abelson and Sussman practice a style of programming that's more about "growing" a language to talk about a problem you're trying to solve. So you start with the base programming language, and then define new terms/syntax that let you talk about your problem at the appropriate level of detail.

So in this case, we've got two functions for conceptually different tasks:

`midpoint`

finds the middle point between two points (indexes in an array in this case).`average2`

finds the average of two values (floating point height values here).

They happen to be implemented in the same way, but I think giving them different names makes them easier to talk about and easier to keep straight in your head.

`(`*defn* mpd-displace [heightmap lx rx by ty spread]
(*let* [cx (midpoint lx rx)
cy (midpoint by ty)
bottom-left (heightmap-get heightmap lx by)
bottom-right (heightmap-get heightmap rx by)
top-left (heightmap-get heightmap lx ty)
top-right (heightmap-get heightmap rx ty)
top (average2 top-left top-right)
left (average2 bottom-left top-left)
bottom (average2 bottom-left bottom-right)
right (average2 bottom-right top-right)]
(heightmap-set! heightmap cx by (jitter bottom spread))
(heightmap-set! heightmap cx ty (jitter top spread))
(heightmap-set! heightmap lx cy (jitter left spread))
(heightmap-set! heightmap rx cy (jitter right spread))))
(*defn* midpoint-displacement-d2 [heightmap]
(mpd-init-corners heightmap)
(mpd-displace heightmap
0 heightmap.last
0 heightmap.last
0.1))

`mpd-displace`

is the meat of the algorithm. It takes a heightmap, the spread
(how much random "jitter" to apply after averaging the values), and four values
to define the square we're working on in the heightmap:

`lx`

: the x coordinate for the left-hand corners`rx`

: the x coordinate for the right-hand corners`by`

: the y coordinate for the bottom corners`ty`

: the y coordinate for the top corners

For the first iteration we want to work on the entire heightmap, so we pass in
`0`

and `heightmap.last`

for the left/right and top/bottom indices.

Run the demo a few times and watch how the midpoints fall between the corner points because they're averaged (with a little bit of jitter thrown in).

### Set the Center

The next step is to set the center element of the square to the average of those midpoints we just set, plus a bit of jitter:

Column 0 1 2 3 4 ┌───┬───┬───┬───┬───┐ 0 │ 2 │ │ 5 │ │ 8 │ ├───┼───┼───┼───┼───┤ 1 │ │ │ │ │ │ R ├───┼───┼───┼───┼───┤ o 2 │ 1 │ │ ◉ │ │5.5│ w ├───┼───┼───┼───┼───┤ 3 │ │ │ │ │ │ ├───┼───┼───┼───┼───┤ 4 │ 0 │ │1.5│ │ 3 │ └───┴───┴───┴───┴───┘ Column 0 1 2 3 4 ┌───┬───┬───┬───┬───┐ 0 │ 2 │ │ 5 │ │ 8 │ ├───┼───┼─┬─┼───┼───┤ 1 │ │ │ │ │ │ │ R ├───┼───┼─▼─┼───┼───┤ o 2 │ 1 ├───▶ ◉ ◀───┤5.5│ w ├───┼───┼─▲─┼───┼───┤ 3 │ │ │ │ │ │ │ ├───┼───┼─┴─┼───┼───┤ 4 │ 0 │ │1.5│ │ 3 │ └───┴───┴───┴───┴───┘ Column 0 1 2 3 4 ┌───┬───┬───┬───┬───┐ 0 │ 2 │ │ 5 │ │ 8 │ ├───┼───┼─┬─┼───┼───┤ 1 │ │ │ │ │ │ │ R ├───┼───┼─▼─┼───┼───┤ o 2 │ 1 ├───▶3.3◀───┤5.5│ w ├───┼───┼─▲─┼───┼───┤ 3 │ │ │ │ │ │ │ ├───┼───┼─┴─┼───┼───┤ 4 │ 0 │ │1.5│ │ 3 │ └───┴───┴───┴───┴───┘

We can add this to the displacement function pretty easily:

`(`*defn* mpd-displace [heightmap lx rx by ty spread]
(*let* [cx (midpoint lx rx)
cy (midpoint by ty)
bottom-left (heightmap-get heightmap lx by)
bottom-right (heightmap-get heightmap rx by)
top-left (heightmap-get heightmap lx ty)
top-right (heightmap-get heightmap rx ty)
top (average2 top-left top-right)
left (average2 bottom-left top-left)
bottom (average2 bottom-left bottom-right)
right (average2 bottom-right top-right)
center (average4 top left bottom right)]
(heightmap-set! heightmap cx by (jitter bottom spread))
(heightmap-set! heightmap cx ty (jitter top spread))
(heightmap-set! heightmap lx cy (jitter left spread))
(heightmap-set! heightmap rx cy (jitter right spread))
(heightmap-set! heightmap cx cy (jitter center spread))))

And now our center is ready to go:

### Iterate

To finish things off, we need to iterate so that we can fill in the rest of the values in progressively smaller squares.

A 3x3 square is the base case, because its corners will be filled in to start, and the edges and midpoint get filled in. We just call our displacement function with the top/bottom/left/right values defining the square and we're done.

A 5x5 starts with one iteration that displaces the whole map. Then it does another "pass" that displaces each of the 3x3 "sub-squares":

╔═══════════╗───┬───┐ ┌───┬───╔═══════════╗ ║ 2 │ │ 5 ║ │ 8 │ │ 2 │ ║ 5 │ │ 8 ║ ║───┼───┼───║───┼───┤ ├───┼───║───┼───┼───║ ║ │ │ ║ │ │ │ │ ║ │ │ ║ ║───┼───┼───║───┼───┤ ├───┼───║───┼───┼───║ ║ 1 │ │3.3║ │5.5│ │ 1 │ ║3.3│ │5.5║ ╚═══════════╝───┼───┤ ├───┼───╚═══════════╝ │ │ │ │ │ │ │ │ │ │ │ │ ├───┼───┼───┼───┼───┤ ├───┼───┼───┼───┼───┤ │ 0 │ │1.5│ │ 3 │ │ 0 │ │1.5│ │ 3 │ └───┴───┴───┴───┴───┘ └───┴───┴───┴───┴───┘ ┌───┬───┬───┬───┬───┐ ┌───┬───┬───┬───┬───┐ │ 2 │ │ 5 │ │ 8 │ │ 2 │ │ 5 │ │ 8 │ ├───┼───┼───┼───┼───┤ ├───┼───┼───┼───┼───┤ │ │ │ │ │ │ │ │ │ │ │ │ ╔═══════════╗───┼───┤ ├───┼───╔═══════════╗ ║ 1 │ │3.3║ │5.5│ │ 1 │ ║3.3│ │5.5║ ║───┼───┼───║───┼───┤ ├───┼───║───┼───┼───║ ║ │ │ ║ │ │ │ │ ║ │ │ ║ ║───┼───┼───║───┼───┤ ├───┼───║───┼───┼───║ ║ 0 │ │1.5║ │ 3 │ │ 0 │ ║1.5│ │ 3 ║ ╚═══════════╝───┴───┘ └───┴───╚═══════════╝

If we have a bigger heightmap (e.g. 9x9) we'll need three passes:

Pass 1 ▩▩▩▩▩▩▩▩▩ ▩▩▩▩▩▩▩▩▩ ▩▩▩▩▩▩▩▩▩ ▩▩▩▩▩▩▩▩▩ ▩▩▩▩▩▩▩▩▩ ▩▩▩▩▩▩▩▩▩ ▩▩▩▩▩▩▩▩▩ ▩▩▩▩▩▩▩▩▩ ▩▩▩▩▩▩▩▩▩ Pass 2 ▩▩▩▩▩▢▢▢▢ ▢▢▢▢▩▩▩▩▩ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▩▩▩▩▩▢▢▢▢ ▢▢▢▢▩▩▩▩▩ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▩▩▩▩▩▢▢▢▢ ▢▢▢▢▩▩▩▩▩ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▩▩▩▩▩▢▢▢▢ ▢▢▢▢▩▩▩▩▩ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▩▩▩▩▩▢▢▢▢ ▢▢▢▢▩▩▩▩▩ ▩▩▩▩▩▢▢▢▢ ▢▢▢▢▩▩▩▩▩ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▩▩▩▩▩▢▢▢▢ ▢▢▢▢▩▩▩▩▩ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▩▩▩▩▩▢▢▢▢ ▢▢▢▢▩▩▩▩▩ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▩▩▩▩▩▢▢▢▢ ▢▢▢▢▩▩▩▩▩ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▩▩▩▩▩▢▢▢▢ ▢▢▢▢▩▩▩▩▩ Pass 3 ▩▩▩▢▢▢▢▢▢ ▢▢▩▩▩▢▢▢▢ ▢▢▢▢▩▩▩▢▢ ▢▢▢▢▢▢▩▩▩ ▩▩▩▢▢▢▢▢▢ ▢▢▩▩▩▢▢▢▢ ▢▢▢▢▩▩▩▢▢ ▢▢▢▢▢▢▩▩▩ ▩▩▩▢▢▢▢▢▢ ▢▢▩▩▩▢▢▢▢ ▢▢▢▢▩▩▩▢▢ ▢▢▢▢▢▢▩▩▩ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▩▩▩▢▢▢▢▢▢ ▢▢▩▩▩▢▢▢▢ ▢▢▢▢▩▩▩▢▢ ▢▢▢▢▢▢▩▩▩ ▩▩▩▢▢▢▢▢▢ ▢▢▩▩▩▢▢▢▢ ▢▢▢▢▩▩▩▢▢ ▢▢▢▢▢▢▩▩▩ ▩▩▩▢▢▢▢▢▢ ▢▢▩▩▩▢▢▢▢ ▢▢▢▢▩▩▩▢▢ ▢▢▢▢▢▢▩▩▩ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▩▩▩▢▢▢▢▢▢ ▢▢▩▩▩▢▢▢▢ ▢▢▢▢▩▩▩▢▢ ▢▢▢▢▢▢▩▩▩ ▩▩▩▢▢▢▢▢▢ ▢▢▩▩▩▢▢▢▢ ▢▢▢▢▩▩▩▢▢ ▢▢▢▢▢▢▩▩▩ ▩▩▩▢▢▢▢▢▢ ▢▢▩▩▩▢▢▢▢ ▢▢▢▢▩▩▩▢▢ ▢▢▢▢▢▢▩▩▩ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▢▢▢▢▢▢▢▢▢ ▩▩▩▢▢▢▢▢▢ ▢▢▩▩▩▢▢▢▢ ▢▢▢▢▩▩▩▢▢ ▢▢▢▢▢▢▩▩▩ ▩▩▩▢▢▢▢▢▢ ▢▢▩▩▩▢▢▢▢ ▢▢▢▢▩▩▩▢▢ ▢▢▢▢▢▢▩▩▩ ▩▩▩▢▢▢▢▢▢ ▢▢▩▩▩▢▢▢▢ ▢▢▢▢▩▩▩▢▢ ▢▢▢▢▢▢▩▩▩

Let's think about what's happening here.

- On pass
`0`

, we break up the grid into a 1x1 set of chunks and displace them (well, "it"). - On pass
`1`

, we break up the grid into a 2x2 set of chunks and displace them. - On pass
`2`

, we break up the grid into a 4x4 set of chunks and displace them.

This 9x9 heightmap was made with an exponent of 3 (\(2^3 + 1 = 9\)) and we need 3 "passes" to finish it. This isn't a coincidence — it's the reason we chose our heightmap resolutions to be \(2^n + 1\).

So at iteration \(i\), we need to displace a \(2^i\)x\(2^i\) collection of squares. Let's go:

`(`*defn* midpoint-displacement [heightmap]
(mpd-init-corners heightmap)
(*loop* [iter 0
spread 0.3]
(when (< iter heightmap.exponent)
(*let* [chunks (Math.pow 2 iter)
chunk-width (/ (- heightmap.resolution 1) chunks)]
(do-nested xchunk ychunk chunks
(*let* [left-x (* chunk-width xchunk)
right-x (+ left-x chunk-width)
bottom-y (* chunk-width ychunk)
top-y (+ bottom-y chunk-width)]
(mpd-displace heightmap
left-x right-x bottom-y top-y
spread))))
(recur (+ 1 iter) (* spread 0.5))))
(normalize heightmap))

That's a lot to take in. Let's break it down:

`(`*defn* midpoint-displacement [heightmap]
(mpd-init-corners heightmap)
(*loop* [iter 0
spread 0.3]
(when (< iter heightmap.exponent)
; Process the chunks at this iteration
(recur (+ 1 iter) (* spread 0.5))))
(normalize heightmap))

We initialize the corners at the beginning and normalize the heightmap at the
end. Then we loop a number of times equal to the exponent we used to make the
heightmap, as described above. Each time through the loop we increment `iter`

and cut the `spread`

for those points in half.

For each pass, we figure out what we're going to need to do:

`(`*let* [chunks (Math.pow 2 iter)
chunk-width (/ (- heightmap.resolution 1) chunks)]
; Actually do it...
)

Let's make a quick "nested for loop" syntax. It'll make things easier to read in the main function and we'll need it for later algorithms too:

`(`*defmacro* do-nested [xname yname width & body]
(*let* [iterations (gensym)]
`(*let* [~iterations ~width]
(do-times ~xname ~iterations
(do-times ~yname ~iterations
~@body)))))

This lets us say something like `(do-nested x y 5 (console.log [x y]))`

to mean:

```
for (x = 0; x < 5; x++) {
for (y = 0; y < 5; y++) {
console.log([x, y]);
}
}
```

And now we can see how we process the chunks:

`(`*let* [chunks (Math.pow 2 iter)
chunk-width (/ (- heightmap.resolution 1) chunks)]
(do-nested xchunk ychunk chunks
(*let* [left-x (* chunk-width xchunk)
right-x (+ left-x chunk-width)
bottom-y (* chunk-width ychunk)
top-y (+ bottom-y chunk-width)]
(mpd-displace heightmap
left-x right-x bottom-y top-y
spread))))

We loop over the indices of the chunks, calculate the bounds for each one, and displace it. Simple enough, but a bit fiddly to get right. And now we've finished displacing all the points!

We could have done this recursively, but I was in an iterative mood when I wrote this. Feel free to play around with it.

## Result

Now that we've got a working algorithm we can play around with the values and see what effect they have on the resulting terrain.

**Warning**: don't set the exponent higher than 8 unless you want to crash this
browser tab. An exponent of `8`

results in a heightmap with `66049`

elements.
Going to `9`

is `263169`

, which most browsers seem to... dislike.

And that's it, we've grown some mountains!

Now that we've got some of the boilerplate down I'm planning on doing a couple more posts like this looking at other noise algorithms, and maybe some erosion algorithms to simulate weathering. Let me know if you've got any specific requests!