In this post we will discuss generating fractal terrain. In the previous post we implemented our own fast Fourier transform in order to simulate an ocean surface. In that post we implemented an unnormalized inverse transform. It seemed logical to employ our FFT object as a means of generating terrain, but in order to do so we will need to add a method to compute the forward transform. This will be a relatively brief post. Most of the foundation has been laid with our ocean simulation, so we will simply retool our FFT object for the purpose of generating terrain.

The process we will implement begins with generating a two-dimensional Gaussian white noise. We then take the Fourier transform of the white noise and apply a \( 1/f^\beta \) filter in the frequency domain. This filter has the effect of curbing the amplitudes of the higher frequencies. We then apply the inverse transform and interpret the result as a height field. Below are two screen captures, one for \( \beta=2.4 \) and the other for \( \beta=2.0 \). Periodicity of the Fourier transform allows us to tile our maps. The images below have been tiled.

In the previous post we found our twiddle factors for the inverse transform to be,

\begin{align*}

T_N^x &= \exp\left(\frac{i2\pi x}{N}\right) \\

\end{align*}

If we follow that analysis for the case of the forward transform, we find the twiddle factors,

\begin{align*}

{T_N^x}' &= \exp\left(\frac{-i2\pi x}{N}\right) \\

\end{align*}

So the only modifications we need to make to our FFT object is to precompute a second set of twiddle factors and to apply those in the case of the forward transform. See `src/fft.h`

and `src/fft.cc`

in the project download.

In order to generate terrain using the method of Fourier synthesis our first step will be to generate a Gaussian white noise (see `src/terrain.h`

and `src/terrain.cc`

for the complete implementation).

// generate noise in the spatial domain for (int i = 0; i < N * N; i++) { c[i] = gaussianRandomVariable(); }

A rendering of the white noise is below.

We then apply the forward transform to transform the white noise into the frequency domain.

// transform to the frequency domain for (int j = 0; j < N; j++) fft.fft_(c, c, 1, j * N); for (int i = 0; i < N; i++) fft.fft_(c, c, N, i);

Below we have mapped the real component of the result of the forward transform to the height of the terrain.

In the third step we apply the filter to curb the higher frequencies.

// apply 1/(f^p) filter double f; for (int j = 0; j < N; j++) { for (int i = 0; i < N; i++) { index = i + j * N; f = sqrt((i-N/2)/float(N) * (i-N/2)/float(N) + (j-N/2)/float(N) * (j-N/2)/float(N)); f = f < 1.0/N ? 1.0/N : f; c[index] = c[index]*(1.0f/ pow(f, power)); } }

The result of applying the low-pass filter is below. Note that we have four spikes because we have tiled our map in the \( x \) and \( z \) directions.

Lastly, we transform back into the spatial domain by applying the inverse transform.

// transform back to spatial domain for (int j = 0; j < N; j++) { fft.fft(c, c, 1, j * N); for (int i = 0; i < N; i++) { c[i + j * N] = c[i + j * N] * (1.0f/N); } } for (int i = 0; i < N; i++) { fft.fft(c, c, N, i); for (int j = 0; j < N; j++) { c[i + j * N] = c[i + j * N] * (1.0f/N); } }

We interpret the components of the result of the inverse transform as our height map yielding the rendering below.

In the `cTerrain`

object we have a simple `color`

method to generate a color for a vertex given the height at that vertex.

vector3 cTerrain::color(double h, double min, double max) { // c0 is the color at minimum height, c0 + c1 at maximum vector3 c0( 66.0/255, 40.0/255, 18.0/255); vector3 c1(189.0/255, 185.0/255, 138.0/255); vector3 c = c0 + c1 * ((h - min) / (max - min)); return c; }

In the project we have also evaluated the normal for each vertex and applied our lighting calculations in the vertex and fragment shaders.

Hopefully, in the next post we will be able to implement a level of detail algorithm in order to generate much larger terrains and display them at interactive frame rates. The source code for this project is available as always. The code is a bit sloppy this time around, but have a look at it. If you have any comments, let me know.

Download this project: terrain.tar.bz2

I have a question about your method. Some algorithms and definitions of Brownian surfaces have a parameter ɑ that determines the variance of two neighboring points. If one point's height is f(x,y), and we let h and k be increments, then the resulting random values f(x+h, y+k) - f(x,y) have a Gaussian distribution with a variance of (h^2 + k^2)^ɑ. The Hausdorff dimension of the resulting surface is then d = 3-ɑ.

How does your method account for this ɑ parameter?