float ld(const float dist) {
    return (2.0 * near) / (far + near - dist * (far - near));

}

float sqr(float x) {
    return x * x;
}

vec2 hash42v2(vec2 x) {

    vec2 x2 = fract(mod(x.xy, PI) * vec2(10.214, 9.637).xy);
    x2.xy += dot(x2.xy, x2.yx + vec2(38.549, 37.759));
    return fract((x2.xx + x2.yy) * x2);
}

float G1V(float dotNV, float k) {

    return 1.0 / (dotNV * (1.0 - k) + k);
}
vec3 GGX2(vec3 n, vec3 v, vec3 l, float r, vec3 F0) {
    r = max(pow(r, 2.5), 0.0001);

    vec3 h = l + v;
    float hn = inversesqrt(dot(h, h));

    float dotLH = clamp(dot(h, l) * hn, 0., 1.);
    float dotNH = clamp(dot(h, n) * hn, 0., 1.);
    float dotNL = clamp(dot(n, l), 0., 1.);
    float dotNHsq = dotNH * dotNH;

    float denom = dotNHsq * r - dotNHsq + 1.;
    float D = r / (PI * denom * denom);
    vec3 F = F0 + (1. - F0) * exp2((-5.55473 * dotLH - 6.98316) * dotLH);
    float k2 = .25 * r;

    return dotNL * D * F / (dotLH * dotLH * (1.0 - k2) + k2);
}

#ifdef LAB_REFLECTIONS
uniform sampler2D colortex4; // reflections

const vec3[] metalF0 = vec3[](vec3(0.78, 0.77, 0.74), // Iron
vec3(1.00, 0.90, 0.61), // Gold
vec3(1.00, 0.98, 1.00), // Aluminum
vec3(0.77, 0.80, 0.79), // Chrome
vec3(1.00, 0.89, 0.73), // Copper
vec3(0.79, 0.87, 0.85), // Lead
vec3(0.92, 0.90, 0.83), // Platinum
vec3(1.00, 1.00, 0.91)  // Silver
);
const vec3[] metalF82 = vec3[](vec3(0.74, 0.76, 0.76), vec3(1.00, 0.93, 0.73), vec3(0.96, 0.97, 0.98), vec3(0.74, 0.79, 0.78), vec3(1.00, 0.90, 0.80), vec3(0.83, 0.80, 0.83), vec3(0.89, 0.87, 0.81), vec3(1.00, 1.00, 0.95));

vec3 fresnelDielectric(float cosTheta, float n) {
    float gSq = sqr(n) + sqr(cosTheta) - 1.0;

    if (gSq < 0.0)
        return vec3(1.0); // Imaginary g => TIR

    float g = sqrt(gSq);
    float a = g - cosTheta;
    float b = g + cosTheta;

    return vec3(0.5 * sqr(a / b) * (1.0 + sqr((b * cosTheta - 1.0) / (a * cosTheta + 1.0))));
}
float linStep(float x, float low, float high) {
    return clamp((x - low) / (high - low), 0, 1);
}
vec3 fresnelLazanyi2019(float cosTheta, vec3 f0, vec3 f82) {
    vec3 a = 17.6513846 * (f0 - f82) + 8.1666666 * (1.0 - f0);
    float m = pow(1.0 - cosTheta, 5);
    return clamp(f0 + (1.0 - f0) * m - a * cosTheta * (m - m * cosTheta), 0, 1);
}
float precalc = far * 1.73205080757;

vec3 rayTrace(vec3 dir, vec3 position, float dither, float fresnel) {
    float quality = mix(1, 16, fresnel);

    vec3 clipPosition = toClipSpace3(position);
    float rayLength = ((dir.z * precalc + position.z) > -near) ? dir.z / (-near - position.z) : precalc;

    vec3 direction = normalize(toClipSpace3(position + dir * rayLength) - clipPosition);

    // get at which length the ray intersects with the edge of the screen
    vec3 maxLengths = (step(0.0, direction) - clipPosition) / direction;
    float mult = min(min(maxLengths.x, maxLengths.y), maxLengths.z);

    vec3 stepv = direction * mult / quality * vec3(RENDER_SCALE, 1.0);
    vec3 spos = clipPosition * vec3(RENDER_SCALE, 1.0) + stepv * dither;
    float tolerance = max(abs(stepv.z) * 3.0, 0.02 / (position.z * position.z)); // from DrDesten <3

    for (int i = 0; i < quality; i++) {
        float sp = texelFetch2D(depthtex0, ivec2(spos.xy / texelSize), 0).r;
        if (sp < spos.z && abs(tolerance - (spos.z - sp)) < tolerance) {
            return vec3(spos.xy, sp);
        }
        spos += stepv;

    }

    return vec3(1.1);
}
float OrenNayar(in vec3 l, in vec3 n, in vec3 v, float r) {

    float r2 = r * r;
    float a = 1.0 - 0.5 * (r2 / (r2 + 0.57));
    float b = 0.45 * (r2 / (r2 + 0.09));

    float nl = dot(n, l);
    float nv = dot(n, v);

    float ga = dot(v - n * nv, n - n * nl);

    return max(0.0, nl) * (a + b * max(0.0, ga) * sqrt((1.0 - nv * nv) * (1.0 - nl * nl)) / max(nl, nv));
}

