Assembling Seamless Meshes

In order to provide a huge variety of trees without requiring each and every tree to be a separate mesh (requiring a huge amount of resources), our preferred approach is - rather than having pre-prepared meshes of full trees - assemble them out of reusable segments which can be combined in a vast number of different ways achieving a huge space of possibilities even within comparatively restrictive memory requirements.

As it's unreasonable to expect anyone to assemble segments seamlessly via simple placing by hand (not to speak of the horrendous inconvenience) - our first component is an automated "segment-analyzer":

Automated Segment Analysis & Categorization

Let's assume, that tree segments are small meshes with multiple openings in the shape of a regular polygon (but with a varying amount of vertices), just like this one:

As should be apparent, this mesh has four openings: One regular heptagon at the bottom and three regular pentagons at the top.

Finding Edges

To get started, let's categorize out vertices in a way where each vert contains a reference to the tris it appears in. We've placed an arbitrary limit of six triangles on this list, as we didn't have any examples where edge-vertices (the only thing relevant for this process), were referenced in more tris than that. In order to identify a point on an opening edge, let's iterate over all triangles it occurs in and note how many times each other vertex is visited by this process. There should be two vertices, which are only visited once - those are the ones that lie on the same edge as the one we started with. If we're unable to find two such neighboring vertices, then move on with the next candidate.

Now, that we've found the edge, we can compute the opening direction (= plane normal), via the cross-product of the vectors to each of the neighboring vertices we selected. As there's always two normalized orthogonal vectors to a plane, let's pick a random point from the triangles we share with any of our selected vertices and compare its direction to the computed opening direction; If it lies in the same hemisphere as our opening direction (= dot product of the direction vs the direction to the shared vertex), we flip our opening direction.

Calculating Openings

Let's take those three points we've selected so far as a starting point and continue further from one end: Again, we iterate over all the triangles this vertex occurs in and mark all vertex occurences - but this time, we ignore all vertices that are already in our current chain (except the first vertex, which we will eventually match to determine that our opening has been completed). Just as before, we pick the one that only appears once - and therefore must lie on an edge; and having found this continuation vertex, we repeat the process until we finally match the first vertex again.

Having calculated the opening direction & corresponding opening plane, we can ensure at every step, that our newly matched vertex actually fits our criteria, and discard it otherwise - potentially even discarding the entire chain in the process.

As we've now found our opening vertices, we average them to determine the center point and take the distance from the first vertex to the center as the radius.

Determining the Rotation Offset

Given a certain center point, radius and plane and number of vertices, there's still infinitely many regular polygons fit this description - we need another factor, a rotation offset from some arbitrary origin, to narrow it down to exactly one unique match, but how can we determine this rotation offset?

There is technically a nice closed-form expression for this conundrum, but unfortunately the inaccuracies appear to explode so quickly, that it isn't of particularly great use to us. What now? - Iteratively narrowing Approximation to the rescue!

Let's first define a transform which should be able to produce our reference point in the plane: (assuming a left-handed coordinate system, with z pointing down)

Now, with the free parameter rotation_offset, we can gather the diff from our observed location vs the values [0 ~ ]τ to find the appropriate value. Rather than blindly trying to iterate over the given value range, we can sample the difference function at a bunch of fixed positions in this range - let's say every ½π - and then grab the two values surrounding the value with the lowest error (n_min - ½π, n_min + ½π, to account for the spherical implications of those factors) and use Brent's Method to find the minimum in the provided range. As there are vertex_count symemtrical variants of each of those matches, we can even reduce the offset further via the following:

Assembling Meshes

Starting with the pre-analyzed segments, we can easily categorize them by the significant opening's number of vertices (we'll pick the largest one - or the one that's closes to [0, 0, 0] as the "significant" or "base"-opening). Now, we start with a randomly picked segment of a chosen category.

To continue, let's first zoom out a bit: Each of the segments is defined in its own segment-space (the direction the significant opening points to, the scale of the opening and its rotation offset). Additionally, there's a local position, rotation & scale at each of the openings of any placed segment in world-space. These two factors have to be considered separately when generating seamless segments. Let's consider those two factors separately, in order to grow our combined mesh in world-space, but correct for the local transforms of each segment, when handling it.

Now, at any given segment, we wanna proceed as follows:

• First, we consider only the category specified by the previous opening (or in the case of the first segment, which has been provided to us via a starting parameter) and select one of the matching openings.

• Now, we calculate the per-segment correction transform via
  (rotation_x(π) * segment_quaternion(segment.base.up, segment.base.rotation_offset)).Inverse()
  (see the implementation of segment_quaternion above)

  Note: Whether we choose to rotate along the x- or y-axis is pretty arbitrary, all we want to achieve is a "flip" along the up/down direction (to have the opening's up direction to pointing down instead, directly into the previous opening's up vector). This could easily be achieved via applying a scale of -1 on th z axis, but as our render instance format doesn't support non-uniform scaling, we've chosen to rotate instead.

• Then we place the mesh at the uninhabited opening (or the start position, in case of the first segment) - using our correction transform to align it with the corresponding opening. next_mesh.scale = previous_opening.radius / next_segment.base.radius;

const vec3f offset = next_segment.base.center_position * next_mesh.scale;
const matrix rotationMat = matrix_from_quaternion(correction_quat * previous_opening.rotation);
const vec3f rotatedOffset = rotationMat.transform(offset);

next_mesh.position = previous_opening.position - rotatedOffset;
next_mesh.rotation = correction_quat * previous_opening.rotation;
• Now, we iterate over all the openings of this new segment and note the "end-point positions" in world space - alongside their radius & orientation. We may here pick from one of the vertex_count possible rotations which produce a seamless transition (via a rotation on the up/down axis of n * τ / vertex_count ∧ n ∈ ℕ, n < vertex_count). const size_t rotation_idx = [0 .. ]opening.vertex_count;
const float variant_rot = (τ / opening.vertex_count) * rotation_idx;

const quaternion corrected_rotation = correction_quat * world_space_rotation;
const matrix rotation_mat = matrix_from_quaternion(corrected_rotation);
const quaternion seg_quat = segment_quaternion(opening.up, variant_rot + opening.rotation_offset);

