Writing in the front

This is a recent formula derivation assignment I did in a training camp of a large factory. This assignment prompted me to seriously learn the derivation of the Cook-Torrance BRDF formula in the PBR rendering equation and found that the water in PBR was much deeper than I expected. When I first started reading it, I was so confused by all the mathematical definitions that I didn’t know what the big guys were talking about. After gritty and careful combing, I gradually understood something, but I was still not sure if I understood the definition of the cubic differential of a microsurface correctly. Then I can understand the formulas of line distribution function and geometric occlusion term. Therefore, in the extra points, I analyzed the F,G,DF, G,DF, D functions of mirror reflection term in more detail when I analyzed the PBR implementation of open source engine. However, in the analysis of the realization of one of the formulas, a parameter was different from what I recognized from the multifarious materials, so I submitted the problem. However, a few days later I found another derivation, which was more clear. After reading, I instantly feel that MY questions are good without level…… But regret is no use already handed in.

Simple rendering equation derivation

Centered on a tiny area of dAdAdA on the surface of an object:

  • Let the unit hemisphere scope be Omega \Omega.

  • Let the Radiance of this area reflected from the outgoing direction WOW_OWO be Lo(ω O)L_o(\ omega_O)Lo(ω O).

  • Set Irradiance to be received within the hemispheric scope to E=LicosθiE = L_i COS \theta_iE=Licosθ I.

  1. The bidirectional reflection distribution function is the ratio of the two. Therefore, the ratio of the differential of the incoming Irradiance in the ω I \ omega_I ω I direction to the differential of the outgoing Radiance in the ω O \ omega_O ω O direction can be derived as follows:

f r ( Omega. i . Omega. o ) = d L o ( Omega. o ) L i ( Omega. i ) c o s Theta. i d Omega. i f_r(\omega_i,\omega_o) = \frac{dL_o(\omega_o)}{L_i(\omega_i){cos\theta_i}d\omega_i}
  1. The differential product of Radiance can be obtained:

L o ( Omega. o ) = ∫ H 2 f r ( Omega. i . Omega. o ) L i ( Omega. i ) c o s Theta. i d Omega. i L_o(\omega_o)=\int_{H^{2}} f_r(\omega_i,\omega_o){L_i(\omega_i){cos\theta_i}d\omega_i}

In Cook-Torrance BRDF, the proportion of light reflected by mirror and refracted is taken into account, and the refracted part is absorbed and the rest is bounced out of the object surface, which can be regarded as diffuse reflection term. Therefore, BRDF in the rendering equation can be divided into two terms:


f r ( Omega. i . Omega. o ) = k d f l a m b e r t + k s f c o o k – t o r r a n c e f_r(\omega_i,\omega_o)=k_df_{lambert}+k_sf_{cook-torrance}

Derivation of diffuse reflection term BRDF

The diffuse reflectance term is calculated by the Lambertian illumination model, in which it is assumed that the incident light is evenly distributed throughout the half sphere, so the conversion process of the incident Irradiance to the outgoing Radiance is independent of the incident direction, that is, both BRDF and LiL_iLi are constant terms:


L o ( Omega. o ) = f r L i ∫ H 2 c o s Theta. i d Omega. i L_o(\omega_o)= f_rL_i\int_{H^{2}}{{cos\theta_i}d\omega_i}

  1. The differential area over a sphere is defined as dA=rsinθdθrdϕdA =rsin \theta d\theta rd \phidA=rsinθdθrdϕ. The differential area divided by r2r^2r2 yields a differential solid Angle that defines dω=dA/r2=sinθdθdϕd\omega =dA/r ^2 =sin \theta d\theta d\phidω=dA/r2=sinθdθdϕ.

  2. The integral of cosθcos\thetacosθ in the hemispherical domain can be calculated by integrating all the differential solid angles to obtain the complete sphere, and the sphere integral formula is transformed into a double integral form. The integral result of cosθcos\thetacosθ in the hemispherical domain is:


∫ Θ d Omega. = ∫ 0 2 PI. ∫ 0 PI. s i n Theta. d Theta. d ϕ \int_{\Theta}d\omega = \int_{0}^{2\pi}\int_{0}^{\pi}{sin\theta} {d\theta} {d\phi}

