// From http://graphics.cs.aueb.gr/graphics/research_illumination.html
// "Real-Time Diffuse Global Illumination Using Radiance Hints"
// paper and shader code

uniform float R_wcs = 10.;            // Rmax: maximum sampling distance (in WCS units)
uniform vec3 extents;
uniform mat4 RHMatrix;
uniform mat4 RSMMatrix;
uniform sampler2D dtex;
uniform sampler2D ctex;
uniform sampler2D ntex;

flat in int slice;
layout (location = 0) out vec4 SHRed;
layout (location = 1) out vec4 SHGreen;
layout (location = 2) out vec4 SHBlue;

vec3 resolution = vec3(32, 16, 32);
#define SAMPLES 16

vec4 SHBasis (const in vec3 dir)
{
    float   L00  = 0.282095;
    float   L1_1 = 0.488603 * dir.y;
    float   L10  = 0.488603 * dir.z;
    float   L11  = 0.488603 * dir.x;
    return vec4 (L11, L1_1, L10, L00);
}

vec4 DirToSh(vec3 dir, float flux)
{
    return SHBasis (dir) * flux;
}

void main(void)
{
    vec3 normalizedRHCenter = 2. * vec3(gl_FragCoord.xy, slice) / resolution - 1.;
    vec3 RHcenter = (RHMatrix * vec4(normalizedRHCenter * extents, 1.)).xyz;

    vec4 ShadowProjectedRH = RSMMatrix * vec4(RHcenter, 1.);

    vec3 RHCellSize = extents / resolution;
    vec2 RHuv = .5 * ShadowProjectedRH.xy / ShadowProjectedRH.w + .5;
    float RHdepth = .5 * ShadowProjectedRH.z / ShadowProjectedRH.w + .5;

    vec4  SHr = vec4(0.);
    vec4  SHg = vec4(0.);
    vec4  SHb = vec4(0.);

    int x = int(gl_FragCoord.x), y = int(gl_FragCoord.y);
    float phi = 30. * (x ^ y) + 10. * x * y;

    for (int i = 0; i < SAMPLES; i++)
    {
        // produce a new sample location on the RSM texture
        float alpha = (i + .5) / SAMPLES;
        float theta = 2. * 3.14 * 7. * alpha;
        float h = alpha;
        vec2 offset = h * vec2(cos(theta), sin(theta));
        vec2 uv = RHuv + offset * 0.01;

        // Get world position and normal from the RSM sample
        float depth = texture2D(dtex, uv).z;
        vec4 RSMPos = inverse(RSMMatrix) * (2. * vec4(uv, depth, 1.) - 1.);
        RSMPos /= RSMPos.w;
        vec3 RSMAlbedo = texture(ctex, uv).xyz;
        vec3 normal = normalize(2. * texture(ntex, uv).xyz - 1.);

        // Sampled location inside the RH cell
        vec3 offset3d = vec3(uv, 0);
        vec3 SamplePos = RHcenter + .5 * offset3d.xzy * RHCellSize;

        // Normalize distance to RSM sample
        float dist = distance(SamplePos, RSMPos.xyz) / R_wcs;
        // Determine the incident direction.
        // Avoid very close samples (and numerical instability problems)
        vec3 RSM_to_RH_dir = (dist <= 0.00) ? vec3(0.) : normalize(SamplePos - RSMPos.xyz);
        float dotprod = max(dot(RSM_to_RH_dir, normal.xyz), 0.);
        float factor = dotprod / (0.1 + dist * dist);

        vec3 color = RSMAlbedo.rgb * factor;

        SHr += DirToSh(RSM_to_RH_dir, color.r);
        SHg += DirToSh(RSM_to_RH_dir, color.g);
        SHb += DirToSh(RSM_to_RH_dir, color.b);
    }

    SHr /= 3.14159 * SAMPLES;
    SHg /= 3.14159 * SAMPLES;
    SHb /= 3.14159 * SAMPLES;

    SHRed = SHr;
    SHGreen = SHg;
    SHBlue = SHb;
}