Rainbow Chromatic Dispersion Algorithm in GLSL

This article details the algorithm for computing rainbow chromatic angular dispersion in GLSL shaders using geometric refractive sphere equations. By combining Snell's law of refraction inside a spherical raindrop, René Descartes' deflection angle model, and Cauchy's dispersion equation for wavelength-dependent refractive indices, real-time graphics pipelines can accurately simulate the optical phase function and color separation of primary and secondary rainbows without relying on precomputed lookup textures.

Geometric Optics of a Refractive Sphere

The path of a light ray interacting with a spherical water droplet is governed by Snell's law at spherical interfaces. When a parallel incident ray hits a droplet at an incidence angle \(i \in [0, \pi/2]\), it refracts into the droplet at an angle \(r\) determined by:

\[\sin(i) = \eta(\lambda) \sin(r)\]

where \(\eta(\lambda)\) is the relative refractive index of water for a specific wavelength \(\lambda\).

For a primary rainbow (one internal reflection, \(p = 1\)), the total deflection angle \(\theta\) relative to the incoming anti-solar vector is expressed as:

\[\theta(i) = 4r - 2i = 4 \arcsin\left(\frac{\sin(i)}{\eta(\lambda)}\right) - 2i\]

The caustic concentration (the rainbow angle) occurs at the minimum deflection where the derivative \(d\theta/di = 0\). Solving this stationary condition yields the critical incidence angle \(i_c\):

\[\cos^2(i_c) = \frac{\eta(\lambda)^2 - 1}{3}\]

Evaluating \(\theta(i_c)\) gives the exact peak scattering angle for a specific wavelength, typically hovering near \(42^\circ\) for red light and \(40.4^\circ\) for violet light.

Chromatic Dispersion via Cauchy's Formula

Water is dispersive; its refractive index decreases as the wavelength increases. GLSL calculates the refractive index per spectral wavelength using Cauchy's empirical equation:

\[\eta(\lambda) = A + \frac{B}{\lambda^2}\]

For liquid water across the visible spectrum (measured in micrometers, \(\mu\text{m}\)):

• \(A \approx 1.3199\)
• \(B \approx 0.0051 \, \mu\text{m}^2\)

Approximating standard RGB primaries at representative wavelengths:

• Red (\(\lambda_R = 0.650 \, \mu\text{m}\)): \(\eta_R \approx 1.332\)
• Green (\(\lambda_G = 0.532 \, \mu\text{m}\)): \(\eta_G \approx 1.338\)
• Blue (\(\lambda_B = 0.450 \, \mu\text{m}\)): \(\eta_B \approx 1.345\)

The GLSL Evaluation Algorithm

In a fragment shader, the algorithm evaluates the scattering angle between the view ray and the anti-solar direction, applies the spherical stationary caustic equation for each wavelength component, and shapes the intensity using an Airy caustic or Gaussian intensity falloff.

// Computes the primary rainbow peak deflection angle for a given refractive index
float descartesRainbowAngle(float eta) {
    float cos_ic = sqrt(clamp((eta * eta - 1.0) / 3.0, 0.0, 1.0));
    float ic = acos(cos_ic);
    float rc = asin(sin(ic) / eta);
    return 4.0 * rc - 2.0 * ic;
}

// Evaluates spectral intensity around the caustic peak
vec3 evaluateRainbowDispersion(vec3 viewDir, vec3 sunDir) {
    // Phase angle relative to anti-solar point
    float cosView = dot(normalize(viewDir), normalize(sunDir));
    float scatterAngle = acos(clamp(cosView, -1.0, 1.0));

    // Wavelengths in micrometers: Red, Green, Blue
    vec3 wavelengths = vec3(0.650, 0.532, 0.450);

    // Cauchy dispersion formula for water
    const float A = 1.3199;
    const float B = 0.0051;
    vec3 eta = A + (B / (wavelengths * wavelengths));

    // Calculate peak angles for RGB
    vec3 peakAngle;
    peakAngle.r = descartesRainbowAngle(eta.r);
    peakAngle.g = descartesRainbowAngle(eta.g);
    peakAngle.b = descartesRainbowAngle(eta.b);

    // Angular width of the caustic intensity profile
    float halfWidth = 0.015; // Radians (~0.85 degrees)

    // Intensity profile: steep falloff outside, softer internal tail (Airy approximation)
    vec3 deltaAngle = scatterAngle - peakAngle;
    vec3 intensity = exp(-pow(max(-deltaAngle, 0.0) / halfWidth, 2.0)) 
                   * exp(-max(deltaAngle, 0.0) / (halfWidth * 2.5));

    return clamp(intensity, 0.0, 1.0);
}

Computational Complexity and Integration

This algorithm runs entirely analytically per-pixel within \(O(1)\) constant time. By evaluating the deflection angle natively inside GLSL, the simulation dynamically accounts for atmospheric humidity, varying droplet refractive indices, and multi-order reflections (such as secondary rainbows by adjusting the deflection function to \(\theta_2(i) = \pi + 2i - 6r\)) without requiring large precomputed scattering phase tables.