Fix: EEVEE: runtime-generated BSDF LUT does not match the precomputed LUT

* `RUNTIME_LUT_CREATION` is disabled by default, but if set to `true`,
the new `bxdf_lut_frag.glsl` should generate the LUT that is submitted
in this commit. This new LUT is very close to the previously stored one.
* Generate the LUT by Monte-Carlo sampling with \(I = 1/N\sum\frac{f}{p}\),
except that the samples are not random, but mapped from a regular grid,
following previous approach to avoid interpolation artefacts caused by
noise.
* Added comments to make the computed quantities and the LUT usage more
clear.
* Glass with `IOR < 1` is now slightly darker due to single-scattering.
The lost energy will be recovered in a future commit.
* The special cases of `IOR == 0`, `roughness == 0` and `roughness == 1`
are handled well during microfacet sampling and evaluation, no need to
clamp.
* When `IOR == 1`, in theory BRDF is zero and BTDF is one, but we do not
make it a special case to allow smooth transition.

Pull Request: https://projects.blender.org/blender/blender/pulls/111632

source/blender/draw/engines/eevee/shaders/bsdf_lut_frag.glsl

@@ -5,55 +5,52 @@
 #pragma BLENDER_REQUIRE(common_utiltex_lib.glsl)
 #pragma BLENDER_REQUIRE(bsdf_sampling_lib.glsl)
 
+/* Generate BRDF LUT following "Real shading in unreal engine 4" by Brian Karis
+ * https://cdn2.unrealengine.com/Resources/files/2013SiggraphPresentationsNotes-26915738.pdf
+ * Parametrizing with `x = roughness` and `y = sqrt(1.0 - cos(theta))`.
+ * The result is interpreted as: `integral = f0 * scale + f90 * bias`. */
 void main()
 {
   /* Make sure coordinates are covering the whole [0..1] range at texel center. */
-  float y = floor(gl_FragCoord.y) / (LUT_SIZE - 1);
   float x = floor(gl_FragCoord.x) / (LUT_SIZE - 1);
+  float y = floor(gl_FragCoord.y) / (LUT_SIZE - 1);
+
+  /* Squaring for perceptually linear roughness, see [Physically Based Shading at Disney]
+   * (https://media.disneyanimation.com/uploads/production/publication_asset/48/asset/s2012_pbs_disney_brdf_notes_v3.pdf)
+   * Section 5.4. */
+  float roughness = x * x;
+  float roughness_sq = roughness * roughness;
 
   float NV = clamp(1.0 - y * y, 1e-4, 0.9999);
-  float a = x * x;
-  float a2 = clamp(a * a, 1e-4, 0.9999);
-
   vec3 V = vec3(sqrt(1.0 - NV * NV), 0.0, NV);
 
   /* Integrating BRDF */
-  float brdf_accum = 0.0;
-  float fresnel_accum = 0.0;
+  float scale = 0.0;
+  float bias = 0.0;
   for (float j = 0.0; j < sampleCount; j++) {
     for (float i = 0.0; i < sampleCount; i++) {
       vec3 Xi = (vec3(i, j, 0.0) + 0.5) / sampleCount;
       Xi.yz = vec2(cos(Xi.y * M_2PI), sin(Xi.y * M_2PI));
 
       /* Microfacet normal */
-      vec3 H = sample_ggx(Xi, a, V);
+      vec3 H = sample_ggx(Xi, roughness, V);
       vec3 L = -reflect(V, H);
       float NL = L.z;
 
       if (NL > 0.0) {
-        float NH = max(H.z, 0.0);
-        float VH = max(dot(V, H), 0.0);
+        /* Assuming sample visible normals, `weight = brdf * NV / (pdf * fresnel).` */
+        float weight = bxdf_ggx_smith_G1(NL, roughness_sq);
 
-        float G1_v = G1_Smith_GGX_opti(NV, a2);
-        float G1_l = G1_Smith_GGX_opti(NL, a2);
-        /* See G1_Smith_GGX_opti for explanations. */
-        float G_smith = 4.0 * NV * NL / (G1_v * G1_l);
+        /* Schlick's Fresnel. */
+        float s = pow(1.0 - saturate(dot(V, H)), 5.0);
 
-        float brdf = (G_smith * VH) / (NH * NV);
-
-        /* Follow maximum specular value for principled bsdf. */
-        const float specular = 1.0;
-        const float eta = (2.0 / (1.0 - sqrt(0.08 * specular))) - 1.0;
-        float fresnel = F_eta(eta, VH);
-        float Fc = F_color_blend(eta, fresnel, vec3(0)).r;
-
-        brdf_accum += (1.0 - Fc) * brdf;
-        fresnel_accum += Fc * brdf;
+        scale += (1.0 - s) * weight;
+        bias += s * weight;
       }
     }
   }
-  brdf_accum /= sampleCount * sampleCount;
-  fresnel_accum /= sampleCount * sampleCount;
+  scale /= sampleCount * sampleCount;
+  bias /= sampleCount * sampleCount;
 
-  FragColor = vec2(brdf_accum, fresnel_accum);
+  FragColor = vec2(scale, bias);
 }