∫ Ω c o s Theta. d Omega. = ∫ 0 2 PI. ∫ 0 PI. 2 c o s Theta. s i n Theta. d Theta. d ϕ = PI. \int_{\Omega}cos\theta d\omega = \int_{0}^{2\pi}\int_{0}^{\frac{\pi}{2}}{cos\theta sin\theta} {d\theta} {d\phi}=\pi

Lo=πfrLiL_o = \ PI f_rL_iLo=πfrLi, assuming that the incident light is uniformly incident in the hemispherical range, LiL_iLi is independent of the incident Angle, but considering the loss of refraction, LiL_iLi and LoL_oLo should have a loss ratio CCC, and the BRDF of diffuse reflection term is fr= C /πf_r = {c}/{\ PI} FR = C /π.

BRDF derivation of specular reflection term

In the Cook-Torrance lighting model, specular reflection on micro-surfaces needs to be considered as follows:

  • Fresnel item Fr (o omega, omega I) F_r (\ omega_o, \ omega_i) Fr (o omega, omega I) : mirror emission proportion and refraction occur.
  • Normal distribution function D(ωh)D(\omega_h)D(ωh) : Probability density function of all microsurfaces with normal direction HHH.
  • Geometric function G (o omega, omega I) G (\ omega_o, \ omega_i) G (o omega, omega I) : eliminating shadows shade in the light path.

Incidence Irradiance

  1. D(ωh)D(\omega_h)D(ωh) is defined as the area of all microsurfaces with normal directions of HHH per unit area, per unit solid Angle. The definition of analogy Radiance L = d2 Φ / dA an omega (d) L = {d ^ 2 \ Phi} / (\ d {omega dA ^ {\ perp})} L = d2 Φ / dA an omega (d), A(ωh)A(ωh)A(ωh)A(ωh)A(ωh) is similar to φ \Phi φ.

d 2 A ( Omega. h ) = D ( Omega. h ) d Omega. h d A d^2A(\omega_h) = D(\omega_h)d\omega_hdA
  1. From the definition of Irradiance, E= D φ /dAE = {d\Phi}/{dA}E= D φ /dA, the radiation flux of incident light on differential area dAdAdA can be obtained as follows:

d Φ = L i c o s Theta. i d A d\Phi = L_icos\theta_idA
  1. For all specular reflection with normal direction HHH, the direction of incidence and the direction of emission are determined. Therefore, for the radiation flux of Irradiance incident on micro-surface, d φ D \Phid φ, the following results can be obtained by differentiating it at the unit solid Angle of incidence direction and again at the unit solid Angle of normal direction HHH:

d 3 Φ i = L i ( Omega. i ) d Omega. i d 2 A an ( Omega. h ) d^3\Phi_i = L_i(\omega_i)d\omega_i d^2A^{\perp}(\omega_h)

= L i ( Omega. i ) d Omega. i c o s Theta. h D ( Omega. h ) d Omega. h d A = L_i (\omega_i) d\omega_i cos\theta_h D(\omega_h) d\omega_h dA

Emergent Irradiance

Ignoring Fresnel term and shadow occlusion term for the time being, it can be seen that there is no change in the radiation flux of incident and outgoing radiation:


d 3 Φ o = d 3 Φ i d^3\Phi_o = d^3\Phi_i
  1. From the definition of Radiance, we know that it is the radiation flux φ \Phi Phi twice the differential of Phi Phi on unit area and unit solid Angle. Then, we can obtain the second differential of Radiance on all units of solid Angle with normal direction HHH:

d L o ( Omega. o ) = d 3 Φ o d Omega. o d A an = L i ( Omega. i ) d Omega. i c o s Theta. h D ( Omega. h ) d Omega. h d A d Omega. o c o s Theta. o d A dLo(\omega_o) = \frac{d^3\Phi_o}{d\omega_o dA^{\perp}}=\frac{ L_i (\omega_i) d\omega_i cos\theta_h D(\omega_h) d\omega_h dA}{d\omega_o cos\theta_o dA }
  1. It can be obtained from the definition of BRDF:

