Terrain Generation: Spectral Synthesis vs. Diamond-Square

A comparative study of two classical terrain-generation techniques — and what happens when you combine them.

Engine: Unity · Language: C# · Team: 2 Module: IS71021E Mathematics for Games and VR/AR (2023–24), Coursework 2 — MSc Computer Games Programming, Goldsmiths University of London

---

The Question

Diamond-Square and frequency-domain spectral synthesis both produce plausible terrain, but they arrive there from opposite directions. Diamond-Square works in space — recursively subdividing a grid and displacing midpoints by a diminishing random offset. Spectral synthesis works in frequency — transforming noise into the frequency domain, attenuating according to a power law, and transforming back.

We implemented both, then asked a third question the brief didn't: what does terrain look like if you mix them?

The project runs all five modes at runtime, switchable from the UI:

public enum meshChoice { FourrierTransform, DiamondSquare, AverageMix, AddMix, Shader };

---

Results

Diamond-Square	Averaged Mix	Additive Mix

The two mixing strategies behave very differently. Averaging the heightmaps ((ds + fft) / 2) preserves the overall silhouette while smoothing the sharp ridge artefacts Diamond-Square produces at subdivision boundaries. Adding them (ds + fft) compounds amplitude instead, keeping Diamond-Square's ridges and layering spectral detail on top — more dramatic, and much easier to push into unusable extremes.

A shader-driven generation mode is included as a fifth comparison point.

---

Spectral Synthesis (FFT.cs)

The pipeline is: generate noise → transform to frequency domain → filter → transform back → displace vertices.

Noise. Three octaves of Perlin noise are sampled at different scales and offsets, averaged, then passed through a smoothstep contrast curve:

combinedSample = combinedSample * combinedSample * (3.0f - 2.0f * combinedSample);

Transform. Each vertex's greyscale intensity is sampled bilinearly through the mesh's UVs and loaded into a System.Numerics.Complex array as a real value, then transformed by the 2D discrete Fourier transform:

$$F(u,v) = \frac{1}{NM}\sum_{x=0}^{N-1}\sum_{y=0}^{M-1} f(x,y), e^{-2\pi i\left(\frac{xu}{N} + \frac{yv}{M}\right)}$$

Filter — this is the part that actually makes it terrain. Each frequency component is attenuated by its distance from the spectrum's centre, raised to an exponent:

float distance = Mathf.Sqrt(Mathf.Pow(u - halfN, 2) + Mathf.Pow(v - halfM, 2));
float f = Mathf.Max(distance, 1);      // avoid division by zero
inputData[u, v] /= Mathf.Pow(f, r);

That is a 1/f^r power-law filter — the classic spectral synthesis approach to fractal terrain. Low frequencies (broad landmasses) pass nearly untouched; high frequencies (fine roughness) are suppressed in proportion to r. Varying the single exponent r moves the output continuously between rolling hills and jagged noise, because it is directly controlling the fractal dimension of the surface.

Then the inverse transform returns the filtered spectrum to a heightmap, which drives vertex Y positions.

On the name: this is a DFT, not an FFT

The class is called FFT, but it implements the direct discrete Fourier transform — a four-deep nested loop evaluating every (u,v) against every (x,y). We documented the cost in the source rather than hiding it:

"Since this is a direct implementation of the DFT algorithm, the time complexity is O(n^4)"

For a 128×128 grid that is ~268 million complex exponentials per transform, which is why generation is not interactive. A Cooley-Tukey decomposition would bring it to O(n² log n) — implementing the naive form first was deliberate, since it maps line-for-line onto the summation formula and made the filter behaviour easy to reason about, but it is the single biggest weakness here.

---

Diamond-Square (diamondSquare.cs)

Diamond-Square seeds the four corners of the grid, then alternates two steps at halving intervals: the diamond step sets each square's centre from its four corners, and the square step sets each diamond's centre from its four neighbours — each with a random offset scaled down by a smoothness exponent every iteration.

The file deliberately retains three implementations:

//performDiamondSquare();                    // Version 1 — self-written
//performDiamondSquareOGLDEV(smoothness);    // adapted from the OGLDEV reference
performDiamondSquareV2();                    // Version 2 — in use

They are kept rather than deleted because the comparison between them was part of the coursework: writing a naive version, studying an established reference implementation, and then rewriting with what that revealed about boundary handling and offset scaling.

---

What I'd Do Differently

• Implement a real FFT. The O(n⁴) DFT is the dominant limitation — it blocks interactive regeneration and caps the practical resolution at 128×128.
• Move the transform off the main thread, or into a compute shader. The whole pipeline is data-parallel and currently runs synchronously in Start().
• Mathf.Pow(u - halfN, 2) for squaring is markedly slower than multiplying the value by itself, and it sits in the innermost filter loop.
• The class name lies. FFT should be DiscreteFourierTransform, or the algorithm should match the name.

---

Credits and Sources

Team: Laurent Klein, Samuel Ding.

The DFT formulation follows the derivation in Dickerson, Terrain Generation (Williams College), cited inline in FFT.cs. One Diamond-Square variant (performDiamondSquareOGLDEV) is adapted from the OGLDEV tutorial implementation and is named accordingly in the source. A YouTube walkthrough referenced during development is linked at the top of FFT.cs.

---

Related

• procedural-dungeon — dungeon generation via branching random walkers and bitmask autotiling.
• cubic-spline-interpolation — Bézier, Hermite, Catmull-Rom and B-spline from scratch, from the same mathematics module.