source/blender/draw/engines/eevee/shaders/btdf_lut_frag.glsl

@@ -5,12 +5,15 @@
 #pragma BLENDER_REQUIRE(common_utiltex_lib.glsl)
 #pragma BLENDER_REQUIRE(bsdf_sampling_lib.glsl)
 
+/* Generate BSDF LUT for `IOR < 1`. Returns the integrated BTDF and BRDF, multiplied by the cosine
+ * foreshortening factor. */
 void main()
 {
+  /* Make sure coordinates are covering the whole [0..1] range at texel center. */
   float x = floor(gl_FragCoord.x) / (LUT_SIZE - 1.0);
   float y = floor(gl_FragCoord.y) / (LUT_SIZE - 1.0);
 
-  float ior = clamp(sqrt(x), 0.05, 0.999);
+  float ior = sqrt(x);
   /* ior is sin of critical angle. */
   float critical_cos = sqrt(1.0 - saturate(ior * ior));
 
@@ -22,67 +25,48 @@ void main()
   y += critical_cos;
   float NV = clamp(y, 1e-4, 0.9999);
 
-  float a = z_factor * z_factor;
-  float a2 = clamp(a * a, 1e-8, 0.9999);
+  /* Squaring for perceptually linear roughness, see [Physically Based Shading at Disney]
+   * (https://media.disneyanimation.com/uploads/production/publication_asset/48/asset/s2012_pbs_disney_brdf_notes_v3.pdf)
+   * Section 5.4. */
+  float roughness = z_factor * z_factor;
+  float roughness_sq = roughness * roughness;
 
   vec3 V = vec3(sqrt(1.0 - NV * NV), 0.0, NV);
 
-  /* Integrating BTDF */
-  float btdf_accum = 0.0;
-  float fresnel_accum = 0.0;
+  /* Integrating BSDF */
+  float btdf = 0.0;
+  float brdf = 0.0;
   for (float j = 0.0; j < sampleCount; j++) {
     for (float i = 0.0; i < sampleCount; i++) {
       vec3 Xi = (vec3(i, j, 0.0) + 0.5) / sampleCount;
       Xi.yz = vec2(cos(Xi.y * M_2PI), sin(Xi.y * M_2PI));
 
       /* Microfacet normal. */
-      vec3 H = sample_ggx(Xi, a2, V);
+      vec3 H = sample_ggx(Xi, roughness, V);
+      float fresnel = F_eta(ior, dot(V, H));
 
-      float VH = dot(V, H);
-
-      /* Check if there is total internal reflections. */
-      float fresnel = F_eta(ior, VH);
-
-      fresnel_accum += fresnel;
-
-      float eta = 1.0 / ior;
-      if (dot(H, V) < 0.0) {
-        H = -H;
-        eta = ior;
+      /* Reflection. */
+      vec3 R = -reflect(V, H);
+      float NR = R.z;
+      if (NR > 0.0) {
+        /* Assuming sample visible normals, accumulating `brdf * NV / pdf.` */
+        brdf += fresnel * bxdf_ggx_smith_G1(NR, roughness_sq);
       }
 
-      vec3 L = refract(-V, H, eta);
-      float NL = -L.z;
-
-      if ((NL > 0.0) && (fresnel < 0.999)) {
-        float LH = dot(L, H);
-
-        /* Balancing the adjustments made in G1_Smith. */
-        float G1_l = NL * 2.0 / G1_Smith_GGX_opti(NL, a2);
-
-        // btdf = abs(VH*LH) * (ior*ior) * D * G(V) * G(L) / (Ht2 * NV)
-        // pdf = (VH * abs(LH)) * (ior*ior) * D * G(V) / (Ht2 * NV)
-        float btdf = G1_l * abs(VH * LH) / (VH * abs(LH));
-
-        btdf_accum += btdf;
+      /* Refraction. */
+      vec3 T = refract(-V, H, ior);
+      float NT = T.z;
+      /* In the case of TIR, `T == vec3(0)`. */
+      if (NT < 0.0) {
+        /* Assuming sample visible normals, accumulating `btdf * NV / pdf.` */
+        btdf += (1.0 - fresnel) * bxdf_ggx_smith_G1(NT, roughness_sq);
       }
     }
   }
-  btdf_accum /= sampleCount * sampleCount;
-  fresnel_accum /= sampleCount * sampleCount;
-
-  if (z_factor == 0.0) {
-    /* Perfect mirror. Increased precision because the roughness is clamped. */
-    fresnel_accum = F_eta(ior, NV);
-  }
-
-  if (x == 0.0) {
-    /* Special case. */
-    fresnel_accum = 1.0;
-    btdf_accum = 0.0;
-  }
+  btdf /= sampleCount * sampleCount;
+  brdf /= sampleCount * sampleCount;
 
   /* There is place to put multi-scatter result (which is a little bit different still)
    * and / or lobe fitting for better sampling of. */
-  FragColor = vec4(btdf_accum, fresnel_accum, 0.0, 1.0);
+  FragColor = vec4(btdf, brdf, 0.0, 1.0);
 }