Ray-Sphere Intersection Quadratic Equation in GLSL

This article explains the mathematical formulation and GLSL implementation of the quadratic equation used to find the intersection points between a ray and a sphere. By combining the vector equations of a 3D line and a spherical surface, we derive the standard quadratic equation \(at^2 + bt + c = 0\), analyze its discriminant to evaluate hit conditions, and demonstrate how to implement the solution efficiently inside a GLSL shader.

Mathematical Formulation

A ray is defined parametrically by an origin point \(\mathbf{O}\) and a normalized direction vector \(\mathbf{D}\), parameterized by distance \(t\):

\[\mathbf{P}(t) = \mathbf{O} + t\mathbf{D}\]

A sphere is defined by its center point \(\mathbf{C}\) and radius \(r\). Any point \(\mathbf{P}\) on the surface of the sphere satisfies:

\[\|\mathbf{P} - \mathbf{C}\|^2 = r^2\]

Substituting the ray equation into the sphere equation yields:

\[\|(\mathbf{O} + t\mathbf{D}) - \mathbf{C}\|^2 = r^2\]

Let the vector from the sphere center to the ray origin be \(\mathbf{V} = \mathbf{O} - \mathbf{C}\). The equation expands via the dot product:

\[(t\mathbf{D} + \mathbf{V}) \cdot (t\mathbf{D} + \mathbf{V}) = r^2\]

\[t^2(\mathbf{D} \cdot \mathbf{D}) + 2t(\mathbf{D} \cdot \mathbf{V}) + (\mathbf{V} \cdot \mathbf{V}) - r^2 = 0\]

The Quadratic Equation

The expansion matches the standard quadratic form \(at^2 + bt + c = 0\), where:

• \(a = \mathbf{D} \cdot \mathbf{D}\) (which simplifies to \(1.0\) if \(\mathbf{D}\) is normalized)
• \(b = 2.0 \cdot (\mathbf{D} \cdot \mathbf{V})\)
• \(c = (\mathbf{V} \cdot \mathbf{V}) - r^2\)

The discriminant \(\Delta = b^2 - 4ac\) determines the nature of the intersection:

• \(\Delta < 0\): The ray misses the sphere (no real roots).
• \(\Delta = 0\): The ray is tangent to the sphere (one real root).
• \(\Delta > 0\): The ray enters and exits the sphere (two real roots).

When \(\Delta \ge 0\), the distances along the ray to the intersection points are:

\[t = \frac{-b \pm \sqrt{\Delta}}{2a}\]

To optimize computations in GLSL, shaders frequently substitute \(b = 2h\), where \(h = \mathbf{D} \cdot \mathbf{V}\). This cancels the factor of 4 in the discriminant and the factor of 2 in the denominator, resulting in:

\[\Delta' = h^2 - ac\] \[t = \frac{-h \pm \sqrt{\Delta'}}{a}\]

GLSL Implementation

Assuming the ray direction rayDir is normalized, \(a = 1.0\), which further simplifies the computation.

// Returns true if an intersection occurs, setting t to the nearest valid hit distance
bool intersectRaySphere(vec3 rayOrigin, vec3 rayDir, vec3 sphereCenter, float radius, out float t) {
    vec3 oc = rayOrigin - sphereCenter;
    
    // Coefficients (with normalized rayDir, a = 1.0)
    float h = dot(rayDir, oc);
    float c = dot(oc, oc) - radius * radius;
    
    float discriminant = h * h - c;
    
    // Ray misses the sphere
    if (discriminant < 0.0) {
        return false;
    }
    
    float sqrtD = sqrt(discriminant);
    
    // Check nearest hit point first
    float t0 = -h - sqrtD;
    float t1 = -h + sqrtD;
    
    if (t0 > 0.0) {
        t = t0;
        return true;
    } else if (t1 > 0.0) {
        t = t1; // Inside the sphere
        return true;
    }
    
    // The sphere is behind the ray origin
    return false;
}

This simplified form minimizes arithmetic operations, branch costs, and register pressure, making it optimal for ray marching and path tracing pipelines within fragment shaders.