From b2dcdbdb8072da3e96d3463a9f9cb8af2407c18b Mon Sep 17 00:00:00 2001 From: Ulrik Guenther Date: Wed, 17 Apr 2024 10:51:21 +0200 Subject: [PATCH] DeferredLighting: Fix energy conservervation for Oren-Nayar model --- .../backends/shaders/DeferredLighting.frag | 114 +++++++++++------- 1 file changed, 71 insertions(+), 43 deletions(-) diff --git a/src/main/resources/graphics/scenery/backends/shaders/DeferredLighting.frag b/src/main/resources/graphics/scenery/backends/shaders/DeferredLighting.frag index 6f6449769..35fffb620 100644 --- a/src/main/resources/graphics/scenery/backends/shaders/DeferredLighting.frag +++ b/src/main/resources/graphics/scenery/backends/shaders/DeferredLighting.frag @@ -70,29 +70,43 @@ layout(set = 6, binding = 0, std140) uniform ShaderParameters { int displayHeight; }; +float chi(float v) { + return v > 0 ? 1.0 : 0.0; +} -float GGXDistribution(vec3 normal, vec3 halfway, float roughness) { - float a = roughness*roughness; - float aSquared = a*a; - float NdotH = abs(dot(normal, halfway)); +float BeckmannDistribution(float NdotH, float roughness) { + float alpha2 = roughness*roughness; float NdotH2 = NdotH*NdotH; - float denom = ((NdotH2 * (aSquared - 1.0) + 1.0)); - return aSquared/(denom*denom*PI); + float denom = NdotH2 * alpha2 + (1.0 - NdotH2); + return (chi(NdotH) * alpha2)/(denom*denom*PI); } float GeometrySchlick(float NdotV, float roughness) { float r = (roughness + 1.0); float k = (r*r) / 8.0; - return (NdotV)/(NdotV*(1.0 - k) + k); + return NdotV/(NdotV*(1.0 - k) + k); +} + +float GGXPartialGeometry(vec3 N, vec3 V, vec3 H, float roughness) { + float VdotH2 = clamp(dot(V, H), 0.0, 1.0); + float Chi = chi(VdotH2/clamp(dot(V, N), 0.0, 1.0)); + + VdotH2 = VdotH2 * VdotH2; + float tan2 = (1.0 - VdotH2)/VdotH2; + + return (Chi * 2.0)/(1 + sqrt(1.0 + roughness * roughness * tan2)); } -float GeometrySmith(vec3 normal, vec3 view, vec3 light, float roughness) { - float NdotV = abs(dot(normal, view)); - float NdotL = abs(dot(normal, light)); +float V1(vec3 N, vec3 X, vec3 H, float roughness) { + float alpha2 = roughness * roughness; + float NdotX = clamp(dot(N, X), 0.0001, 1.0); + return 1.0 / (NdotX + sqrt(alpha2 + (1.0 - alpha2) * NdotX * NdotX)); +} - return GeometrySchlick(NdotV, roughness) * GeometrySchlick(NdotL, roughness); +float CookTorranceVisibility(vec3 N, vec3 V, vec3 H, vec3 L, float roughness) { + return V1(N, V, H, roughness) * V1(N, L, H, roughness); } vec3 FresnelSchlick(float cosTheta, vec3 F0) { @@ -174,7 +188,7 @@ vec3 viewFromDepth(float depth, vec2 texcoord, const mat4 invProjection) { float SinPhi(float Ndotw, vec3 w) { - float SinTheta = sqrt(max(0, 1 - Ndotw)); + float SinTheta = sqrt(max(0, 1 - Ndotw*Ndotw)); if(Ndotw < 0.0001) { return 0.0; @@ -184,7 +198,7 @@ float SinPhi(float Ndotw, vec3 w) { } float CosPhi(float Ndotw, vec3 w) { - float SinTheta = sqrt(max(0, 1 - Ndotw)); + float SinTheta = sqrt(max(0, 1 - Ndotw*Ndotw)); if(Ndotw < 0.0001) { return 1.0; @@ -197,9 +211,9 @@ vec2 alphabeta(float NdotL, float NdotV) { vec2 ab = vec2(0.0); if(abs(NdotL) > abs(NdotV)) { - ab = vec2(sqrt(max(0.0, 1 - NdotV)), sqrt(max(0.0, 1 - NdotL) / max(abs(NdotL), 0.0000001))); + ab = vec2(sqrt(max(0.0, 1 - NdotV*NdotV)), sqrt(max(0.0, 1 - NdotL*NdotL) / (NdotL*NdotL))); } else { - ab = vec2(sqrt(max(0.0, 1 - NdotL)), sqrt(max(0.0, 1 - NdotV) / max(abs(NdotV), 0.0000001))); + ab = vec2(sqrt(max(0.0, 1 - NdotL*NdotL)), sqrt(max(0.0, 1 - NdotV*NdotV) / (NdotV*NdotV))); } return ab; @@ -559,48 +573,62 @@ void main() // Oren-Nayar model for diffuse and Cook-Torrance for Specular else if(reflectanceModel == 0) { - float roughness = MaterialParams.r * PI / 2.0; + vec3 kD = vec3(1.0f); + const float roughness = MaterialParams.r * PI / 2.0; + const float metallic = MaterialParams.g; - float LdotV = max(dot(L, V), 0.0); - float NdotL = max(dot(L, N), 0.0); - float NdotV = max(dot(N, V), 0.0); + const float LdotV = max(dot(L, V), 0.0); + const float NdotL = max(dot(N, L), 0.0); + const float NdotV = max(dot(N, V), 0.0); + const float NdotH = max(dot(N, H), 0.0); - float sigma2 = roughness * roughness; - float A = 1.0 - sigma2 / (2.0 * (sigma2 + 0.33)); - float B = 0.45 * sigma2 / (sigma2 + 0.09); + if(metallic < 0.99f) { + // Oren-Nayar has sigma parameter in radians in [0, pi/2], we use [0, 1]. - vec2 ab = alphabeta(NdotL, NdotV); - float m = max(0, CosPhi(NdotL, L) * CosPhi(NdotV, V) + SinPhi(NdotL, L) * CosPhi(NdotV, V)); + const float sigma2 = roughness * roughness; + const float A = 1.0 - sigma2 / (2.0 * (sigma2 + 0.33)); + const float B = 0.45 * sigma2 / (sigma2 + 0.09); - float L1 = NdotL / PI * (A + B * m * ab.x * ab.y); + const vec2 cosTheta = vec2(clamp(dot(N, L), 0.0, 1.0), clamp(dot(N, V), 0.00001, 1.0)); + const vec2 cosThetaSq = cosTheta * cosTheta; - vec3 inputColor = intensity * emissionColor.rgb * Albedo.rgb * lightOcclusion; + const float sinTheta = sqrt((1 - cosThetaSq.x)*(1 - cosThetaSq.y)); - diffuse = inputColor * L1; + vec3 LP = normalize(L - cosTheta.x * N); + vec3 VP = normalize(V - cosTheta.y * N); + float cosPhi = clamp(dot(LP, VP), 0.0, 1.0); - if(Specular > 0.0 || MaterialParams.g > 0.0) { - float metallic = MaterialParams.g; - vec3 F0 = vec3(0.04); - F0 = mix(F0, Albedo.rgb, metallic); + float C = cosPhi * sinTheta/max(cosTheta.x, cosTheta.y); + float orenNayar = cosTheta.x * (A + B * C) / PI; - float roughness = MaterialParams.r; + vec3 Li = intensity * emissionColor.rgb * Albedo.rgb * lightOcclusion; - float NDF = GGXDistribution(N, H, roughness); - float G = GeometrySmith(N, V, L, roughness); - vec3 F = FresnelSchlick(max(dot(H, V), 0.0), F0); + diffuse = Li * orenNayar; + } - vec3 BRDF = (NDF * G * F)/max(4.0 * abs(NdotL) * abs(NdotV), 0.001); + // Cook-Torrance + // TODO: Have a material property for index of refraction + // float ior = 1.0; + // vec3 F0 = vec3(abs((1.0 - ior)/(1.0 + ior))); + // F0 *= F0; + vec3 F0 = vec3(0.04); + F0 = mix(F0, Albedo.rgb, metallic); - vec3 kS = F; - vec3 kD = (vec3(1.0) - kS); - kD *= 1.0 - metallic; + const float roughnessCookTorrance = clamp(MaterialParams.r, 0.001, 1.0); - vec3 radiance = intensity * emissionColor.rgb; - specular = (kD * Albedo.rgb / PI + BRDF) * radiance * NdotL; - } - } + vec3 fresnel = FresnelSchlick(max(dot(H, V), 0.0), F0); + float visibility = CookTorranceVisibility(N, V, H, L, roughnessCookTorrance); + float NDF = BeckmannDistribution(NdotH, roughnessCookTorrance); + vec3 BRDF = fresnel * visibility * NDF * NdotL; + vec3 kS = fresnel; + kD = (vec3(1.0) - kS) * (1.0 - metallic); + + vec3 radiance = intensity * emissionColor.rgb; + // combine Cook-Torrance and Oren-Nayar, conserving energy + lighting = (kD * diffuse + BRDF * radiance) * lightAttenuation; + } if(debugLights > 0) { if(debugLights == 3) {