const vec3f uncorrected_offset = (opening.center - base.center) * current_mesh.scale;
const vec3f end_point_offset = rotation_mat.transform(uncorrected_offset);

endPoint.radius = opening.radius * current_mesh.scale;
endPoint.position = previous_opening.position + end_point_offset;
endPoint.up = rotation_mat.transform(opening.up);
endPoint.rotation = seg_quat * correction_quat * world_space_rotation;

We repeat this process until no open "end-points" remain. This concludes our mesh assembly procedure.

Vertex + Fragment Shader Hybrid Movement Vectors

Most modern deferred-shaded rendering engines use temporal Anti-Aliasing to accumulate subpixel detail over multiple frames by reprojecting information from the previous frame into the current one. To reproject the positions from the previous frame, we need information for each moved pixel (be it via camera movement, by object movement or a combination of the two) of where it resided in the previous frame, a so-called "velocity buffer".

As the motion vectors usually gradually shift over flat surfaces - especially singular triangles - it's in many cases very accurate to just calculate the movement direction in the vertex or geometry shader, to spare the fragment shader from costly per-pixel calculations of information identical or near-identical to what would've been achieved via non-perspective-corrected barycentric interpolation. This works really well for most cases as perspective distortion is usually negligable, however, there's one notable exception:

For exceptionally large triangles, either current or past frame vertices may land significantly off screen - potentially even beyond of the clipping plane. This makes the resulting interpolated values utterly incorrect as they massively overestimate the motion vectors of corresponding pixels.

Vertex Shader Aspect

Now, rather than always calculating movement vectors per pixel to address this, we can build a detection mechanism in the vertex shader, which allows us to conditionally enable per-pixel calculations if necessary:

Fragment Shader Aspect

For the fragment shader, the conditional calculation is now trivial:

Smooth LOD Blending

For modern traditionally-rendered open-world games LOD (level of detail) switching is an essential mechanic by which objects which take up less space on screen are represented by lower-poly versions than if the same model were to be rendered up-close. For many objects the LOD transition (switching out one LOD for another, as the camera moves throughout the world and objects become larger or smaller on screen) can be done somewhat seamlessly, without any fancy techniques (as long as the low poly mesh accurately represents the higher-fidelity version from further away) - or even if not, the resulting "pop-in" effect usually occurs scattered across the screen wherever a mesh happens to pass the megic boundary into a different LOD, however, for dense vegetation (e.g. grass, flowers) - as this transition will continually occur at a consistent distance / position on screen, simply switching out LODs will quickly result in a very apparent transition-point, or a very high number of LOD transitions is needed to prevent this.

In order to tackle this, we can instead simply fade smoothly between different LODs. With a deferred rendering setup, this means that we need to dither between LODs, smoothly shifting the resulting likelyhood of pixel-rejection across the LOD boundary.

The implication of wanting to be able to fade between LODs is that both levels of detail of a mesh need to be loaded & (at least partially) rendered (wherever fades occur), so, simple LOD switching aims towards linear effort to render a model at any given distance (with the complexity of the mesh increasing the more space it takes up), but with fading LODs, the effort to render a model would naively be highest when transitioning between displayed levels of detail.

There is, however, some trickery we can perform, in order to at least reduce the complexity of rendering models in transition...

Distance-Aware Geometry Shader Triangle-Rejection