// Calculates the Fresnel term using Schlick's approximation.
float FresnelSchlick(float cosTheta, float F0) {
    return F0 + (1.0 - F0) * pow(1.0 - cosTheta, 5.0);
}

// GGX Distribution function, part of the microfacet model for rough surfaces.
float GGXDistribution(vec3 N, vec3 H, float alpha) {
    float alphaSq = alpha * alpha;
    float NdotH = max(dot(N, H), 0.0);
    float denom = (NdotH * NdotH) * (alphaSq - 1.0) + 1.0;
    return alphaSq / (PI * denom * denom);
}
// Schlick-GGX geometry function.
float GeometrySchlickGGX(float NdotV, float alpha) {
    float alphaSq = alpha * alpha;
    float k = alphaSq / 2.0;
    return NdotV / (NdotV * (1.0 - k) + k);
}

// Geometry Smith function with GGX geometry function, accounts for shadowing and masking.
float GeometrySmith(vec3 N, vec3 V, vec3 L, float alpha) {
    float NdotV = max(dot(N, V), 0.0);
    float NdotL = max(dot(N, L), 0.0);
    float ggx1 = GeometrySchlickGGX(NdotV, alpha);
    float ggx2 = GeometrySchlickGGX(NdotL, alpha);
    return ggx1 * ggx2;
}

// Cook-Torrance BRDF function.
vec3 CookTorranceBRDF(vec3 N, vec3 V, vec3 L, vec3 F0, float roughness) {
    vec3 H = normalize(V + L);
    float D = GGXDistribution(N, H, roughness * roughness); // Convert roughness to alpha
    float G = GeometrySmith(N, V, L, roughness * roughness);
    float F = FresnelSchlick(max(dot(H, V), 0.0), dot(F0, vec3(0.333))); // Assuming F0 is the same for RGB

    // Specular reflection component
    float numerator = D * G * F;
    float denominator = 4.0 * max(dot(N, V), 0.0) * max(dot(N, L), 0.0) + 0.001; // Adding a small value to prevent division by zero
    vec3 specular = vec3(numerator / denominator);

    // Assuming no diffuse component for this example
    return specular;
}

vec3 computeReflections(vec4 spectex, float lightmapY, float roughness, vec3 f0, vec3 normal, vec3 fragpos, vec3 albedo, vec3 ambient) {
    vec4 reflection;
    vec3 reflections;
    float porosity = (spectex.z * 255 <= 64.5) ? step(0.0, spectex.z) : 0.0;
    float rainMult = sqrt(lightmapY) * wetness * (1.0 - (porosity * porosity));

    roughness = mix(roughness, 0.01, rainMult);
    f0 = mix(f0, vec3(0.02), rainMult);
    vec3 f82 = f0;

    if (f0.y > 229.5 * recip) {
        int index = int(int(spectex.y * 255) - 229);
        f0 = metalF0[index];
        f82 = metalF82[index];
    }

    float cosTheta = max(dot(normal, fragpos), 0.0);
    vec3 fresnel = fresnelDielectric(cosTheta, spectex.y);
    if (spectex.y > 229.5 * recip)
        fresnel = fresnelLazanyi2019(cosTheta, f0, f82);
    if (spectex.x > clamp(REFLECTION_CUTOFF, 0.1, 1.0)) {
        vec3 reflectedVector = vec3(0.0);
        vec3 scatteredVector = vec3(0.0);
        vec4 noise3 = blueNoise3(gl_FragCoord.xy);
        const int quality = 2;
        for (int i = 0; i < quality; i++) {
            int seed = (frameCounter % 40000) * quality + i;

            vec2 ij = fract(R2_samples3(seed) + noise3.rg);
            vec3 rayDir = normalize(cosineHemisphereSample2(ij));

            rayDir = TangentToWorld2(normal, rayDir);
            rayDir = mat3(gbufferModelView) * rayDir;
            reflectedVector = reflect(normalize(fragpos.xyz), normal);
            scatteredVector = reflect(normalize(rayDir), normal);
            vec3 roughReflectedWithRoughness = normalize(mix(reflectedVector, scatteredVector, roughness));
            //vec3 roughReflectedWithRoughness = CookTorranceBRDF(normal, sunPosition, rayDir, fresnel, roughness);

            vec3 rtPos = rayTrace(roughReflectedWithRoughness, fragpos.xyz, fract(seed / 1.6180339887 + noise3.r), 1);

            if (rtPos.z < 1) {
                vec3 previousPosition = reproject(mat3(gbufferModelViewInverse) * toScreenSpace(rtPos) + gbufferModelViewInverse[3].xyz, false);

                if (clamp(previousPosition.xy, 0, 1) == previousPosition.xy) {
                    reflection += vec4(texture(colortex4, previousPosition.xy).rgb, 1);
                }
            } else {

                reflection.rgb += ambient;
            }
        }
        reflection /= quality;
        reflection = mix(vec4(ambient, 1), reflection, reflection.a);
        reflections = (((reflection.rgb) * (fresnel * albedo + albedo * 0.5)) * clamp(fresnel, 0, 1));

    } else {
        reflections = ((ambient) * (fresnel * 0.01)) * clamp(fresnel, 0, 1);
    }
    return reflections;
}
#endif
