GLSL Caustics from Overlapping Wave Normal Maps

Animated caustic patterns simulate the focusing of light rays through turbulent water surfaces onto underwater geometry. In real-time rendering with GLSL, this effect is simulated by animating and blending multiple scrolling wave normal maps, calculating how incoming directional light refracts through the combined surface, and deriving light intensity concentrations from the divergence or spatial overlap of these refracted rays.

1. Animating and Blending Normal Maps

Water waves consist of multiple interfering frequencies. To replicate this without complex physics simulations, two distinct tiling normal maps (or one map sampled at different scales) are translated across the water surface using time offsets:

vec2 uv1 = worldUV * scale1 + time * direction1;
vec2 uv2 = worldUV * scale2 + time * direction2;

vec3 normal1 = texture(normalMap, uv1).rgb * 2.0 - 1.0;
vec3 normal2 = texture(normalMap, uv2).rgb * 2.0 - 1.0;

To preserve surface slopes and avoid flattening high-frequency details, combine these tangent-space normals using Reoriented Normal Mapping (RNM) or partial derivative blending instead of simple linear interpolation:

vec3 t = normal1 + vec3(0.0, 0.0, 1.0);
vec3 u = normal2 * vec3(-1.0, -1.0, 1.0);
vec3 blendedNormal = normalize(t * dot(t, u) - u * t.z);

2. Computing Refraction Vectors

Light bends at the boundary between air and water according to Snell’s Law. GLSL provides the native refract function to calculate the direction of the light ray entering the water. Using an Index of Refraction (IOR) of approximately \(1.333\) for water (\(1.0 / 1.333 \approx 0.75\) relative index):

vec3 lightDir = normalize(sunDirection);
vec3 refractedLight = refract(lightDir, blendedNormal, 1.0 / 1.333);

This refracted vector defines the trajectory along which the light propagates toward the seabed or submerged surfaces.

3. Estimating Convergence and Caustic Brightness

Real caustics are envelopes of light rays converging onto small surface areas. Because real-time forward or deferred rendering paths cannot trace millions of photons per frame, GLSL approximations infer convergence directly from the refracted field.

Method A: Grid Ray Displacements and Screen Derivatives

Displace the underwater coordinates using the lateral components (xy) of the refracted vector scaled by the water depth:

vec2 displacedUV = receiverUV + refractedLight.xy * (waterDepth / max(0.1, -refractedLight.z));

The concentration of light at a given point is inversely proportional to the area change of the displaced surface mapping. This is evaluated using GLSL screen-space derivative functions (dFdx, dFdy):

vec2 ddx = dFdx(displacedUV);
vec2 ddy = dFdy(displacedUV);
float area = abs(ddx.x * ddy.y - ddx.y * ddy.x);
float causticIntensity = clamp(1.0 / (area * sensitivity), 0.0, maxBrightness);

Where rays diverge, the area increases, causing causticIntensity to drop. Where rays converge, the area approaches zero, producing bright focal lines.

Method B: Overlapping Dual-Layer Trajectories

A performance-friendly alternative computes the relative convergence directly between two independently refracted sample passes:

vec2 coordA = receiverUV + refract(lightDir, normal1, 0.75).xy * depth;
vec2 coordB = receiverUV + refract(lightDir, normal2, 0.75).xy * depth;

float caustic = dot(texture(noiseMap, coordA).r, texture(noiseMap, coordB).r);
caustic = pow(caustic, sharpnessExponent);

4. Enhancements: Chromatic Aberration and Depth Attenuation

Water disperses light across different wavelengths at slightly different angles. Simulating this in GLSL involves slightly altering the IOR per color channel before deriving light brightness:

vec3 eta = vec3(0.750, 0.752, 0.754); // Red, Green, Blue refraction ratios
vec3 rR = refract(lightDir, blendedNormal, eta.r);
vec3 rG = refract(lightDir, blendedNormal, eta.g);
vec3 rB = refract(lightDir, blendedNormal, eta.b);

Finally, caustic visibility diminishes as water depth increases due to absorption and scattering. Modulate the computed caustic pattern using an exponential extinction factor based on the path length through the water column:

float depthFade = exp(-waterDepth * absorptionCoefficient);
vec3 finalCaustic = causticRGB * depthFade * sunColor;