Our Vegtation-LOD architecture provides each object with a minimum and maximum renderable distance, expressed as a LOD-transition point index - hence, they always occur at pre-calculated distances which can be retrieved from a set of lookup-tables for the fade-in & fade-out distances of any given LOD-transition-point (they actually contain the squared distances in our implementation to speed up the comparison, but that's besides the point).

Now, when we arrive in the geometry shader for any given triangle, we can easily determine if any of the vertices are within the set bounds of our corresponding LOD-regions (or some are behind the visible LOD area and some are in front - to handle triangles spanning multiple LOD regions) - if they fail this check, we can simply reject the triangle outright (by returning without any calls to EmitVertex or EndPrimitive), excluding it from any further processing.

This leaves us with the triangles that are either fully visible or within the blending region. For each vertex of those triangles we can now calculate the fadeInFactor of the first LOD region and fadeOutFactor of the last LOD region. Those fade factors are now valid to be interpolated barycentricly into the fragment shader.

As meshes may span multiple LOD-barriers we have to consider even & odd transition boundaries separately, as - depending on which side of the LOD-barrier we are, when entering a blending region - we need to either keep or reject the fragment with a gradual increase (or decrease) of rejection likelyhood the further we step into the transition region.

This leaves us with nicely transitioning vegetation meshes. The same principle can easily be applied when meshes are rendered into light depth buffers to achieve matching smoothly interpolating LOD shadows (the camera position to calculate the world-space vertex distances from does not become the light position here, but remains as the player camera position).

Shared Texture-Streaming Architecture

Adding vegetation objects as a separate entity type to render objects into the engine - as they require different types of attributes to be associated with the models - leaves us with a challenge:

Our existing texture-streaming architecture is currently targeted towards gathering only render object textures (and updating the corresponding gpu resources pointing to them), but it wouldn't make sense to entirely duplicate this functionality for vegetation objects - and unifying the behavior across different types also scales much better with adding even more different types of rendered objects in the future (it also allows for using a single texture across different types of objects).

Resource Loading Journey

To illustrate the process, let's have a look at the different stages which need to be completed before a texture can be rendered to a certain type of rendered object:

• A Chunk with different types of render objects is added to the renderer
• Basic information about the associated meshes is retrieved (Extents, Mesh File ID, Texture File ID, etc.)
• The data-chunk with generic information is uploaded to the GPU
• When the upload is completed, the chunk can be forwarded to the culling compute shader
• If a certain bounding box is considered visible by the culling compute shader, it's associated object's visible counter will be incremented (for the relevant LOD, if applicable)
• On the next frame the CPU-side iterates over all elements contained in the pools of all rendered object types, finds the ones that are unloaded, but have been referenced on the GPU side and enqueues them into the "Please Load" queue. The objects are now marked as "requested"
• If load-capacity is available the request for the corresponding model & texture is enqueued to the asset-query system
• The model or texture is asynchronously loaded from disk, decompressed & converted to a usable format
• The asset-query handle is checked for completion, if completed, the resources are collected by the renderer-streamer
• The required resources are allocated for storing the asset on the GPU
• The Asset-Data upload to the allocated resource is enqueued, and the corresponding upload-handle is associated with the rendered object / texture data
• When iterating the objects at the beginning of the frame for unloaded visible objects & the "requested" flag is set, we check if there's an associated upload handle and check for upload completion
• If the resource is a texture and the upload is completed, we register the new texture (with the associated texture sampler) in our bindless texture array
• We then update the object's reference to no longer point to unloaded texture / mesh data, but the newly loaded resource.
• Now, when the cull compute-shader finds those instances to still be visible (whose objects now point to loaded resources), it creates draw calls for each of the matching instances
• Finally, the meshes are drawn to the screen with the freshly loaded mesh & texture data

Request States

For textures specifically, this means that - from the perspective of the renderer-streamer - we can split the asset loading process into two phases with various sub-states each:

CPU Load Request

• Awaiting an Asset Query Slot
• Awaiting Asset Data from a running Asset Query
• Alternatively: Awaiting Completion of a different CPU Load Request (created by a different object, but using the same texture)
• Ready (the request will continue as a GPU Update Request)

GPU Update Request

• Allocating the Texture
• Enqueueing the Texture Upload
• Awaiting completion of the Texture Upload
• Setting the texture to the Bindless Texture Array
• Setting the Texture Array Index to the Objects (making it available to the renderer)

State Handling

With this rough overview of the procedure of loading textures in mind, let's consider the following cases and how they're being handled now in a bit more detail: (Steps performed as part of the asset-query subsystem are displayed with reduced focus)

A newly loaded object requires an unloaded texture

In this case, we simply follow the steps detailed above:

• Enqueue Texture asset request → dequeue → forward request, store handle → process request → load asynchronously → mark ready
• Allocate Texture → request GPU upload & store upload-handle → upload completes
• Set texture reference to bindless texture array
• Update entry in the Render- / Vegetation-Object Pool with the bindless texture array index

A newly loaded object requires a texture which is already requested

• Store a dependent texture request with the type & index of the corresponding object
• Every frame, we check if the underlying object has advanced its state from "Waiting for Asset Request"
• If the Object has been advanced, we enqueue our GPU Update Request pointing to the upload handle
• When the upload has been completed, if the texture hasn't been set to the Bindless Texture Array yet, we perform this ourselves (can happen if uploads complete whilst iterating GPU Update Requests)
• Finally, we update the corresponding entity in the associated GPU object pool to point to the index in the Texture Array

A newly loaded object requires a texture which is already uploading

• We merely need to insert a GPU Update Request for the corresponding object
• When the associated upload handle has been marked completed, if the texture hasn't been set to the Bindless Texture Array yet, we perform this ourselves
• Finally, we update the entity in the associated GPU object pool to point to it

A newly loaded object requires a texture which is already available

• Here, we only need to update the object's entity in the associated GPU object pool to point to the corresponding index in the Bindless Texture Array

Rendering Soft Shadows

Achieving Soft shadows with traditional rendering techniques is far from being an accurate depiction of such phenomena in ray-traced scenes, but still significantly improves the overall look of the final output. Whereas soft shadowing in ray-tracing is an inherent property of tracing towards the larger surface of the light fixture, with a traditional rendering pipeline, there is no direct equivalent of this.

Experimenting with sampling the light fixture area over multiple frames and accumulating the result in screen-space has proven tricky, so usually simply accumulating the shadowed-ness of neighboring shadow-map texels is being used for a smoother overall look - but it doesn't come with all of the "real" properties of soft shadows, such as softness increasing with distance from the shadow-casting object. These things can be mimicked though by manually increasing the sampling radius depending on the depth difference between world-space and shadow map.

For sub-texel blending, we can simply subtract the fractional sub-texel position from the sample offsets and throw those positions into a weighting function, then accumulate the weighted samples and divide the sum by the combined weights.

In order to decrease the number of required samples per pixel, we begin by casting a sparse but wide net of shadow samples. If the accumulated result is partially shadowed (so neither all, nor none of the samples are fully in shadow), we collect more samples to resolve the full quality soft shadow.

To determine the smoothing range to be used for distance-based shadow penumbra sizes, we gather the average depth difference from the initial net of shadow samples, use that to interpolate between the minimum and maximum smoothing radius and use the inverse of the resulting value as bounds for our weighting formula.

Temporal Global-Illumination Accumulation via Reprojection

Anti-Aliasing is commonly implemented via temporal multi-sampling and reprojecting previous frames into the current view-space in modern rendering engines. However, the same techniques can be adapted to combat the rough surface ray-tracing noise by accumulating and bluring ray-traced Global-Illumination samples over multiple frames. As we're tracing into a voxelized version of our rendered scene, we can additionally apply a morton sequence to our voxelization offset positions, to further improve the visual quality of the reflected geometry.

As there's quite a bit of overlap, it's very convenient to combine the GI accumulation with the application of Temporal Anti-Aliasing, as e.g. the color and depth differences observed in one, can aid the blending in the other. Unlike the TAA-pass, however, GI is meant to blur on a much larger range in both time and space, as the voxel-traced samples are much more unreliable.

That's why we take a blurred version (restricted by depth differences) of the previous frame as the starting point of the current one and also allow neighboring samples from the current frame to be included in the final color of the current frame's GI as well. That way, GI can slowly "bleed" across surfaces, merely inhibited by stark differences in depth curvature (= edges) between the surrounding samples.

There is, however, one complication to consider here: Global-Illumination sample color is inherently linked to the albedo surface color of which it had been reflected, so without accounting for this surface textures will get very washed out and blurry, as their final color is indirectly influenced by the neighboring texels via their contributions to the blurred GI sample colors.

To combat this, we remove some of the textural information from the GI data, by simply dividing the accumulated GI color by the current frame's albedo color (offset by some constant, to limit the negative impact of incorrectly attributed samples) before storing the final accumulated sample. When referencing this accumulated GI texel during the next frame, we'll obviously have to re-apply the albedo color, before combining it with the newly acquired sample.

It's important to not fall back to the shaded output color here, but rather use the raw albedo color, as otherwise directional differences in shading would have an enormous side-effect on the accumulated Global-Illumination samples.

Another point to consider is that GI samples inherently cannot be interpolated between texels, as otherwise surface information from one surface may bleed into an unrelated other surface through such means.

The Troubled Seas of Voxel-Field Tracing

In contrast to ray tracing against lower poly versions of the rendered scene, where flat surfaces can trivially be represented accurately, there are additional complications with such areas as soon as such surfaces are rasterized to a grid that doesn't exactly conform the directional preferences of the surface plane.

First of all, when resolving from the world-space position in the rendered output to the voxel grid position, we may land inside or potentially even behind the voxelized equivalent.

Secondly, as flat surfaces that don't align neatly with the voxel grid, will inhibit at least some stair-stepping, we could - especially at extremely shallow viewing angles - accidentally encounter a ray intersection with the same surface, by running into such a stair-step.

Lastly, when tracing staring from a voxelized corner position, we may accidentally run into the voxelized remains of a close surface that lies in an irrelevant direction, which occupies our path by simply being too close to the origin surface.

To combat these points, we firstly offset our voxelized position in direction of the current surface normal. If we're still stuck inside a voxel at this point, we allow the voxel tracer a couple of steps through solid voxels, until we'd give up and assume the target voxel to be occupied.

This helps to achieve much more pleasing results than simply doing one or the other, as those factors combine nicely to combat the issues mentioned above pretty convincingly in all but exceptional edge-cases.

Cascaded Shadow Maps

Rendering sun shadows usually comes with a bit of a trade-off: As texture-resolution is limited, shadows may either be high quality, but fade out very quickly, or they fade out much futher away, but are much less detailed.

As a solution to this conundrum, which stays within the traditional rendering pipeline, we can simply render the world from multiple increasingly large orthographic views (biased towards the look direction). This results in having a view with high quality shadows for close-up detail, but also far-reaching shadow map coverage (with lesser quality where it matters less).

When applying the shadows to our scene, we can simply iterate through the list of our so-called cascades, until we end up in-bounds of the light matrix. This is handled as usual, we simply project the world-space position of our fragment to the view-projection matrix of the corresponding shadow-map and end up with a view-space coordinate into the associated texture.

In order to obscure the transitions between the different shadow-maps, we can employ a simple edge-fade, smoothly transitioning from one cascade to another.

Renderer Architecture: Occlusion Culling

As we're generating multi-draw-indirect-buffer for rendered chunks from a compute shader, we'd like to reduce the number of objects being rendered to only consist of what's actually visible on screen.

One such approach is Frustum Culling - a simple technique, which attempts to reject objects behind or next to the View Frustum. It's fairly cheap and easy to implement - simply requiring a bounding box to be projected to the screen, but it's also inherently limited, as it can't combat overdraw.

A more intricate and powerful technique is called Occlusion Culling. Here, we attempt to compare the depth of the projected object or bounding box against a reprojected depth buffer from the previous frame (which may additionally include some static close objects being re-rendered on top, to combat gaps in the reprojected depth buffer).

To make those lookups cheaper, we can construct a mip-chain from the depth buffer, which contains both the minimum and maximum depth value of the represented depth-buffer values (as a best- & worst-case depth). To make using this mip-chain much more straightforward (as we can't just assume that it's valid to split a certain texel down the middle and end up in the right place of the next mip), we can simply restrict the reprojected depth map size to the highest power of two in x and y less or equal to the depth buffer resolution.

Now, with the help of this depth buffer mip chain, we can easily check if the worst case depth of our object - so the closest it could possibly be, derived from its bounding box - can even occupy a certain area on screen (as opposed to being behind some other static object from the reprojected depth).

If our object's worst-case depth lands between the minimum and maximum depth represented by a certain mip texel and is at the edge of the bounding box, we can recursively check if the data stored in the underlying mip's texels makes this decision more obvious. Otherwise, we have to draw the object on screen.

Renderer Architecture: Objects & LODs

To achieve fully GPU driven render object streaming and for some significant throughput advancements, we'd only like to stream full chunks of static render objects to the GPU as a simple list of render object instances to generate a lodded multi-draw-indirect-buffer.

New chunks are simply supplied to the renderer as a list of such render object descriptors (position, rotation, scale). As the descriptors need to be able to look up information about the corresponding render objects, common to all render object instances (bounding boxes, vertex & buffer offsets, etc.), we want this lookup to be a simple index lookup on the GPU side, we allocate space for each distinct render object into a pool of Render Object Information.

Now, for every chunk the compute shader which generates the multi-draw-indirect-buffer is called once, with the offset and size of the corresponding list of render objects.

Here, we can now do all sorts of culling operations and determine the correct LOD to be used for this draw-call.

GPU driven LOD streaming

So far, we've simply assumed, that all LODs have been loaded and are available to the compute shader. However, we most likely don't even wanna stream in all of those LODs - at least in most cases.

By simply adding a requestCount to our Render Object Information for each object, we can simply atomically increment (can also be weighted by importance / size) whenever an LOD is being requested or drawn.

After each frame the CPU-side iterates over the occupied slots of the pool, finds the missing LODs with the most requests and forwards it to the asset-streaming subsystem, which will asynchronously read, prepare and upload the LOD's mesh data to the GPU. Prioritization is handled via the requestCount.

As even loaded LODs are marked this way, the unloading of LODs can be handled in a similar fashion. Obviously, the requestCount needs to be reset after every frame.

Binary Serialization via Type-Seggregated Buffer Compression

In order to serialize generated assets into a custom model format, we need a minimal, fast & adaptable (de-)serialization API.

As discussed by gaffer on games, in order to achieve single-function de-/serialization per object, we introduce a similar API with a custom serialization-context input parameter.

Our serialization context now collects values of different types into separate arrays and when reaching the end of the stream, we can simply RLE-compress, Reference & Integer-Encode (where applicable) all the collected values per type to get a significantly reduced output size. By collecting symbols on-the-fly into a recent, approximate common and most-frequent buffer (upgrading values inside those buffers as we collect them), we can further improve our ease of referencing common values (although, this process is fairly costly upon decompression as well). Rather than just immediately encoding the symbols in a bitstream (which isn't even possible with rANS), we can gather symbol-frequency statistics over the generated symbols into histograms & use those to entropy-encode the symbol stream.

Afterwards, the collected symbols can be written to a file-/memory-/network-bytestream of our choosing.

Further improvements may be achieved, by meta-compressing the type specific symbol streams, in order to possibly combine commonly succeeding symbols into a separate meta-symbol, etc.

One of the limitations of this approach, is that all values need to be collected first (need to reside in RAM at once, rather than being able to be streamed in & out concurrently), and can only be compressed afterwards. However, following a data-oriented approach, this is usually a faster way to do this (compared to the constant context switching of compression & structure descent) and presumably results in a much simpler implementation.

On Noise Patterns

When having a look at our original global-illumination implementation, there was a large amount of sporadic noise, that didn't result in visually consistent and seemingly smooth lighting. The noise source for this implementation had been a white-noise texture, that was being read at locations tied to the screen and world position & was being offset by a random vector for each frame.

Fully random noise (white noise) may lead to long-term convergence to the average, but not necessarily to a stable look over the short term, which is why usually blue (or even custom spacio-temporal) noise is used for such purposes.

This may lead to a more "patterned" looking noise texture, but also to a much higher rate of convergence to a visually consistent state.

Voxel Traced Reflections

As global-illumination rays are cast in the deferred section of the render loop, we have access to roughness information associated with the surface of the first bounce. When encountering shiny surfaces - if we're fully respecting and considering surface properties - this can lead to some nice voxel traced reflections. On a per-frame basis, these reflections clearly display the resolution limitations of the voxel grid, although by jiggeling the voxel grid offset slightly over time, one could achieve a much smoother look (as we can observe in a future article). Unfortunately the limited (but also time-dithered) color accuracy is quite obvious here, as we haven't overhauled the noise pipeline on that end yet.

Comparing Triangle-Collapse and Edge-Collapse

Sadly, Vertex-Elimination is too limited to really produce any interesting results with reasonably high-poly meshes, as it's decimation floor is relatively high for any given mesh. So, let's only have a look at how Edge-Collapse & Triangle-Collapse fare on some more complex models and what their characteristics are when pushed to their limits.

Triangle-Collapse

Edge-Collapse

With the bottom plate we can clearly see how Triangle-Collapse quickly starts eating away even at large overall shapes. In order to get this under control, more "globalist" error metrics would need to be employed - but Edge-Collapse would most certainly similarly benefit from this. On the other hand, this is actually the only case where Edge-Collapse also ends up damaging normals (by flipping the triangle vertex-order). Still, it's far less pronounced as Triangle-Collapse. Ignoring all the issues, I think I prefer Edge-Collapse here over Triangle-Collapse, as it resulted in much smoother looking decimation. I expected the thin paper geometry to cause issues, but this fortunately didn't occur.

Triangle-Collapse

Edge-Collapse

Here we can observe how Triangle-Collapse slowly arrives at an increasing number of problematic cases, where normals are accidentally flipped, whilst Edge-Collapse remains much more stable. I generally prefer the look of Triangle-Collapse here though - that is until it starts cutting of relevant areas entirely due to the more rapid increase of introduced error with every decimation step.

Triangle-Collapse

Edge-Collapse

Lastly, this is one of the cases where Triangle-Collapse (apart from the cut-off podium) produces a much nicer looking output - especially when focusing on the wings and hair of the statue. It somehow is able to achieve a somewhat painterly look on a lot of models.

Performance

As already discussed in the Edge-Collapse entry, Triangle-Collapse is quite a lot faster in its current implementation than Edge-Collapse. For a more direct head-to-head comparison, let's take a look at performance metrics from the bird model shown in the Vertex-Elimination post. Please note, that Vertex-Elimination isn't even reaching anything close to 1/8 the triangle count, as it pollutes its own geometry with too many unfavourable polygons during the decimation process.

[Model of the Cover & second Image: Eber by noe-3d.at is licensed under Creative Commons Attribution-NonCommercial]

Mesh Decimation by Edge Collapse

Triangle-Collapse may be a terrific Mesh Decimation methodology for some scenarios, but for our specific requirements (Terrain LODs with certain immutable vertices at the edge) it isn't optimal as its caveats lie exactly in the areas that hurt the most. So, let's take our "Collapse"-Strategy one step back and - rather than collapsing entire triangles into a single vertex, merely collapse a single edge at a time.

Given the similarity to Triangle-Collapse, the implementation follows along the same lines: Pick an edge, update all triangles that touched the removed vertex to point to the remaining vertex instead and optionally update the remaining vertex position to a new position. We can choose to either collapse onto the center position of the edge, or pick the position of the vertex we consider to be more relevant (in our case the latter strategy is only used if one of the vertices is immutable). If a triangle contains both of the vertices, it can be removed entirely.

Performance

Edge-Collapse is unfortunately yet again a bit slower than Triangle-Collapse and Vertex-Elimination. To collapse the upper body of a high-resolution antique statue scan to 1/256 the vertex count we get the following metrics:

This is unfortunately quite a lot slower for large models, but for our use case of decimating relatively small and low-complexity terrain geometry-fragments, this should still be acceptable.

Complications

Edge-Collapse isn't really more flexible, but more robust than Triangle-Collapse. Like Triangle-Collapse the strategy necessarily introduces averaged vertex positions into the mesh, but this didn't really cause any trouble with vertex-order or massive differences in any relevant scenarios yet. Sure, it's much slower - especially than the blazingly fast Vertex-Elimination, but it's robust and doesn't catch on fire if the source geometry isn't exactly to the specific tastes of the algorithm.

[Model of the Cover Image: Relief Depiction of Ma'at by IPCH Digitization Lab is licensed under Creative Commons Attribution]

Mesh Decimation by Triangle Collapse

As our Vertex-Elimination Decimation had too many caveats to be used as a general purpose decimation algorithm - and even ran into problems with our target of terrain decimation - we'll have to come up with a different solution. So, let's try to address the largest problems we ran into: Dependency on "nice" geometry, that we'd end up ruining eventually anyways & needing to do extra steps to fix vertex-ordering in triangles, if requested.

So how about if we don't start with vertices or lines - which may be part of complex polygons, but rather with triangles: Simply picking a triangle and collapsing it onto its barycenter should allow us to decimate even the most complex of meshes - with the downside, that with every collapse operation, we necessarily introduce a certain amount of inaccuracy, which may carry over into future collapses.

We can simply treat all of a triangle's vertices as references to the same remaining (now position-adjusted) result-vertex and ensure that all of their distinct neighbors reference this replacement vertex instead. If a neighbor triangle shares more than one vertex with our current triangle it can also be removed. We would assume triangle vertex order to stay consistent, as we're simply collapsing three vertices into one and don't really create any new vertices, but as no vertex position is "safe" from our decimation, it may occur, that triangles just happen to be collapsed into new positions that result in a flip of a triangle's vertex-order direction.

Our goal of restricting which vertices remain in place is slightly more difficult here, as we need to either not collapse any triangles that feature vertices that may not be moved or - if only one vertex may not be moved - collapse onto that vertex. In practice both of those strategies aren't flawless, especially since our mesh-generation scenario has to deal with a ton of open edges that may not be moved and need to line up exactly with another mesh's vertices - and the better of the two strategies happens to be the "not touching triangles with immutable vertices" one, which ends up severely limiting the level of decimation possible from a given input set of partially locked vertices.

Error Calculation

The error metric to prioritize triangles that collapse with the least amount of assumed visual difference can be fairly straightforward: Simply summing up the area of the 3d pyramid between the collapsed and non-collapsed vertex positions for each of the affected triangles (except the collapsed one) suffices.

Performance

Triangle-Collapse is unfortunately quite a lot slower than Vertex-Elimination, but still pretty fast. To collapse the upper body of a high-resolution antique statue scan to 1/256 the vertex count we get the following metrics:

The performance of Triangle-Collapse also benefits from the same chunked list data-structure, as implemented for Vertex-Elimination. For a more direct comparison, here's the metrics from running a light 1/8 decimation on the same model as in the Vertex-Elimination entry.

Complications

As already discussed, although much more flexible than Vertex-Elimination, Triangle-Collapse isn't without it's flaws. However, these are mostly limited to triangle vertex-ordering and necessarily introducing new vertex positions with every triangle collapse, so all-in-all it may be very much suitable for some use-cases. Triangle-Collapse can produce almost painterly features and is sometimes able to preserve sharp details much longer than other competing strategies, but unfortunately it wasn't a perfect match for ours.

[Model in the Image below: Abecedary tablet by IPCH Digitization Lab is licensed under Creative Commons Attribution-NonCommercial]

Mesh Decimation by Vertex Elimination

Since Border-Alternating Terrain LODs require full control over which triangles are decimated, a custom mesh decimation algorithm is needed - and even if not, it'd be very convenient to have a fast, reliable & flexible way of decimation at our finger tips.

So decimating a mesh by eliminating the least-relevant vertices seemed like a novel & fast solution to the problem. But first, how do we want to eliminate those least-relevant vertices?

Let's assume, we have this geometry-fragment, and we've determined, that the tip of this shape is the vertex we'd like to eliminate. In order to remove said vertex from our object, we have to ensure, that all neighboring vertices that previously connected to the tip do now connect to one another. So, by observing the other vertices in all triangles connected with our chosen vertex, we can construct the polygon(s) of the resulting area. After sorting neighboring triangles into coherent polygons by finding common (= overlapping) vertices (which can be a bit tricky with "open" meshes or at the rims of mesh fragments, as there may not be any further overlapping vertices after a certain point, but still some triangles left to process), the polygon(s) can simply be triangulated and the chosen vertex has been eliminated.

As simple as this may seem, it sadly isn't that easy to implement in practice: Polygon triangulation in 3D isn't really a thing if the vertices don't all lie on a flat plane and even 2D polygon triangulation isn't trivial at all. As we're merely evaluating the suitability of this approach, we'll simply ignore those cases for now, and disallow any vertices from being eliminated if this would result in such a horrible case.

Regarding our goal of being able to leave singular vertices untouched, to line up with higher-fidelity LODs, this can be achieved very easily here, as the decimation procedure happens per-vertex anyways.

Error Calculation

In order to prioritize vertices that result in minimal change to the mesh, we need to define an error function. Fortunately, this is rather easy: The error of removing the vertex at the tip can be formulated as the area between the re-triangulated base and the tip (= sum of the areas of the tetrahedron between the generated triangle and tip). As there's no singular valid triangulation of a polygon, this calculation unfortunately has to be repeated for all vertices and can't be blindly carried over to neighbors, as they may end up with a different triangulation, that would result in a different error value. To demonstrate, that not only one valid triangulation exists, let's consider the case of the following opening and its two different triangulations:

Performance

This approach is actually blazingly fast, even when forced to re-triangulate the resulting areas from the perspective of each vertex. To decimate the Statue pictured at the top of this entry (left) to the point where all vertex eliminations would result in horrific polygons (right), this is what we end up with:

This required an incredibly fast chunked list data-structure to be written, which allows for blazingly fast random removals & inserts (removing old error values in a sorted list, inserting new error values at the correct position) and a bunch of other little tricks.

Complications

What should now be somewhat apparent from the Performance Output above is, that many of the vertices needed to be skipped, as it would've been too complex to triangulate the resulting geometry. This usually isn't as terrible with smaller models, but heavily impacts the flexibility of this algorithm. What can be somewhat worked around, but isn't particularly pleasing about this approach is, that triangle vertex-ordering (which is sometimes used to determine normal directions and allows for front- or backface culling) isn't automatically preserved and needs to be fixed during the decimation process, if this behavior is relevant.

[Model of the Cover Image: Grabfigur by noe-3d.at is licensed under Creative Commons Attribution-NonCommercial]
[Model in the Image below: Abecedary tablet by IPCH Digitization Lab is licensed under Creative Commons Attribution-NonCommercial]

Border-Alternating Terrain LODs

With connected geometry - like terrain geometry - simply dissecting the model into tiles and providing various LODs for each tile (with larger tile sizes for lower quality LODs) unfortunately leads to problems whenever different LOD layers meet, as with a naive approach LOD levels would disagree at their borders about the level of fidelity. Forcing higher quality rims onto lower quality levels leads to an escalation, where even the lowest quality LODs would need to provide the highest quality borders, even though they would never actually border any of the highest quality levels, but as every lower level needs to be able to connect to a higher and lower quality level seamlessly, this forces highest quality borders onto all layers, which will increase the geometric detail of the lowest-quality LODs so much, that it entirely defeats the purpose of having LOD levels in the first place, as all of the detail is spent on irrelevant high-fidelity borders.

In order to resolve this, LOD levels can't share borders across all layers. This makes streaming LODs and generating LODs a lot more complex, but completely eliminates the issue at hand: Let's assume we have this lowest-quality LOD0 of our terrain.

If we want to load a segment of LOD0, we need to split it into its four children, as otherwise we'd end up rendering some segments twice: Once with LOD0 and more with LOD1.

Now, to address the issue discussed above: We need to choose a different boundary, at which higher-resolution LODs line up with LOD1. So, let's lock the vertices on the dark-turquoise lines to match the fidelity of LOD2.

Let's assume the worst case: We're interested in more detail in the selected tile of LOD2. Simply loading in this tile, we would disagree with the neighboring tiles of LOD1 on the level-of-fidelity at the two borders highlighted in white.

To address this, loading the tile highlighted in white forces the turquoise tiles to be loaded at the same LOD, as they share the same dark-turquoise square.

In order to again address the issue of sharing boundaries between LOD levels, we force the vertices on the dark-green lines to match the vertices at LOD3.

Let's again assume one of the more interesting cases, and let's say we want to load a higher quality LOD at the tile highlighted in green. Similar as before, this leads to a disagreement on the level of fidelity on the borders highlighed in white.

As discussed before, loading in the tile highlighted in white forces the tiles in green to be loaded as well. Actually any of these 4 tiles being loaded would force the remaining 3 tiles to be loaded in as well.

As two of the parent nodes of those nodes haven't been loaded in already, this also forces those parent nodes to be loaded in at a higher LOD as well, namely LOD2.

Implementing this intra-tile dependency isn't actually as tricky as it may seem, as it's always very straightforward to determine, which of the neighboring tiles need to be loaded in as well. A simple look-up-table suffices to store the vectors to the required neighboring ("sibling") tiles from the perspective of each child-index of a parent tile:

Now one can simply recursively iterate over tiles by starting at the root, splitting every tile that is sufficiently close and loading in all sibling tiles that need to be loaded alongside a given child.

Handling Procedural Geometry & Mesh Fragments

As we always need to know the vertices higher-resolution LOD layers in advance, it's highly recommended to generate such LODs with a bottom-up approach, as otherwise it's much more cumbersome to force our current LOD's vertices to match up with higher-resolution ones we haven't even seen yet. When doing so for partially loaded Meshes (especially common when procedurally generating or manipulating the underlying geometry), generating higher LOD layers with matching outer vertices (yes, not the inner vertices, that match up with higher LODs, but ones that need to match up with neighboring tiles of the same LOD level) can be tricky to handle. This is because neighboring LODs will be decimated separately, but we need to ensure that both we and them pick the same vertex positions on the rims of our tiles, as we'd end up with gaps between triangles otherwise. An easy strategy to handle this, is to manually Edge-Collapse every second vertex on the rim onto its neighbor, then mark the remaining rim vertices as immutable, before the entire tile is decimated to the desired level. This doesn't necessarily result in the most pleasing geometry, but it ensures consistency between the tiles and prevents any gaps from forming.

Clustered Deferred Lighting

To significantly reduce the number of required shadow map lookups in a fragment shader, the view frustum can be subdivided into many small frustums in a pre-process (which we approximate by their bounding sphere). We can now check for the intersection of the light-view-frustum and view-frustum-segment-bounding-sphere and reference the relevant light indices of a given segment in a small lookup table.

Retrieving the corresponding cluster index in the fragment shader is very straightforward, as the only information required to resolve the fragment to a cluster is the screenspace texture coordinate and linearized depth ratio between the near and far plane. All of these inputs are already commonplace in deferred shading setups.

Additionally, such a lookup can easily be utilized to implement a basic volumetric lighting shader.

Procedurally generated PBR textures

To fill our procedurally generated worlds with unique looking environments, we render the corresponding colors for a given region directly into the texture atlas.

As material indexes are directly passed to the shader, alongside the world-position of the corresponding terrain, we can select the correct type of texture and procedurally generate the values to the atlas. This results in a lot of variety, even if the same material is used:

Blending is done on by simply applying a noise function to the barycentric coordinates that are used to fetch the per-vertex material index in the fragment shader from a texture. Sadly it's currently not possible to use different noise patterns for different kinds of material combinations, but the current noise pattern results in very pleasing material transitions already (see video in the header)

This way we can create a wide variety of procedural textures and selectively apply them to the underlying terrain. These are all of the material types we currently offer:

Improving and Optimizing the Texture Atlas

The original implementation was already capable of processing hundreds of thousands of triangles per second, but the utilization of the texture area wasn't particularly good, as the coverage quadtree initially preferred placing triangle clusters at high-layer node positions. This left arbitrarily sized holes in the quadtree that are neither particularly easy to find, when looking for available space for future clusters nor necessarily good spots to place clusters into, as the total available area is quickly fragmented by the first couple of clusters.

We'd much prefer clusters to be placed closer to existing clusters, so we'll modify the atlas slightly: After the final dimensions have been found, we allow the quadtree to attempt to move the current bounds to the top left by one node per layer if possible. This achieves much better atlas coverage - still not optimal, but improved to the point that models with double the unit-size can consistently fit into the atlas now. Apart from the coverage benefits, for some models this ends up positively impacting the runtime performance as well, as it's much faster to find empty spots with the improved coverage.

After a bunch of smaller optimization (mostly bloom-filters to avoid costlier randomish-memory-access-dependent-conditionals) we're somewhere between handling 500k - 1M Triangles per Second on a single thread.

Stanford Dragon Texture Coordinates & Model Render

To actually get a decent performance benchmark of the texture atlas, I fed it the entire stanford dragon model (871k Triangles), as well as a bunch of other large standard models (like a simplified version of Stanford Lucy).

Although, with the simplified Stanford Lucy model however (28k Triangles), the old method appears still to be faster - probably because of the larger variance in terms of triangle area:

Fast Texture Atlas for Procedural Geometry

Since world chunks are going to be represented with consistent, but procedural geometry, we need an efficient way of generating texture information for chunks to make manual changes and unique-looking locations possible. To sample and store the corresponding texture, we need an efficient way of generating texture coordinates for our procedural geometry.

As many existing texture atlas implementations aren't particularly fast or assume that the geometry already comes with UV partitioning and texture coordinates that only need to be assembled into one big atlas, I chose to implement a new texture atlas system for the engine.

The foundation of the texture atlas is a coverage-quadtree where nodes are either empty, partially-covered or full. Nodes are stored in a large buffer as there's no need to stream information on demand. The last layer is merely represented as a large bitmap, rather than the actual nodes containing information about the pixel contents. Utilizing this coverage-quadtree it's fairly straightforward to check if a certain amount of space is currently available in the atlas.

The atlas procedure consists of two steps:

Setup

Add all triangles with 3d position and vertex indexes to the atlas. The atlas calculates the side lengths and angles of the triangle and takes note of existing triangles that share an edge with the current one.

Texture Coordinate Generation

After all triangles have been added, they are first sorted by area. Then chunks of triangles are attempted to be added to the atlas branching out from the initial large triangle's edges.

There's some special cases that have to be handled here, like when all three vertexes of a triangle are already present in the atlas and connecting them into a new triangle would result in a reasonable approximation of sizes and angles in respect to the original 3d geometry. Likewise, in order to ensure that the resulting areas don't get out of hand with curved geometry, as angles are used to construct the 2d triangles in the texture atlas and the resulting area may slowly creep up or down in size compared to their 3d counterparts, mismatching triangles are rejected from the current batch and will either be added separately or as part of another triangle-cluster.

Whenever no more connecting triangles would fit into the area available in the texture atlas, triangle bounds are marked as covered and we move on to the next possible collection of triangles - starting with the next largest triangle that hasn't been placed yet.

Obviously this doesn't result in optimal utililzation of the texture atlas, especially as the bounding rectangle of a triangle is a terrible estimate of it's actual coverage, but this comes with a large performance benefit, as no triangle rasterization is required and the texture generation part can be offloaded to the GPU.

Terrain Chunk Streaming & Atmospheric Effects

In order to traverse such a large scale terrain, LOD & streaming need to be considered immediately. Chunks of different sizes are streamed into fixed-size textures and are dynamically loaded when needed.

Whenever a large chunks would contain a chunk with vertexes that are closer than a set distance from the camera, the chunk is split and the contained chunks of higher fidelity are streamed in.

To reduce the pop-in artifacts usually seen in even tripple-A productions, we use local maxima - rather than average terrain height - to construct LODs.

Obviously this leaves visible gaps in the landscape wherever multiple chunk-types meet, but we'll deal with specific meshes to solve this problem at a later date.

A low LOD texture atlas

Atmospheric Effects

The sun position is now 'rotated' into the coordinate space mapped around the 'flat' planet and camera position using a vector alignment matrix:

The way the sun passes through the atmosphere is approximated with an iterative shader that draws an equirectangular view and approximates rayleigh and mie scattering with a fixed number of samples. The scale of the planetary system and atmosphere have been adjusted to achieve earth-like atmospheric conditions.

Micro-Biome Extraction

After the terrain has been generated and eroded, regions need to be classified into biomes based on the physical properties of their surroundings. To calculate these biomes like various grass types wherever the top soil is grass a noise map can be sampled. However this isn't entirely accurate as especially regional & global elevation changes and water proximity can influence biome parameters significantly.

Therefore we compute a signed distance field of water and one that maps terrain height, blur the height map over a large area and then compute the local difference to an areas neighbors. This leaves us with the water map and the map shown above for local terrain elevation changes.

Using these two maps and a noise buffer, we can derive micro-biome features like mossy stones near water or in crests, small puddles can be turned to mud, large elevation change inhibit the growth of trees or tall grass, ...

Wind Erosion

With our newly developed surface wind simulation, we can now erode the surface terrain. Terrain that's close to the wind height will be eroded more than terrain that's further from it and wind speeds influence how far particles are dragged when eroded. The various terrain layers have different properties of how easily they can be erroded by wind & water.

As the eroded particles are dragged through the wind, the terrain is iteratively checked for matching height as the particles slowly descend based on their weight.

Similarly to the wind simulation, eroded particles are added onto a separate buffer and are applied to the new terrain after all cells have completed their erosion procedure.

Orbital Dynamics

In order to have semi-accurate climatic conditions producing terrain with only one hot and one cold zone, an somewhat excentric elliptic orbital path with fine tuned rotational velocity and direction is required. Solving orbital dynamics formulas - even with a fixed star that isn't affected by the satellites mass - proved quite difficult however.

In Polar Coordinates, the orbital distance from the origin at angle theta with orbital parameters a and b can be expressed like this:

This alone doesn't help much however, as the velocity at each position can easily be calculated, but woulnd't translate to a closed form solution and I didn't want to resort to step-wise simulation that'd be hard to rely on or accurately predict with different time steps. So, I was able to construct a function angle(time) to derive the radius by pluggin the result into the radius function:

In order to arrive at a combined mass that actually produces a period of 2 * pi, we have to calculate the mass from orbital parameter a:

As the combined mass of the two orbital bodies doesn't change with the system, it can be calculated on initialization when setting the other orbital parameters.

Planets, Average Solar Energy on the Surface, Current Solar Energy

Layer Based Terrain

Terrain is built in layers, the bottom layer is non-erodable bedrock. Depending on a signed-distance field, a noise layer and two cubic slopes a certain height of terrain of a given layer will be added on top of the existing terrain. The signed distance field chooses the min and max value from the two cubic slopes (basically a min- and max height) and the noise value (between 0 and 1) is used to smoothly interpolate between those two values. The generated terrain is rough at first, but erosion will massively improve and shape those layers later.