Volumetric Light Shaft Rendering in GLSL
Volumetric light shafts, commonly known as sunbeams or god rays, are rendered in real-time graphics by simulating participating media along a camera's viewing ray. In GLSL, this is achieved through ray marching, where the fragment shader samples points along a ray traveling from the camera into the scene, tests each point against a shadow map for visibility, and accumulates in-scattered light based on phase functions and atmospheric transmittance. This article explains the step-by-step mechanics of calculating, attenuating, and accumulating sunbeam intensity across discrete ray samples directly within a GLSL shader.
The Ray Marching Model
To evaluate the light arriving at the camera through participating media, the rendering pipeline computes the Radiative Transfer Equation (RTE) numerically. A viewing ray is defined by the camera position as the origin \(\vec{O}\) and a normalized viewing direction \(\vec{D}\):
\[\vec{P}(t) = \vec{O} + t \vec{D}\]
The ray march begins at a near distance \(t_{\text{start}}\) (typically the near plane) and ends at \(t_{\text{end}}\), which is either the scene geometry depth fetched from the depth buffer or an artificial maximum distance for skybox rendering. The ray is partitioned into \(N\) discrete steps of length \(\Delta t = \frac{t_{\text{end}} - t_{\text{start}}}{N}\).
Per-Sample Visibility and In-Scattering
For every step \(i\) along the ray, the world-space sample coordinate \(\vec{P}_i\) is computed. The shader must determine two key factors at this position:
- Visibility (\(V\)): The point \(\vec{P}_i\) is transformed into light projection space to sample the directional shadow map. If the sampled depth is less than the coordinate's light-space depth, the point is in shadow (\(V = 0\)). Otherwise, the sunbeam hits the media (\(V = 1\)).
- Phase Function (\(P(\theta)\)): Atmospheric dust, aerosols, and moisture do not scatter light uniformly in all directions (Mie scattering). The Henyey-Greenstein phase function or the Cornette-Shanks approximation is evaluated using the cosine of the angle \(\theta\) between the light direction \(\vec{L}\) and the view ray \(\vec{D}\):
\[P(\theta) = \frac{1 - g^2}{4\pi (1 + g^2 - 2g \cos\theta)^{3/2}}\]
The asymmetry factor \(g \in (-1, 1)\) dictates the forward or backward scattering bias. For prominent sunbeams, a forward-scattering value between \(0.6\) and \(0.85\) is standard.
Transmittance and Beer-Lambert Law
Light traveling through media suffers extinction (absorption and out-scattering). The proportion of light transmitted over a distance is governed by the Beer-Lambert law:
\[T(d) = e^{-\sigma_t \cdot d}\]
where \(\sigma_t\) is the extinction coefficient. In uniform media, this can be tracked as a running transmittance value:
- Transmittance to Eye: Accounts for light scattered at \(\vec{P}_i\) that survives traveling back to the camera.
- Transmittance from Light: Accounts for sunbeam light reaching \(\vec{P}_i\) from the top of the atmosphere.
Intensity Accumulation Loop in GLSL
The accumulation process combines visibility, phase scattering, and attenuation into a single iterative loop within the fragment shader:
uniform sampler2D shadowMap;
uniform mat4 lightSpaceMatrix;
uniform vec3 lightDir;
uniform vec3 lightColor;
uniform float scatteringCoeff;
uniform float extinctionCoeff;
uniform float gFactor; // Henyey-Greenstein asymmetry
float henyeyGreenstein(float cosTheta, float g) {
float g2 = g * g;
return (1.0 - g2) / (4.0 * 3.14159265 * pow(1.0 + g2 - 2.0 * g * cosTheta, 1.5));
}
vec3 accumulateSunbeams(vec3 rayOrigin, vec3 rayDir, float maxDist, int steps, float dither) {
float stepSize = maxDist / float(steps);
float cosTheta = dot(rayDir, -lightDir);
float phase = henyeyGreenstein(cosTheta, gFactor);
vec3 accumulatedLight = vec3(0.0);
float currentTransmittance = 1.0;
// Offset the starting point with noise to reduce visual banding
float t = stepSize * dither;
for (int i = 0; i < steps; ++i) {
vec3 samplePos = rayOrigin + rayDir * t;
// Transform to light space
vec4 lightSpacePos = lightSpaceMatrix * vec4(samplePos, 1.0);
vec3 projCoords = lightSpacePos.xyz / lightSpacePos.w;
projCoords = projCoords * 0.5 + 0.5;
// Check if sample is within light frustum
if (projCoords.z <= 1.0 && projCoords.x >= 0.0 && projCoords.x <= 1.0 &&
projCoords.y >= 0.0 && projCoords.y <= 1.0) {
float shadowDepth = texture(shadowMap, projCoords.xy).r;
float shadowBias = 0.002;
float visibility = (projCoords.z - shadowBias < shadowDepth) ? 1.0 : 0.0;
if (visibility > 0.0) {
// In-scattered light at this step
float stepExtinction = exp(-extinctionCoeff * stepSize);
vec3 inScattered = lightColor * phase * scatteringCoeff * visibility;
// Accumulate light integrated over step length
accumulatedLight += currentTransmittance * inScattered * stepSize;
currentTransmittance *= stepExtinction;
}
}
t += stepSize;
if (currentTransmittance < 0.01) break; // Early exit if ray is fully extinguished
}
return accumulatedLight;
}Eliminating Banding with Jittering
Discrete sampling produces visible wood-grain slicing artifacts across the screen. To eliminate these patterns without massively increasing step counts:
- Screen-Space Noise: A high-frequency dither value (such as interleaved gradient noise or blue noise) is sampled using screen coordinates and mapped to the range \([0, 1)\).
- Ray Offset: The initial step distance \(t\) is multiplied by this noise value. This jitters the sample locations along the ray per-pixel, converting structural banding into fine, high-frequency noise.
- Temporal Filtering or Blur: The resulting noise is easily resolved using a temporal anti-aliasing (TAA) pass or an edge-aware bilateral blur before compositing the volumetric light over the final scene color buffer.