Volumetric Raymarched Clouds — Realistic Cumulus with Silver Linings and Light Scattering

Volumetric clouds remain one of the most visually rewarding — and technically demanding — effects in real-time rendering. Titles like Horizon Zero Dawn, Microsoft Flight Simulator, and Red Dead Redemption 2 have raised the bar for what players expect from skies. At the heart of all these systems is the same core idea: march a ray through a 3D volume, sample a density function, and accumulate light scattered toward the camera. In this post we build a complete, self-contained GLSL fragment shader that produces animated, sunlit cumulus-style clouds over a gradient sky — no textures required.

Complete Standalone Shader

Paste this directly into any WebGL1 fragment-shader preview (like the shader player on this forum). It expects only iResolution and iTime. Everything — noise, density modeling, lighting, phase function — is procedural.

precision mediump float;

uniform vec2 iResolution;
uniform float iTime;

// Hash functions for noise
float hash(float n) {
    return fract(sin(n) * 43758.5453123);
}

float hash3(vec3 p) {
    p = fract(p * 0.3183099 + 0.1);
    p *= 17.0;
    return fract(p.x * p.y * p.z * (p.x + p.y + p.z));
}

// 3D value noise
float noise3D(vec3 p) {
    vec3 i = floor(p);
    vec3 f = fract(p);
    f = f * f * (3.0 - 2.0 * f);

    float n = i.x + i.y * 157.0 + 113.0 * i.z;

    float a = hash(n + 0.0);
    float b = hash(n + 1.0);
    float c = hash(n + 157.0);
    float d = hash(n + 158.0);
    float e = hash(n + 113.0);
    float g = hash(n + 114.0);
    float h = hash(n + 270.0);
    float k = hash(n + 271.0);

    return mix(
        mix(mix(a, b, f.x), mix(c, d, f.x), f.y),
        mix(mix(e, g, f.x), mix(h, k, f.x), f.y),
        f.z
    );
}

// FBM with 5 octaves for cloud density
float fbm(vec3 p) {
    float f = 0.0;
    float amp = 0.5;
    float freq = 1.0;
    for (int i = 0; i < 5; i++) {
        f += amp * noise3D(p * freq);
        freq *= 2.02;
        amp *= 0.5;
    }
    return f;
}

// Low frequency FBM for large cloud shapes
float fbmLow(vec3 p) {
    float f = 0.0;
    float amp = 0.5;
    float freq = 1.0;
    for (int i = 0; i < 3; i++) {
        f += amp * noise3D(p * freq);
        freq *= 2.0;
        amp *= 0.5;
    }
    return f;
}

// Cloud density function
float cloudDensity(vec3 p) {
    // Cloud slab from y=3 to y=8
    float cloudBase = 3.0;
    float cloudTop = 8.0;
    float cloudMid = 0.5 * (cloudBase + cloudTop);
    float cloudHalfH = 0.5 * (cloudTop - cloudBase);

    // Height within cloud layer [0,1]
    float hNorm = (p.y - cloudBase) / (cloudTop - cloudBase);
    if (hNorm < 0.0 || hNorm > 1.0) return 0.0;

    // Height-based density envelope: thicker in middle, rounded falloff
    // Peaks around 30-40% height (cumulus shape)
    float heightEnv = smoothstep(0.0, 0.15, hNorm) * smoothstep(1.0, 0.6, hNorm);
    // Make bottom flatter and top more rounded
    heightEnv *= mix(1.0, 0.7, pow(hNorm, 0.5));

    // Wind animation
    vec3 wind = vec3(iTime * 0.3, 0.0, iTime * 0.1);

    // Large scale shape noise (creates distinct cloud puffs)
    float largeScale = fbmLow((p + wind) * 0.15);
    // Create distinct puffs by using a coverage threshold
    float coverage = 0.45;
    float baseShape = smoothstep(coverage, coverage + 0.25, largeScale);

    if (baseShape < 0.01) return 0.0;

    // Medium detail noise
    float detail = fbm((p + wind) * 0.35);

    // Combine: base shape modulated by detail
    float density = baseShape * heightEnv;
    // Erode edges with detail noise
    density -= (1.0 - detail) * 0.35;
    density = max(density, 0.0);

    // Add wispy detail at edges
    float fineDetail = noise3D((p + wind * 1.5) * 1.8);
    density += density * fineDetail * 0.3;

    return density * 1.2;
}

// Henyey-Greenstein phase function
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));
}

// Light marching for self-shadowing
float lightMarch(vec3 pos, vec3 sunDir) {
    float density = 0.0;
    float stepSize = 0.6;
    vec3 p = pos;
    for (int i = 0; i < 6; i++) {
        p += sunDir * stepSize;
        density += max(cloudDensity(p), 0.0) * stepSize;
        stepSize *= 1.2; // Increase step size for efficiency
    }
    return density;
}

