The procedural generation of terrain is mostly a straightforward process due to the fact that most terrain can be modeled as a function over 2D space, i.e. f(x,z) = y. This function can be generated with fractals or noise or some other method, but the end result is always a height map. Due to the functional nature of this representation, it is impossible for there to be multiple y values at any given (x,z) coordinate. This is a acceptable restriction for typical terrain because gravity and erosion usually remove any structures or spaces that are stacked on top of one another. However, complex 3D formations such as caves cannot be represented this way. The goal of this project is to procedurally generate and render a plausible cave-like structure.
Since, as noted above, a function over (x,z) is insufficient, we will rectify this by adding a degree of freedom and creating a function over (x,y,z). Instead of a function to generate a height field in 2D, this function will create a density field in 3D. We will define this density function as a summation over a list of individual structures that define specific cave elements. Structures are themselves density functions, but over a limited area. The surface of the cave can then be extracted by specifying a threshold that separates "air" from "stone" and rendering the isosurface that results.
The generated caves are made up of four basic structures. The first is a ball, which will be familiar if you have seen or implemented metaballs. The density function of the ball is (1 - r^2)^2, where r is the distance from the center of the ball. This function was selected for its computational simplicity and finite support. Finite support allows us to assume a density of zero outside of a specified region, which is useful for an optimization I'll mention later. All of the structure's density functions are based on this drop-off curve.
Another structure is the line, which is effectively two balls connected by a cylinder. The density along the line portion is determined by the closest distance to the line. The most complex structure is the stalactite, which is based on the line in shape, but has a falloff curve along the line's axis as well. The resulting profile curve matches the shape of a quadratic. While it would have been possible to simply use a cone, this shape better captures the form of real stalactites. The last structure (not show below) is the ground plane, whose density function is simply the distance from a plane. The basic structures can be seen below with their normals rendered for clarity.
If we were to use the basic structures as they are, the resulting cave would be unnaturally smooth due to the regular nature of the density functions. To add a bit more variety, we will add some noise to each structure's density function. Specifically, we will use Perlin noise, which conveniently is also defined over three dimensions. Adding noise to the density function has the effect of a displacement map. The surface will be perturbed in a smooth, pseudorandom fashion.
The overall layout of the cave is specified by a graph. Each node represents a "cavern" and is made up of a cluster of one or more ball structures. Node location is determined by subdividing the given space into a given number of buckets and randomly choosing a point within each bucket. Nodes are connected by line structures. The lines are segmented into smaller pieces and perturbed to give their path more variation. Edges in the graph are added by connecting the closest two nodes until the graph is fully connected. The highest node in the graph is then connected to the ground. The graph structure of the cave can clearly be seen in the image below. This overall layout can be obfuscated by adjusting the various parameters, but it is not as apparent from within the caves.
Stalactites are placed by picking a random point within one of the ball structures that make up a cavern. The x and z coordinates of this point are then used to look up a value in the Perlin noise function. The absolute value of this noise is used as the probability that the given stalactite will be placed. This creates clusters of stalactites with varying density. If the location has been chosen, then we step up through the density function until we find the isosurface, which should be the roof of the cavern. Once we have placed the stalactite, a stalagmite is placed below it with a given probability. If the length of the stalactite and stalagmite is long enough, then they will form a column.
The isosurface of the cave is discretized into triangles by running the marching cubes algorithm. A normal is also associated with each vertex, which is computed from the gradient of the density function. While the "cotton-candy" look of the normal visualization is informative, I also added the ability to view the caves in a more realistic manner. The primary source of illumination in this method is from a spotlight that is fixed to the viewing direction. The intensity of the light falls off radially and through space. The surface of the cave is rendered with a procedural bump map that samples a 3D texture that holds a few octaves of Perlin noise. Color is also generated from this texture by blending between two colors based on the noise at a given point. An attempt was made to emphasize the translucent nature of the stalactites and stalagmites with a pseudo subsurface scattering, but the results were not that impressive.
The naive implementation of marching cubes requires time proportional to the number of voxels in the specified region. This is extremely wasteful since the surface of the cave only takes up a fraction of the total space. I achieved a dramatic acceleration by implementing an incremental version of marching cubes that found a voxel on the isosurface and then only evaluated neighboring voxels that were also on the isosurface. The sampling of the density function can also be sped up. In the naive version, one must sum over all structures to obtain a density value. However, since the structures have finite support they will not affect the density outside a given range. By grouping the structures into a bounding volumn hierarchy, I obtained another significant speed up.
Using a density function is not without its flaws. A few oddities can be found, especially at lower voxel resolutions. Some stalactites suffer from bad normals at their tips, which suggests there may be a discontinuity in the stalactite density field. Occasionally you will see objects floating in the air. This is caused by a combination of low voxel resolution and too much noise. You may also see a rather rough transition from the cave tunnels to the ground, but this can also be alleviated by increasing the voxel density.
The method of creating an isosurface from a structure-summed density function produced reasonably good results. While one can imagine a method that tries to connect discrete parts into a similar structure, it would ultimately be less organic and have less variety since the pieces would have to fit together in predefined ways. Certain cave structures such as columns arise naturally from the interaction of individual structures. While the cave does tend to look manufactured when viewed from the outside, the inside is full of interesting detail.