f r ( Omega. o . Omega. i ) = d L o ( Omega. o ) d E i ( Omega. i ) = L i ( Omega. i ) d Omega. i c o s Theta. h D ( Omega. h ) d Omega. h d A L i ( Omega. i ) c o s Theta. i d Omega. i d Omega. o c o s Theta. o d A = c o s Theta. h D ( Omega. h ) d Omega. h c o s Theta. i d Omega. o c o s Theta. o f_r(\omega_o,\omega_i)= \frac{dL_o{(\omega_o)}}{dE_i(\omega_i)} =\frac{ L_i (\omega_i) d\omega_i cos\theta_h D(\omega_h) d\omega_h dA}{L_i(\omega_i){cos\theta_i}d\omega_id\omega_o cos\theta_o dA }= \frac{cos\theta_h D(\omega_h) d\omega_h}{cos\theta_id\omega_o cos\theta_o}
  1. We know that the Angle of the exit direction is twice the Angle of the half vector direction θo=2θh\theta_o=2\theta_hθo=2θh, The ratio of ϕd ω=sinθdθdϕd\omega =sin \theta d\theta d\phidω=sinθdθdϕ can be deduced from the definition of differential solid Angle:

d Omega. o = s i n 2 Theta. h d 2 Theta. h d ϕ = 4 s i n Theta. h c o s Theta. h d Theta. h d ϕ d\omega_o = sin2\theta_h d2\theta_h d\phi = 4sin\theta_h cos\theta_h d\theta_h d\phi

d Omega. h d Omega. o = 1 4 c o s Theta. h \frac{d\omega_h}{d\omega_o} = \frac{1}{4cos\theta_h}
  1. After the ratio result is substituted into BRDF, the following formula can be obtained:

f r ( Omega. o . Omega. i ) = D ( Omega. h ) 4 c o s Theta. i c o s Theta. o f_r(\omega_o,\omega_i)=\frac{ D(\omega_h)}{4cos\theta_i cos\theta_o}
  1. At this time to join Fresnel item Fr (o) omega F_r (\ omega_o) Fr (o) omega, geometric visible item G (o omega, omega I) G (, \ \ omega_o omega_i) G (o omega, omega I) and BRDF is available:

f r ( Omega. o . Omega. i ) = D ( Omega. h ) F r ( Omega. o ) G ( Omega. o . Omega. i ) 4 c o s Theta. o c o s Theta. i f_r(\omega_o,\omega_i)=\frac{D(\omega_h)F_r(\omega_o)G(\omega_o,\omega_i)}{4 cos\theta_o cos\theta_i}

The complete integral formula

Considering the spontaneous light of the object, after linearly mixing the BRDF of diffuse and specular reflection terms, multiply the corresponding diffuse and specular reflection coefficients kd, KSK_d, k_SKd and ks, the complete rendering equation can be obtained as follows:


L o ( p . Omega. o ) = L e ( p . Omega. r ) + ∫ Ω ( k d c PI. + k s D ( Omega. h ) F r ( Omega. o ) G ( Omega. o . Omega. i ) 4 c o s Theta. o c o s Theta. i ) L i ( p . Omega. i ) c o s Theta. i d Omega. i L_o(p,\omega_o)=L_e(p,\omega_r)+\int_{\Omega} {(k_d \frac{c}{\pi}+k_s\frac{D(\omega_h)F_r(\omega_o)G(\omega_o,\omega_i)}{4 cos\theta_o cos\theta_i})} {L_i(p,\omega_i){cos\theta_i}d\omega_i}

PBR implementation of open source engine MadDream3D

The open source engine MadDream3D implemented the SHAder code of PBR integral formula in the maddream.js file of the project:

void calLight(inout vec3 lcolor, inout vec3 spec, vec3 wo, vec3 dir, vec3 n,Light lt){
    vec3 wi = 1.0*vec3(dir);
    vec3 h = normalize(wi+wo);
    float nv = abs(dot(n,wo))+1e-5;
    float nl = clamp(dot(n,wi),0.0.1.0);
    / / float vh = clamp (dot settlement (h), 0.0, 1.0);
    float lh = clamp(dot(h,wi),0.0.1.0);
    float nh =clamp(dot(h,n),0.0.1.0);
    float roughFix = rough*rough;
    vec3 F = fresnelReflect(F0,lh);
    float D = distributionGGX(nh,roughFix);
    float G = geometryGGX(nv,nl,roughFix);
    vec3 kd = vec3(1.0) - F;
    kd *=1.0 - metal; 
    vec3 Fr = (D*G)*F;
    vec3 Fd = kd/PI;
    lcolor += Fd*nl*lt.diffuse;
    spec += D*G*F*nl*lt.diffuse;
}
Copy the code