// Sky color
vec3 skyColor(vec3 rd, vec3 sunDir) {
    // Blue sky gradient
    float t = max(rd.y, 0.0);
    vec3 sky = mix(vec3(0.6, 0.75, 0.95), vec3(0.25, 0.45, 0.85), t);

    // Sun glow
    float sunDot = max(dot(rd, sunDir), 0.0);
    sky += vec3(1.0, 0.9, 0.7) * pow(sunDot, 64.0) * 0.8;
    sky += vec3(1.0, 0.8, 0.5) * pow(sunDot, 8.0) * 0.15;

    // Horizon haze
    float horizon = 1.0 - abs(rd.y);
    sky = mix(sky, vec3(0.8, 0.85, 0.9), pow(horizon, 8.0) * 0.5);

    return sky;
}

void main() {
    vec2 uv = (gl_FragCoord.xy - 0.5 * iResolution.xy) / iResolution.y;

    // Camera setup
    vec3 ro = vec3(0.0, 2.0, 0.0); // Camera at y=2, looking at clouds
    vec3 rd = normalize(vec3(uv.x, uv.y + 0.3, 1.0)); // Slight upward tilt

    // Sun direction (high and to the right)
    vec3 sunDir = normalize(vec3(0.6, 0.55, 0.4));
    vec3 sunColor = vec3(1.0, 0.95, 0.85);

    // Background sky
    vec3 col = skyColor(rd, sunDir);

    // Raymarch through cloud slab [y=3, y=8]
    float cloudBase = 3.0;
    float cloudTop = 8.0;

    // Only march if ray points upward and can hit cloud layer
    if (rd.y > 0.01) {
        float tMin = (cloudBase - ro.y) / rd.y;
        float tMax = (cloudTop - ro.y) / rd.y;

        if (tMin < 0.0) tMin = 0.0;
        if (tMax < tMin) {
            // Ray doesn't hit slab
            gl_FragColor = vec4(col, 1.0);
            return;
        }

        // Limit max distance for performance
        tMax = min(tMax, 80.0);

        float numSteps = 50.0;
        float stepSize = (tMax - tMin) / numSteps;

        // Accumulated color and transmittance
        vec3 cloudCol = vec3(0.0);
        float transmittance = 1.0;

        float cosAngle = dot(rd, sunDir);
        // Dual-lobe phase function for silver lining + backscatter
        float phase = henyeyGreenstein(cosAngle, 0.6) * 0.7
                    + henyeyGreenstein(cosAngle, -0.2) * 0.3;
        // Normalize roughly
        phase = mix(1.0, phase, 0.7);

        float t = tMin;

        // Add small jitter to reduce banding
        t += stepSize * hash(dot(gl_FragCoord.xy, vec2(12.9898, 78.233)));

        for (int i = 0; i < 50; i++) {
            if (transmittance < 0.01) break;

            vec3 pos = ro + rd * t;

            float density = cloudDensity(pos);

            if (density > 0.001) {
                // Light marching for shadows
                float lightDensity = lightMarch(pos, sunDir);

                // Beer-Lambert absorption
                float lightTransmittance = exp(-lightDensity * 1.5);

                // Powder/sugar effect (fluffy bright edges)
                float hNorm = (pos.y - cloudBase) / (cloudTop - cloudBase);
                float powder = 1.0 - exp(-density * 4.0);
                powder = mix(powder, 1.0, 0.5); // Don't let it go too dark

                // Combine lighting
                vec3 ambient = vec3(0.55, 0.6, 0.75) * 0.4; // Sky ambient
                vec3 directLight = sunColor * lightTransmittance * phase * powder * 1.8;

                // Add extra brightness for cloud tops
                float topBoost = smoothstep(0.4, 0.9, hNorm);
                directLight += sunColor * topBoost * lightTransmittance * 0.3;

                // Dark bottoms with slight warm/cool variation
                float bottomDark = smoothstep(0.3, 0.0, hNorm);
                ambient *= (1.0 - bottomDark * 0.4);

                vec3 luminance = ambient + directLight;

                // Cloud albedo (slightly warm white)
                vec3 cloudAlbedo = vec3(1.0, 0.99, 0.97);
                luminance *= cloudAlbedo;

                // Beer-Lambert for this step
                float stepTransmittance = exp(-density * stepSize * 2.5);
                float absorption = 1.0 - stepTransmittance;

                cloudCol += transmittance * absorption * luminance;
                transmittance *= stepTransmittance;
            }

            t += stepSize;
        }

        // Blend clouds with sky
        col = cloudCol + transmittance * col;
    }

    // Tone mapping (simple Reinhard)
    col = col / (1.0 + col);

    // Gamma correction
    col = pow(col, vec3(1.0 / 2.2));

    // Very slight contrast boost
    col = smoothstep(0.0, 1.0, col);

    gl_FragColor = vec4(col, 1.0);
}

How Volumetric Raymarching Works

Unlike billboard or flat-plane cloud techniques, volumetric clouds treat the sky as a participating medium — a 3D region filled with particles that absorb, scatter, and emit light. The renderer casts a ray from the camera through each pixel, and for every small step along that ray inside the cloud layer it asks two questions: "Is there cloud here?" and "How much light reaches this point from the sun?" The answers are accumulated to produce the final color.

