Computing CIE xyY from Tristimulus Values in GLSL

This article explains how to compute the CIE xyY color representation from CIE XYZ tristimulus values using the OpenGL Shading Language (GLSL). Converting XYZ to xyY separates the chromaticity coordinates (\(x, y\)) from the luminance component (\(Y\)), which is essential for gamut mapping, color analysis, and tone mapping in rendering pipelines. Below is the direct mathematical formulation and an optimized, robust GLSL implementation that handles potential edge cases such as division by zero.

Mathematical Formulation

The CIE XYZ color space represents colors via three tristimulus values: \(X\), \(Y\), and \(Z\). The CIE xyY representation normalizes these values so that chromaticity is independent of luminance:

• \(x = \frac{X}{X + Y + Z}\)
• \(y = \frac{Y}{X + Y + Z}\)
• \(Y = Y\)

The third coordinate, \(z = \frac{Z}{X + Y + Z}\), is redundant because \(x + y + z = 1\) (meaning \(z = 1 - x - y\)). Therefore, only \(x\), \(y\), and the original luminance value \(Y\) are stored.

GLSL Implementation

When implementing the conversion in GLSL, total darkness (\(X = Y = Z = 0\)) causes a division by zero. To prevent NaN artifacts, a small epsilon must guard the denominator, or a conditional check can substitute a fallback chromaticity (typically the reference white point, such as D65 where \(x \approx 0.3127, y \approx 0.3290\)).

// Converts CIE XYZ tristimulus values to CIE xyY
vec3 XYZ_to_xyY(vec3 XYZ) {
    float sum = XYZ.x + XYZ.y + XYZ.z;
    
    // Avoid division by zero for black pixels
    if (sum < 1e-6) {
        // Fallback to D65 white point chromaticity with 0 luminance
        return vec3(0.3127, 0.3290, 0.0);
    }
    
    float invSum = 1.0 / sum;
    float x = XYZ.x * invSum;
    float y = XYZ.y * invSum;
    float Y = XYZ.y;
    
    return vec3(x, y, Y);
}

Branchless Alternative

In high-throughput fragment shaders where branching may impact performance, an epsilon clamp provides a branchless alternative:

vec3 XYZ_to_xyY_branchless(vec3 XYZ) {
    float sum = max(XYZ.x + XYZ.y + XYZ.z, 1e-6);
    return vec3(XYZ.x / sum, XYZ.y / sum, XYZ.y);
}

Inverting the Transformation

To return to tristimulus values for final screen rendering, the inverse transformation calculates \(X\) and \(Z\) from \(x\), \(y\), and \(Y\):

vec3 xyY_to_XYZ(vec3 xyY) {
    if (xyY.y < 1e-6) {
        return vec3(0.0);
    }
    
    float X = (xyY.x * xyY.z) / xyY.y;
    float Y = xyY.z;
    float Z = ((1.0 - xyY.x - xyY.y) * xyY.z) / xyY.y;
    
    return vec3(X, Y, Z);
}

Ensure that input color values are linearized and converted to the CIE XYZ color space via the appropriate \(3 \times 3\) color matrix (such as sRGB-to-XYZ) before executing these functions.