Fresnel item

The Realization of Fresnel term FFF adopts Schlick approximation:


F = F 0 + ( 1 – F 0 ) ( 1 – c o s Theta. ) F = F_0+(1-F_0)(1-cos\theta)

F 0 = ( n 1 – n 2 n 1 + n 2 ) 2 F_0=(\frac{n_1-n_2}{n_1+n_2})^2

The code directly defines a fixed value for F0F_0F0, and does not calculate the refractive index of light in different media:

vec3 F0 = vec3(0.02.0.02.0.02);
/ /...
vec3 fresnelReflect(vec3 f0,float VoH){
    return f0 + (vec3(1.0) - f0) * pow(1.0 - VoH, 5.0);
}
Copy the code

Normal distribution function

The implementation of normal distribution function DDD adopts GGX distribution:


D = Alpha. 2 PI. ( ( n ⋅ m ) 2 ( Alpha. 2 – 1 ) + 1 ) 2 D = \frac{\alpha^2}{\pi((n \cdot m)^2(\alpha^2-1) +1)^2}

As for the value of α\alphaα, it is recommended to take the square of surface roughness α= R2 α= R ^2α= R2 in Disney principle coloring model by referring to MAO Xingyun’s predecessors. In the engine implementation, the author made a certain transformation of the denominator of the formula by factoring:

float distributionGGX(float nh,float a){
    float a2 = a*a;
    float f = (nh * a2 - nh) * nh + 1.0;
    float denom = PI*f*f;
    return a2/denom;
}
Copy the code

Geometry item

Schlick-ggx approximation model is adopted to realize the geometric term GGG. The problem of occlusion exists in both the observation direction VVV and the light source direction LLL, so it should be dealt with in the two directions respectively:


G = G 1 G 2 G = G_1G_2

G 1 = n ⋅ v ( n ⋅ v ) ( 1 – k ) + k . G 2 = n ⋅ l ( n ⋅ l ) ( 1 – k ) + k G_1=\frac{n \cdot v}{(n \cdot v)(1-k)+k},G_2=\frac{n \cdot l}{(n \cdot l)(1-k)+k}

k = ( R + 1 ) 2 8 k=\frac{(R+1)^2}{8}

A complete implementation of the Schlick-GGX approximation model is made in the engine:

float geometryGGX(float NoV,float NoL,float rough){
    float ggx2 = GeometrySchlickGGX(NoV, roughness);
    float ggx1 = GeometrySchlickGGX(NoL, roughness);
    return ggx1 * ggx2;
}
float GeometrySchlickGGX(float NdotV, float roughness)
{
    float r = (roughness + 1.0);
    float k = (r*r) / 8.0;
    float nom   = NdotV;
    float denom = NdotV * (1.0 - k) + k;
    return nom / denom;
}
Copy the code

Other items

Kd, KSK_D, k_SKD and ks in cook-Torrance BRDF formula are diffuse and specular reflection coefficients, but Fresnel term FrF_rFr itself means the probability of specular reflection.

According to the explanation obtained by consulting the data, for the unreflected light refracted into the object, part of it was consumed and part of it would bounce out again, which is the reason for the occurrence of diffuse reflection. The consumed part is the reflection of reflectivity.

Therefore, in the implementation of the engine, the specular reflection coefficient can be replaced by the Fresnel term FFF, and the diffuse reflection coefficient can be replaced by 1−F1-F1−F.

After checking, the diffuse reflection coefficient is multiplied by 1− Metal1-metal1 − Metal. This is because metal does not produce diffuse reflection, so the proportion of mirror reflection taken by metal should be subtracted.

⋅ Wi =cosθin \cdot w_i =cos \theta_in =cosθ I, Times the incident light Li(ω I)L_i(\omega_i)Li(ω I).

vec3 kd = vec3(1.0) - F;
kd *=1.0 - metal; 
vec3 Fr = (D*G)*F;
vec3 Fd = kd/PI;
lcolor += Fd*nl*lt.diffuse;
spec += D*G*F*nl*lt.diffuse;
Copy the code