The outer loop is the primary raymarch. We limit it to the altitude slab between CLOUD_LOWER and CLOUD_UPPER so we never waste steps in empty sky. At each step we evaluate the density function, and if density is above a threshold we fire a secondary ray toward the sun (the light march) to figure out how much the cloud self-shadows at that point.

// Primary raymarch pseudocode:
for each step along the view ray {
    sample density at current position
    if density > 0 {
        light_reaching_point = lightMarch(position, sunDirection)
        scattered += transmittance * light * (1 - e^(-density * stepSize))
        transmittance *= e^(-density * stepSize)   // Beer-Lambert
    }
}

The Density Function — Building Clouds from Noise

Real cloud renderers (like Guerrilla Games' system for Horizon Zero Dawn) sample 3D textures — a low-frequency Perlin-Worley shape texture and a high-frequency Worley detail texture. Because we cannot use external textures here, we replicate this two-tier approach with procedural value noise run through FBM.

The key insight is erosion: we generate a smooth base shape with low-frequency FBM, then subtract a higher-frequency detail noise to carve wispy edges into the cloud. A height-dependent mask shapes the vertical profile — dense in the middle of the layer, thinning toward the top and bottom — imitating how real cumulus clouds are sculpted by temperature gradients.

// Density construction:
float base   = fbm(position * 0.0003);       // large-scale shape
float detail = fbm(position * 0.0009);        // fine erosion
float height = smoothstep(0.0, 0.15, hf)     // bottom ramp-on
             * smoothstep(1.0, 0.6, hf);      // top ramp-off

float density = smoothstep(0.35, 0.7, base) * height;
density -= detail * 0.25;                     // carve edges
density = max(density, 0.0);

Beer-Lambert Absorption

As light travels through a medium, its intensity decreases exponentially with the optical depth traversed. This is the Beer-Lambert law:

// T = transmittance (fraction of light surviving)
// sigma = absorption coefficient (density * factor)
// ds = step distance

T *= exp(-sigma * ds);

We apply this law twice: once along the view ray (to fade the sky behind thick clouds) and once along the light ray (to compute self-shadowing). The self-shadow pass is what gives clouds their dark, moody undersides contrasted against bright sunlit tops.

Henyey-Greenstein Phase Function

In a real cloud, photons do not scatter equally in all directions. Water droplets are strongly forward-scattering — looking toward the sun through thin cloud edges produces a bright silver lining. The Henyey-Greenstein function models this anisotropy with a single parameter g that controls the asymmetry:

float henyeyGreenstein(float cosTheta, float g) {
    float g2 = g * g;
    return (1.0 - g2) /
           (4.0 * 3.14159 * pow(1.0 + g2 - 2.0 * g * cosTheta, 1.5));
}

// Dual-lobe blend for realism:
// g = 0.6  → strong forward scattering (silver linings)
// g = -0.3 → mild back scattering (soft glow when sun is behind you)
float phase = hg(cosTheta, 0.6) * 0.7
            + hg(cosTheta, -0.3) * 0.3;

A positive g (close to 1) concentrates scattering forward — this creates the signature silver lining effect. A small negative lobe adds gentle back-illumination. Blending the two lobes is a well-known trick from film VFX (first popularized in production by Wrenninge et al.) and is used in virtually every AAA cloud renderer.

Silver Linings and Edge Glow

In addition to the phase function, we add an explicit edge-detection term. When the sampled density is very low (meaning the ray is near the cloud boundary), we boost the luminance contribution. This enhances that unmistakable bright halo you see when the sun sits just behind a cumulus tower:

// Edge factor: bright at cloud boundaries, zero inside
float edgeFactor = smoothstep(0.02, 0.0, density) * 2.0;

// Added to the luminance accumulation:
luminance += sunLight * edgeFactor * lightReaching * 0.6;

Performance Considerations for Games

The shader above runs around 48 primary steps and 6 light steps per primary hit — workable for a fullscreen demo but too heavy for a shipped game at full resolution. Production engines use several tricks to bring the cost down:

Quarter-resolution rendering — Clouds are raymarched into a buffer at 1/4 screen resolution, then temporally upsampled using reprojection. Horizon Zero Dawn and Nubis (Guerrilla's second-generation system) rely heavily on this.

Blue-noise jittering — Instead of regular step intervals, each pixel's starting offset is jittered with a blue-noise pattern. Combined with temporal accumulation this hides banding and allows fewer steps.

Early-out with transmittance threshold — Once the accumulated transmittance drops below ~1%, the cloud is effectively opaque and we can stop marching. This saves work in dense cloud cores.

LOD-based step scaling — Distant clouds use larger step sizes and fewer noise octaves since fine detail is invisible at range.

With these optimizations, modern games achieve beautiful volumetric skies at sub-millisecond GPU cost per frame — a dramatic improvement over the multi-millisecond budgets of early implementations.

Going Further

To evolve this shader toward production quality, consider adding: weather map-driven coverage (a 2D texture controlling where clouds appear), Worley noise for more cellular, cauliflower-like shapes, atmospheric perspective (fade distant clouds toward the horizon color), and temporal reprojection to amortize cost across frames. Each of these is a natural extension of the architecture shown above.