In this blog post I will describe how my planet terrain erosion simulator handles rock layers and deformation of layers. Here is a picture of the current results:

You can clearly see the multiple layers and their deformation. Each type of layer has its own erosion intensity, which influences how the final terrain looks after erosion (differential erosion). This is generated all in real time (10ms per 64x64 tile) on the CPU without the use of traditional voxels. Voxels would be much too slow.

Virtual Layers

To make it feasible to have many rock layers without actually storing them, I use One Weird Trick™. The main idea is that since the terrain top surface is a 2D heightmap, we only need to generate the layers that intersect the top surface, as shown in the following image:

Here we have started with a predefined set of global layers, each with its own thickness and material type. I have drawn a dark line which represents the raw terrain elevation which is generated by applying math functions to fractal noise (FBM type). You can see that at each horizontal position, we only have 1 or 2 layers exposed to the surface. Therefore, all we need to do to have as many layers as we want (thousands) is to determine the layer index that intersects the elevation surface at each position. In practice, since we are applying erosion after this step, we need to also include a layer or two below that one as a margin in case the erosion removes the top layer.

To determine the layer index, I precalculate the elevation of the top of each layer in an array. Then, at each position, I use a linear SIMD-accelerated scan through the elevations to find the right layer (from top to bottom). Binary search could be used here but isn't faster for <1000 layers. The layer search can be sped up significantly by determining for each tile the minimum and maximum layer indices, which are used to limit the search space within a tile.

Layer Spatial Variation

So far, I have assumed that we have the same set of layers everywhere. However, in reality the rock strata vary across a planet's surface. The question is how can we do something similar without requiring expensive computations for material layer generation?

It's not feasible to spatially vary layer thickness, since that would require determining the top elevation of each layer for each pixel, which involves first generating the thicknesses for all overlying layers, then summing them up, then finally doing the layer search for each position.

The approach I have settled on is to have N different sets of layers on a planet, and to switch between those layer sets based on fractal noise. Each layer set has the name number of layers (e.g. 2048), but different thicknesses and material types. In the next image there are 2 different layer sets which are randomly chosen based on the location.

This approach allows for the layers to be different in different places. For instance, we could have many tiny layers in an area populated with sedimentary rocks, but could have no layers in an area with only intrusive igneous rocks. The number of possible layer sets could be unlimited, but I'm aiming for a few hundred different layer sets on a planet.

Layer Material Spawning

Once the layer index (and layer set) is determined at each position, I use a CDF-based sampling method to determine the actual material at each position within a given layer. This means that each layer can have multiple materials that spatially vary across the planet surface. Each material is given a spawning probability, and this is used to precompute for each layer the CDF of the material probabilities. Then, binary search is used to find the material for each position based on some input fractal noise. This produces materials that are spawned with the right frequency. The layer index is used to hash the material CDF sampling noise, so that each layer can have a different material spawned.

Layer Deformation

The final thing to implement is the deformation of the layers in order to produce nice results like I showed in the first image. The key thing to recognize is that everything we have done so far is relative to the top of the topmost layer, which is currently flat. When we search for the layer index to use at each position, we can apply a spatially-varying offset to this top layer surface so that it is not flat and has deformations similar to those produced by metamorphic processes on Earth:

To make this work, we first need to calculate the depth of the terrain surface below the top of the topmost layer, assuming it was flat. This is just the sum of the layer thicknesses above the terrain. Then, a random noise offset is added to the terrain surface depth value. This final depth value is used to look up the layer index at each position. This causes the layers to modulate relative to the terrain surface for almost no extra compute.

In total, this layer system spawns materials for a 97x97 tile in only 0.2 ms, with 2048 layers. Doing something similar with traditional voxels would be infeasible in real time.

Here are more results of this system: