前言
因为前边讲了在Unity中如何实现一个BRDF——简单来说就是把公式抄到Shader里,于是这篇文章则主要从原理角度来讲一讲基于物理的BRDF公式为什么长那个样子。本篇文章主要是整理一下去年(2022)十二月底写的关于基于微平面理论的BRDF的一些东西,主要依据的是两篇经典论文:一篇是Cook-Torrance模型的A Reflectance Model for Computer Graphics[1982],另一篇则是GGX模型Microfacet Models for Refraction through Rough Surfaces[2007](这篇论文的作者里也有我的导师…真要说的话本篇文章内容的由来还有一些故事)。
原本计划中应该还有两三篇文章,但这可能是这几天(或者很长一段时间)的最后一篇了,连更几篇确实把我写累了(倒),而且因为一直忙着整理以前的东西,导致我这几天的作业也还没写,导师的项目也几天没动了(其实也是因为还有一些进度存货) 。
BRDF
双向反射分布函数(Bidirectional Reflectance Distribution Function,BRDF)描述了物体表面入射光与反射光的关系,具体从定义上来说是出射方向Radiance与入射方向Irradiance的比值,这里也简单补上辐射度量学里的相关定义:
- 辐射能量(Radiant Energy):电磁辐射的能量, Q Q Q
- 辐射通量/功率(Radiant Flux/Power):单位时间内通过某一区域的能量, Φ = d Q d t \Phi = \frac{dQ}{dt} Φ=dtdQ
- 辐射强度(Radiant Intensity):单位立体角内的辐射通量, I ( ω ) = d Φ d ω I(\omega) = \frac{d\Phi}{d\omega} I(ω)=dωdΦ
- 辐射照度(Irradiance):单位面积内的辐射通量, E ( x ) = d Φ ( x ) d A E(x) = \frac{d\Phi(x)}{dA} E(x)=dAdΦ(x)
- 辐射亮度(Radiance):单位立体角、单位投影面积内的辐射通量, L ( p , ω ) = d 2 Φ ( p , ω ) d ω d A c o s θ L\left( {p,\omega} \right) = \frac{d^{2}\Phi(p,\omega)}{d\omega dAcos\theta} L(p,ω)=dωdAcosθd2Φ(p,ω)
这时候我们再来看渲染方程,
L
i
(
p
,
ω
i
)
n
⋅
ω
i
d
ω
i
L_{i}\left( {p,\omega_{i}} \right)n \cdot \omega_{i}\mathbb{d}\omega_{i}
Li(p,ωi)n⋅ωidωi其实算的就是入射方向的Irradiance
d
E
(
ω
i
)
{d}E(\omega_{i})
dE(ωi),中间的点乘是因为入射的角度不同接收到的能量就不同。
而 f r ( p , ω i , ω o ) L i ( p , ω i ) n ⋅ ω i d ω i f_{r}\left( {p,\omega_{i},\omega_{o}} \right)L_{i}\left( {p,\omega_{i}} \right)n \cdot \omega_{i}\mathbb{d}\omega_{i} fr(p,ωi,ωo)Li(p,ωi)n⋅ωidωi就是从这个方向入射的Radiance经过物体表面后又有多少被分配至出射方向,于是这么一积分再加上自身自发光的Radiance便是我们要的出射方向的Radiance,这也就是BRDF的作用。
L o ( p , ω o ) = L e ( p , ω o ) + ∫ Ω f r ( p , ω i , ω o ) L i ( p , ω i ) n ⋅ ω i d ω i L_{o}\left( {p,\omega_{o}} \right) = L_{e}\left( {p,\omega_{o}} \right) + \int_{\Omega}^{~}{f_{r}\left( {p,\omega_{i},\omega_{o}} \right)L_{i}\left( {p,\omega_{i}} \right)n \cdot \omega_{i}\mathbb{d}\omega_{i}} Lo(p,ωo)=Le(p,ωo)+∫Ω fr(p,ωi,ωo)Li(p,ωi)n⋅ωidωi
微平面理论与Cook-Torrance
回到BRDF本身,基于物理的BRDF基本理论便是微平面理论,其认为物体的宏观表面由一个个微小的表面组成,光线在微平面上发生理想镜面反射或折射,微平面分布越粗糙,物体高光越模糊:
一个经典的基于微平面理论提出的BRDF模型就是Cook-Torrance模型,其BRDF项包含了漫反射
R
d
R_{d}
Rd和镜面反射
R
s
R_{s}
Rs两部分,可以看到公式结构与之前实现里提到的GGX公式基本相同,不同的地方在于
R
s
R_{s}
Rs项Cook-Torrance分母中为
π
\pi
π而GGX为
4
4
4,之后普遍认为这里取
4
4
4更加合理:
R
=
s
R
s
+
d
R
d
R = sR_{s} + dR_{d}
R=sRs+dRd
R s = F π D G ( N ⋅ L ) ( N ⋅ V ) R_{s} = \frac{F}{\pi}\frac{DG}{({\mathbf{N}} \cdot {\mathbf{L}})({\mathbf{N}} \cdot {\mathbf{V}})} Rs=πF(N⋅L)(N⋅V)DG
其中镜面反射描述了从物体表面反射出去的光,而漫反射本质上与次表面散射相同(在Cook-Torrance的论文中,认为漫反射也包含了当表面足够粗糙时光在物体微表面的多次弹射,因为从下文就可以看到,镜面反射中没有考虑到多次弹射的反射现象。而现在似乎有单独对其进行的研究,粗糙度改变时能量应该守恒),都是在进入物体表面后散射出去的光,不同的是相对于渲染像素的距离:如果光线在同一像素内进入并离开则被认为是漫反射,不同像素则认为是次表面散射。
对于漫反射,通常都采用归一化(归一化是为了保证能量守恒,当漫反射强度系数为1时,进入的能量等于出去的能量)的Lambert漫反射模型,模型虽然简单同时也没有考虑菲涅尔反射系数对入射光数量的影响,但效果很好,其公式如下:
R d = 1 π R_{d} = \frac{1}{\pi} Rd=π1
菲涅尔项F
菲涅尔项描述了一个现象,视线与物体表面夹角越小(越接近掠射角)反射越强,如近处的海面可以看到水面以下的东西,而远处的海面则反射出天空的颜色(实质上菲涅尔反射系数是描述了对于一份光,有多少进入了物体表面,又有多少在物体表面被直接反射)。对于同一种材质,在不同波长下,其强度可能不同,这也是为什么金属的高光会有颜色。
通常我们采用一个近似公式来进行描述,其中
F
0
F_{0}
F0为入射角为0时的菲涅尔反射率,通常金属较高,而非金属较低。注意宏观表面的菲涅尔反射效应是由微观表面菲涅尔反射产生的,在公式的表现上就是公式描述了微观表面(理想表面)下的菲涅尔反射(取
m
\mathbf{m}
m的原因),并通过推导成为了宏观BRDF公式的一部分:
F
S
c
h
l
i
c
k
(
i
,
m
)
=
F
0
+
(
1
−
F
0
)
(
1
−
(
i
⋅
m
)
)
5
F_{Schlick}\left( {{\mathbf{i}},{\mathbf{m}}} \right) = F_{0} + \left( 1 - F_{0} \right)\left( 1 - ({\mathbf{i}} \cdot {\mathbf{m}}) \right)^{5}
FSchlick(i,m)=F0+(1−F0)(1−(i⋅m))5
法线分布函数D
法线分布函数描述了微平面法线 m \mathbf{m} m的统计分布。 D ( m ) D({\mathbf{m}}) D(m)是一个密度函数,其定义如下:以 d ω m \mathbb{d}\omega_{m} dωm表示以m为中心的无限小的立体角, d A \mathbb{d}A dA表示一个无限小的宏观面积,则 D ( m ) d ω m d A D({\mathbf{m}})d\omega_{m}dA D(m)dωmdA表示法线位于立体角内的对应微表面的总面积。下边是一些性质:
- 0 ≤ D ( m ) ≤ ∞ 0 \leq D({\mathbf{m}}) \leq \infty 0≤D(m)≤∞
- 1 ≤ ∫ D ( m ) d ω m 1 \leq \int{D({\mathbf{m}})d\omega_{m}} 1≤∫D(m)dωm,微表面面积大于等于对应的宏观表面面积
- ( v ⋅ n ) = ∫ D ( m ) ( v ⋅ m ) d ω m \left( {{\mathbf{v}} \cdot {\mathbf{n}}} \right) = \int{D({\mathbf{m}})\left( {{\mathbf{v}} \cdot {\mathbf{m}}} \right)d\omega_{m}} (v⋅n)=∫D(m)(v⋅m)dωm,两者投影面积相同
- ∫ D ( m ) ( n ⋅ m ) d ω m = 1 \int{D({\mathbf{m}})\left( {{\mathbf{n}} \cdot {\mathbf{m}}} \right)d\omega_{m}} = 1 ∫D(m)(n⋅m)dωm=1,即上式中 v \mathbf{v} v与 n \mathbf{n} n重合的特殊情况,也就是归一化,保证能量守恒
常见如Phong分布、Beckmann分布、GGX(Trowbridge-Reitz)分布,其实都是在尝试用函数来近似拟合物体表面的法线分布情况,其公式分别如下(均做了归一化处理,主要就是为了保证能量守恒,其中Cook-Torranc采用的便是未归一化的Beckmann分布):
- Phong分布: D ( m ) = α + 2 2 π c o s α θ D({\mathbf{m}}) = \frac{\alpha + 2}{2\pi}{cos}^{\alpha}\theta D(m)=2πα+2cosαθ, c o s θ = ( n ⋅ m ) cos\theta = ({\mathbf{n}} \cdot {\mathbf{m}}) cosθ=(n⋅m)
- Beckmann分布: D ( m ) = 1 π α 2 c o s 4 θ e − ( t a n θ α ) 2 D({\mathbf{m}}) = \frac{1}{\pi\alpha^{2}{cos}^{4}\theta}e^{- {(\frac{tan\theta}{\alpha})}^{2}} D(m)=πα2cos4θ1e−(αtanθ)2, c o s θ = ( n ⋅ m ) cos\theta = ({\mathbf{n}} \cdot {\mathbf{m}}) cosθ=(n⋅m)
- GGX分布: D ( m ) = α 2 π c o s 4 θ ( α 2 + t a n 2 θ ) 2 D({\mathbf{m}}) = \frac{\alpha^{2}}{\pi{cos}^{4}\theta\left( \alpha^{2} + {tan}^{2}\theta \right)^{2}} D(m)=πcos4θ(α2+tan2θ)2α2, c o s θ = ( n ⋅ m ) cos\theta = ({\mathbf{n}} \cdot {\mathbf{m}}) cosθ=(n⋅m)
阴影遮蔽函数G
阴影遮蔽项描述了光线被物体微平面遮挡的情况。阴影指入射光方向的遮挡,遮蔽指出射光方向的遮挡,其具有的一些性质如下:
- 0 ≤ G ( i , o , m ) ≤ 1 0 \leq G({\mathbf{i}},{\mathbf{o}},{\mathbf{m}}) \leq 1 0≤G(i,o,m)≤1,小于等于1同样是保证能量守恒
- G ( i , o , m ) = G ( o , i , m ) G\left( {{\mathbf{i}},{\mathbf{o}},{\mathbf{m}}} \right) = G({\mathbf{o}},{\mathbf{i}},{\mathbf{m}}) G(i,o,m)=G(o,i,m),交换入射与出射方向遮挡系数不变
- G ( i , o , m ) = 0 G\left( {{\mathbf{i}},{\mathbf{o}},{\mathbf{m}}} \right) = 0 G(i,o,m)=0,当 ( i ⋅ m ) ( i ⋅ n ) ≤ 0 ({\mathbf{i}} \cdot {\mathbf{m}})({\mathbf{i}} \cdot {\mathbf{n}}) \leq 0 (i⋅m)(i⋅n)≤0或 ( o ⋅ m ) ( o ⋅ n ) ≤ 0 ({\mathbf{o}} \cdot {\mathbf{m}})({\mathbf{o}} \cdot {\mathbf{n}}) \leq 0 (o⋅m)(o⋅n)≤0,即微表面的背面在宏观表面正方向上不可见
GGX中采用了Smith阴影遮蔽函数,将双向阴影遮蔽
G
G
G近似为两个单向阴影项
G
1
G_1
G1的可分离乘积:
G
≈
G
1
(
i
,
m
)
G
1
(
o
,
m
)
G \approx G_{1}({\mathbf{i}},{\mathbf{m}})G_{1}({\mathbf{o}},{\mathbf{m}})
G≈G1(i,m)G1(o,m)
G
1
G_1
G1是从法线分布函数D推导得出的,GGX分布推导出Smith
G
1
G_1
G1如下:
G
1
(
v
,
m
)
=
2
1
+
1
+
α
2
t
a
n
2
θ
,
c
o
s
θ
=
(
v
⋅
n
)
G_{1}\left( {{\mathbf{v}},{\mathbf{m}}} \right) = \frac{2}{1 + \sqrt{1 + \alpha^{2}{tan}^{2}\theta}},cos\theta = ({\mathbf{v}} \cdot {\mathbf{n}})
G1(v,m)=1+1+α2tan2θ2,cosθ=(v⋅n)
GGX公式推导
有了上边这些项的定义,就可以利用微平面的理想镜面反射推导宏观粗糙表面BRDF,通过积分对应微面元的贡献计算宏观BRDF(当然,GGX其实提出的是一个完整的BSDF模型,这里只讲其BRDF部分):
f
r
(
i
,
o
,
n
)
=
∫
∣
i
⋅
m
i
⋅
n
∣
f
r
m
(
i
,
o
,
m
)
∣
o
⋅
m
o
⋅
n
∣
G
(
i
,
o
,
m
)
D
(
m
)
d
ω
m
f_{r}\left( {{\mathbf{i}},{\mathbf{o}},{\mathbf{n}}} \right) = \int{\left| \frac{{\mathbf{i}} \cdot {\mathbf{m}}}{{\mathbf{i}} \cdot {\mathbf{n}}} \right|f_{r}^{m}\left( {{\mathbf{i}},{\mathbf{o}},{\mathbf{m}}} \right)\left| \frac{{\mathbf{o}} \cdot {\mathbf{m}}}{{\mathbf{o}} \cdot {\mathbf{n}}} \right|G\left( {{\mathbf{i}},{\mathbf{o}},{\mathbf{m}}} \right)D({\mathbf{m}})d\omega_{m}}
fr(i,o,n)=∫
i⋅ni⋅m
frm(i,o,m)
o⋅no⋅m
G(i,o,m)D(m)dωm
其中
∣
i
⋅
m
i
⋅
n
∣
\left| \frac{{\mathbf{i}} \cdot {\mathbf{m}}}{{\mathbf{i}} \cdot {\mathbf{n}}} \right|
i⋅ni⋅m
,
∣
o
⋅
m
o
⋅
n
∣
\left| \frac{{\mathbf{o}} \cdot {\mathbf{m}}}{{\mathbf{o}} \cdot {\mathbf{n}}} \right|
o⋅no⋅m
为宏观表面与微面元之间转换的系数,这与之前提到的辐射度量学的相关定义有关,因为微平面的法线与宏观表面的法线可能并不相同。其中:
f
r
m
(
i
,
o
,
m
)
=
F
(
i
,
m
)
δ
ω
m
(
h
,
m
)
4
(
i
⋅
h
)
2
f_{r}^{m}\left( {{\mathbf{i}},{\mathbf{o}},{\mathbf{m}}} \right) = F({\mathbf{i}},{\mathbf{m}})\frac{\delta_{\omega_{m}}\left( {\mathbf{h}},{\mathbf{m}} \right)}{4\left( {\mathbf{i}} \cdot {\mathbf{h}} \right)^{2}}
frm(i,o,m)=F(i,m)4(i⋅h)2δωm(h,m)
δ ω m ( h r , m ) \delta_{\omega_{m}}\left( {\mathbf{h}}_{r},{\mathbf{m}} \right) δωm(hr,m)为狄拉克函数,当 h = m {\mathbf{h}} = {\mathbf{m}} h=m时为无穷大,否则为0,也就是只有镜面反射方向不为0,或者说只有 m \mathbf{m} m与 h \mathbf{h} h重合的微平面才会对宏观BRDF产生贡献(也能用来理解为什么 m \mathbf{m} m最后被换成了 h \mathbf{h} h),其性质如下(g为任意函数):
∫ Ω g ( o ) δ ω o ( s , o ) d ω o = { g ( s ) , i f s ∈ Ω 0 , 其他 \int_{\Omega}^{~}{g({\mathbf{o}})\delta_{\omega_{o}}\left( {{\mathbf{s}},{\mathbf{o}}} \right)d\omega_{o} = \left\{ \begin{matrix}{g({\mathbf{s}}),~if{\mathbf{s}} \in \Omega} \\ {0,其他} \\ \end{matrix} \right.} ∫Ω g(o)δωo(s,o)dωo={g(s), ifs∈Ω0,其他
这个公式实际上是经过一定推导得来的,因为考虑镜面反射,所以 i ⋅ h \mathbf{i} \cdot \mathbf{h} i⋅h和 o ⋅ h \mathbf{o} \cdot \mathbf{h} o⋅h相同,以及狄拉克函数为了进行积分然后换元时产生的雅可比行列式,于是才有了分母中的式子,对详细的数学推导感兴趣的话可以看原论文。
现在我们将刚刚提到的式子代入上方进行积分,就得到了我们常见的这个BRDF公式:
f
r
(
i
,
o
,
n
)
=
F
(
i
,
h
)
G
(
i
,
o
,
h
)
D
(
h
)
4
∣
i
⋅
n
∣
|
o
⋅
n
∣
f_{r}\left( {{\mathbf{i}},{\mathbf{o}},{\mathbf{n}}} \right) = \frac{F\left( {{\mathbf{i}},{\mathbf{h}}} \right)G\left( {{\mathbf{i}},{\mathbf{o}},{\mathbf{h}}} \right)D({\mathbf{h}})}{\left. 4\left| {{\mathbf{i}} \cdot {\mathbf{n}}} \right| \middle| {\mathbf{o}} \cdot {\mathbf{n}} \right|}
fr(i,o,n)=4∣i⋅n∣∣o⋅n∣F(i,h)G(i,o,h)D(h)