#if !defined HAMMON_DIFFUSE_LIB
#define HAMMON_DIFFUSE_LIB
    vec3 hammonDiffuse(in vec3 albedo, in vec3 n, in float nDotV, in float nDotL, in float nDotH, in float lDotV, in float roughness) {
        nDotH = abs(nDotH) + 1e-5;
        float facing = 0.5 + 0.5 * lDotV;
        float rough = nDotH <= 0.0 ? 0.0 : facing * (0.9 - 0.4 * facing) * ((0.5 + nDotH) / max(nDotH, 0.05));
        vec3 fresnel_v = 1.0 - fresnelFunction(nDotV, n, vec3(0.0));
        vec3 fresnel_l = 1.0 - fresnelFunction(nDotL, n, vec3(0.0));
        vec3 energyConservationFactor = vec3(1.0);
        uint N = HAMMON_FRESNEL_HEMISPHERE_DIRECTIONS;
        for(uint i = 0u; i < N; ++i) {
            //Probably would be better if I didn't use a loop, but it doesn't harm performance much at all so who cares.
            vec2 hammersley = hammersley2d(i, N);
            vec3 direction = hemisphereSample_uniform(hammersley.x, hammersley.y);
            energyConservationFactor -= fresnelFunction(direction.z, n, vec3(0.0)) / N;
        }
        vec3 smooth_v = (fresnel_l * fresnel_v) / energyConservationFactor;
        vec3 single = mix(smooth_v, vec3(rough), roughness) * rcp(pi);
        float multi = 0.1159 * roughness;

        return max((albedo * single + albedo * multi) * nDotL, 0.0);
    }
#endif