/***********************************************/
/*           Copyright (c) 2025 Belmu          */
/*             All Rights Reserved             */
/***********************************************/

/*
    [Credits]:
        Jessie - Spectrum to XYZ function, CIE binary search function and CMF weights method (https://github.com/Jessie-LC)

    [References]:
        Measured CIE data (http://www.cvrl.org/)
*/

uniform sampler1D CIE_CMF_2012;
uniform sampler1D CIE_BT709_Basis;

bool inVisibleRange(float wavelength) { 
    return clamp(wavelength, 0.0, 440.0) == wavelength;
}

const float scale      = 683.368;
const float cieAverage = 113.042;

vec3 spectrumToXYZExactCIE2012(float spectrum, float wavelength) {
    float n = wavelength - 390.0;
    
    if (!inVisibleRange(n)) return vec3(0.0);

    int n0 = min(int(n) - 1, int(n));
    int n1 = min(int(n) - 1, n0 + 1);

    vec3 lo = texelFetch(CIE_CMF_2012, n0, 0).rgb;
    vec3 hi = texelFetch(CIE_CMF_2012, n1, 0).rgb;

    vec3 xyz = mix(lo, hi, fract(n)) * spectrum * scale;
        
    //Normalize the conversion
    return xyz * 441.0 / cieAverage;
}

vec3 spectrumToLinear(float spectrum, float wavelength, vec3 cmfWeights) {
    float wavelengthPdf = dot(cmfWeights, texelFetch(CIE_CMF_2012, int(wavelength - 390.0), 0).rgb / cieAverage);
    return fromXYZ(spectrumToXYZExactCIE2012(spectrum, wavelength) / wavelengthPdf);
}

// Based on Jessie's code, adapted it to my coding style
float binarySearchCIE(float lowBound, float highBound, float target, vec3 cmfWeights) {
    float midIndex = (lowBound + highBound) * 0.5;

    int lookupRes = 441 * 2;
    for (int x = 0; x < int(log2(float(lookupRes))); x++) {
        float coordinate = midIndex;
        int coordTweaked = int(coordinate * 441.0);

        float value = saturate(dot(cmfWeights, texelFetch(CIE_CMF_2012, coordTweaked + 441, 0).rgb));

             if (value < target) lowBound  = midIndex;
        else if (value > target) highBound = midIndex;
        else return midIndex;

        midIndex = (lowBound + highBound) * 0.5;
    }
    return midIndex;
}

float srgbToSpectrum(vec3 color, float wavelength) {
    return dot(color, texelFetch(CIE_BT709_Basis, int(clamp(wavelength - 390.0, 0.0, 390.0)), 0).rgb);
}

vec3 estimateCMFWeights(sampler2D tex, vec2 coords) {
    vec3 history   = vec3(0.0);
    const int size = 4;

    for (int x = -size; x <= size; x++) {
        for (int y = -size; y <= size; y++) {
            history += texture(tex, coords + vec2(x, y) * texelSize).rgb;
        }
    }
    
    return mix(saturate(history / dot(vec3(1.0), history)), vec3(1.0 / 3.0), history == vec3(0.0) ? 1.0 : 0.2);
}
