环境光照

上帝说要有光，于是就有了光，在渲染中，也是一样，物体要被看见，必须要有光。

上图是计算机图形学领域里一篇经典论文的“预告图(teaser image)”，图中展示了真实渲染需要解决的问题：

直接光照：直接从光源处发射过来的光源，然后经过物体反射，最终被观察者看到
间接光照：物体接收到光源后，也会反射一部分光源
间接阴影：物体反射出来的间接光照产生的阴影

如下图左图所示，light 是直接光照打到某个物体后发射过来的光线，view 是观察方向，这个间接光打到点 p 后，会向附近反射光线，然后反射到观察方向的光线就会被我们看到，真实世界中会有无数根这样的光照，因此环境光的问题就是找到所有这些光线，并求出这些光线经过 p 点后反射到观察方向的光线。

渲染方程

\[L_o(\boldsymbol{x},\omega_o) = L_e(\boldsymbol{x}, \omega_o) + \int_{H^2}f_r(\boldsymbol{x}, \omega_o,\omega_i)L_i(\boldsymbol{x}, \omega_i)cos\theta_i d\omega_i \]

其中:

\[\begin{aligned} \overbrace{L_0(\boldsymbol{x},\omega_o)}^{outgoing} &= \overbrace{L_e(\boldsymbol{x}, \omega_0)}^{emitted} \\ &+ \overbrace{\int_{H^2}f_r(\boldsymbol{x}, \omega_o,\omega_i)L_i(\boldsymbol{x}, \omega_i)cos\theta_i d\omega_i}^{reflected} \end{aligned} \]

\(\omega_o\)：观察方向
\(\boldsymbol{x}\)：渲染点

IBL（Image-Based Lighting）

我们经常会用到 Cubemap 来做环境光照，贴图里存储了环境光照数据，如下图：

那么我们计算每一个像素的环境光照时，环境光照来自上半球，因此需要对上半球方向上的环境贴图进行随机采样。

这样实时渲染时计算量就非常大，那 IBL 的思想是将一部分耗时的计算提前预计算

兰伯特模型（ Lambertian）

首先考虑简单的光照模型，一束光线打到物体上，假设物体是漫反射（Diffse），则光线会向四周均匀发射，然后我们从观察方向看到该方向上被发射的光线，在兰伯特模型中，漫反射的光线只跟入射方向有光，跟观察方向无关。

兰伯特光照模型是一种基于观察得到的经验模型，不符合真实物理情况。

首先定义一些变量：

\(\boldsymbol{v}\)：观察方向
\(\boldsymbol{n}\)：法线方向
\(\boldsymbol{l}\)：光线方向
其他变量：color 等

入射光照能量跟法线方向公式：

\[cos\theta = \boldsymbol{l} \cdot \boldsymbol{n} \]

垂直的时候能量最多，\(\theta\) 角度越大，能量越少

光线能量跟距离的关系，光线从发射点出发，呈球状分散开来，能量会均匀的分布在球面上，因此能量跟距离的平方成反比

\(Intensity\) : 光照强度

最后得出计算公式：

\[L_d = k_d \left(\frac{I}{r^2}\right)max(0, \boldsymbol{n} \cdot \boldsymbol{l}) \]

\(L_d\)：观察发现漫反射光线
\(K_d\)：漫反射系数，实际应用时是物体纹理贴图
\(max(0, \boldsymbol{n} \cdot \boldsymbol{l})\)：\(\boldsymbol{n} \cdot \boldsymbol{l}\) 可能为负数，

渲染结果：

Blinn-Phong 模型

实际上，光线反射除了漫反射外，有些物体（例如金属）会有高光的效果，反射光线几种在某一部分区域，（如下图黄色区域），因此 Phong 在 Lambort 模型基础上增加了对高光的计算。

\[L_s = k_s \left(\frac{I}{r^2}\right)max(0, \boldsymbol{R} \cdot \boldsymbol{l})^p \]

光高项使用了一个 \(pow\) 项来模拟，观察方向 \(\boldsymbol{v}\) 跟反射方向 \(\boldsymbol{R}\) 越近，高光越强。

Bling-Phong 改进点是使用 半程向量（\(\boldsymbol{h}\)） 来计算

Phong 模型中需要计算反射向量 \(\boldsymbol{R}\)，Bling-Phong 直接通过光照方向跟观察方向得出半程向量 \(\boldsymbol{h}\)。

\[L_s = k_s \left(\frac{I}{r^2}\right)max(0, \boldsymbol{n} \cdot \boldsymbol{h})^p \]

\(k_s\)：高光系数
\(p\)：高光指数参数

最后增加一个 Ambient 环境光项

\[L_a = k_a I_a \]

将所有光线累加起来，得出 Blinn-Phong 光照公式：

\[\begin{aligned} L &= L_a + L_d + L_s \\ &= k_a I_a + k_d \left(\frac{I}{r^2}\right)max(0, \boldsymbol{n} \cdot \boldsymbol{l}) + k_s \left(\frac{I}{r^2}\right)max(0, \boldsymbol{n} \cdot \boldsymbol{h})^p \end{aligned} \]

上面的方法都是理想化的经验模型，忽略了各种能量损失，因此如果需要进一步模拟真实的光照，就需要基于物理的光照方法了。

辐射度量学（Radiometry）

我们首先要把光照通过数学表达出来，这里引入辐射度量学的一些概念，首先介绍立体角

立体角

弧度：与圆心角与其对应的圆弧长度的关系

\[\theta = \frac{l}{r} \]

整个圆有 \(2\pi\) 的弧度

立体角：用来衡量三维球体的圆心角（宏观符号 \(\Omega\)，微分符号 \(\omega\)）

\[\Omega = \frac{A}{r^2} \]

整个球有 \(4\pi\) 的立体角，球面积为： \(S = 4\pi r^2\)，则立体角为 \(\Omega = \frac{S}{r^2} = 4 \pi\)

则在球上，对立体角的微分公式有：

\[\begin{aligned} dA &= (rd\theta)(rsin\theta d \omega \phi) \\ &= r^2sin\theta d\theta d\phi \end{aligned} \]

\[d\omega = \frac{dA}{r^2} = sin\theta d\theta d\phi \]

辐射通量 Radiant Flux

光照其实是一种能量（energy），能量的单位是焦耳（Q），但是图形学中我们一般用辐射通量或者辐射功率来描述光照，它是表示单位时间内的能量（W）

\[\Phi = \frac{dQ}{dt} \]

辐射强度 Radiant Intensity

辐射强度是表示：单位立体角的通过的辐射通量，用 \(I\) 表示，单位是坎德拉（candela）

\[I = \frac{d\Phi}{dw} \]

dw: 单位立体角

辐照度（Irradiance）

到达单位面积的辐射通量

\[E = \frac{d\Phi}{dA^{\perp}} \]

要投影到面积法线方向

有了辐照度后，我们就可以应用辐照度来计算光照能量了：

以及之前提到的距离跟光能量的关系：

随着距离 \(r\) 变大，\(I\) 始终保持不变（单位立体角不变），但是因为球面积增大，辐照度 \(E\) 减小

辐射率(Radiance)

通过单位面积单位立体角的辐射通量，光线从 dA 面往 \(w\) 方向辐射的光线能量

表面发出光线的辐射率（Exiting Radiance）

\[\begin{aligned} L &= \frac{dI}{dA^{\perp}}=\frac{d\frac{d\Phi}{d\omega}}{dAcos\theta}=\frac{d^2\Phi}{d\omega dAcos\theta} \\ &= \frac{dM}{d\omega} \end{aligned} \]

光源射向表面的辐射率（Incident Radiance）

\[L = \frac{dE}{d\omega} = \frac{d\frac{d\Phi}{dA^{\perp}}}{dw}=\frac{d^2\Phi}{d\omega dAcos\theta} \]

光线从 \(\omega\) 方向射入或者射出，计算辐射率时，都要投影到该表面 \(dA\) 的法线方向

名称	符号	单位	公式	解释
辐射能量(Radiant energy)	\(Q\)	焦耳(\(J\))	-	电磁辐射能量
辐射通量(Radiant Flux)	\(\Phi\)	瓦(\(W\))	\(\Phi = \frac{dQ}{dt}\)	单位时间辐射的能量，也叫辐射功率(Radiant Power)或通量(Flux)
辐照度(Irradiance)	\(E\)	瓦/平方米(\(W/m^2\))	\(E = \frac{d\Phi}{dA ^{\perp}}\)	到达单位面积的辐射通量
辐射度(Radiosity)	\(M\)	瓦/平方米(\(W/m^2\))	\(M = \frac{d\Phi}{dA^{\perp}}\)	离开单位面积的辐射通量，也叫辐出度、辐射出射度（Radiant Existance）
辐射强度(Radiant Intensity)	\(I\)	瓦/立体弧度(\(W/sr\))	\(I = \frac{d\Phi}{dw}\)	通过单位立体角的辐射通量
辐射率(Radiance)	\(L\)	瓦/平方米立体弧度(\(W/m^2sr\))	\(L = \frac{d^2\Phi}{dwdA^{\perp}}\)	通过单位面积单位立体角的辐射通量

\(dA^{\perp}=dAcos\theta\)：\(dA\) 法线方向上的投影量

能量积分

一根角度为 \(\omega\) 的光线射向点 \(p\) 的光线能量为：

\[\begin{aligned} L_i(p, \omega) &= \frac{dE(p, \omega)}{cos \theta d \omega} \\ \\ dE(p, \omega) &= L_i(p, \omega)cos \theta d \omega \end{aligned} \]

则对整个半球求积分，即可得出 \(p\) 点获得的所有能量

\[E(p) = \int_{H^2}L_i(p, \omega)cos \theta d \omega \]

BRDF

BRDF 描述的是光线从一个立体角 \(\omega_i\) 打到某个点后，反射到另外一个立体角 \(\omega_r\) 上的能量之间的比例关系。（也可以理解为光线打到某个点后，被物体吸收，然后物体再发出一些光线）

\[f_r(\omega_i \rightarrow \omega_r) = \frac{反射能量_{观察方向}}{入射能量}= \frac{dL_r(\omega_r)}{dE_i(\omega_i)} = \frac{dL_r(\omega_r)}{L_i(\omega_i)cos\theta_i d\omega_i} \]

\[L_r({p, \omega_r}) = \int_{H^2}f_r(p,\omega_i \rightarrow \omega_r)L_i(p, \omega)cos \theta d \omega \]

兰伯特

兰伯特模型中假设没有光照能量损失，BRDF 公式如下：(Lambertian BRDF 推导过程)

\[f_{lambert} = \frac{\rho}{\pi} \]

\(\rho\)：表示 albedo 颜色

Cook-Torrance 模型

BRDF 有很多种实现模型，实时渲染管线中一般使用的是 Cook-Torrance BRDF 模型，它将反射拆成了漫反射和镜面反射两个部分。

\[f_r = k_d f_{lambert} + k_s f_{cook-torrance} \]

\(k_d\)：这里的表示入射光线中被折射部分的能量所占的比率
\(k_s\)：是被反射的部分的比率

镜面反射部分：

\[f_{cook-torrance}=\frac{DFG}{4(\boldsymbol{\omega_o} \cdot \boldsymbol{n})(\boldsymbol{\omega_i} \cdot \boldsymbol{n})} \]

\(D\)：法线分布函数，估算在受到表面粗糙度的影响下，朝向方向与半程向量一致的微平面的数量。这是用来估算微平面的主要函数。
\(F\)：菲涅尔函数，描述的是在不同的表面角下表面所反射的光线所占的比率。
\(G\)：几何函数，描述了微平面自成阴影的属性。当一个平面相对比较粗糙的时候，平面表面上的微平面有可能挡住其他的微平面从而减少表面所反射的光线。

对于上面的每一种函数，根据实现物理机制，每种函数都有不同的实现方式。

菲涅尔方程 (F)

Fresnel-Schlick近似法

\[\begin{aligned} F_{Schilick}(\boldsymbol{h},\boldsymbol{v},F_0) &= F_0 + (1 - F_0)(1 - (\boldsymbol{h}\cdot \boldsymbol{v}))^5 \\ &= F_0 + (1 - F_0)(1 - cos \theta_h)^5 \end{aligned} \]

\(\boldsymbol{h}\) : 半程向量

法线分布函数 (D)

发现分布函数跟粗糙度 \(\alpha\) 有关系：

Trowbridge-Reitz GGX

\[\begin{aligned} D_{GGX}(\boldsymbol{n},\boldsymbol{h},\alpha) &= \frac{\alpha^2}{\pi ((\boldsymbol{n}\cdot \boldsymbol{h})^2(\alpha^2 - 1) + 1)^2} \\ &= \frac{\alpha^2}{\pi ((cos\theta_h)^2(\alpha^2 - 1) + 1)^2} \end{aligned} \]

Beckmann distribution

\[D(\boldsymbol{h}) = \frac{e^{-\frac{tan^2\theta_h}{\alpha^2}}}{\pi \alpha^2 cos^4\theta_h} \]

IBL 预计算原理

不考虑自发光的情况下，渲染方程式：

\[L_o(\boldsymbol{x},\omega_o) =\int_{H^2}f_r(\boldsymbol{x}, \omega_o,\omega_i)L_i(\boldsymbol{x}, \omega_i)cos\theta_i d\omega_i \]

带入之前的 BRDF 函数：

\[\begin{aligned} L_o(\boldsymbol{x},\omega_o) &= \int_{H^2} (k_d f_{lambert} + k_s f_{cook-torrance}) L_i(\boldsymbol{x}, \omega_i)cos\theta_i d\omega_i \\ &= \int_{H^2} k_d f_{lambert} L_i(\boldsymbol{x}, \omega_i)cos\theta_i d\omega_i + \int_{H^2} k_s f_{cook-torrance} L_i(\boldsymbol{x}, \omega_i)cos\theta_i d\omega_i \\ &= L_d(\boldsymbol{x}, \omega_o) + L_s(\boldsymbol{x}, \omega_o) \end{aligned} \]

辐照度图预计算

对于漫反射部分(把常数项移除积分后得到)：

\[\begin{aligned} L_o(\boldsymbol{x},\omega_o) &= \int_{H^2} k_d \frac{c}{\pi} L_i(\boldsymbol{x}, \omega_i)cos\theta_i d\omega_i \\ &= k_d \frac{c}{\pi} \int_{H^2} L_i(\boldsymbol{x}, \omega_i)cos\theta_i d\omega_i \end{aligned} \]

公式中的积分项只依赖 \(\omega_i\)，这部分就可以通过预先计算，存储到对应的一张新的立体贴图中，图中存储了 \(\omega_o\) 方向上，所有表面间接反射光的总和（也就是积分项中的值）

简而言之，以前需要在环境贴图上采样大量 \(\omega_i\) 积分计算的操作，简化为只需要在新生成的 辐照度图 中采样一次 \(\omega_o\) 方向上的值。

这里不去讨论怎么将三维的向量转换到立体贴图的二维坐标空间采样贴图。

镜面反射部分

现在考虑镜面反射部分：

\[L_s(\boldsymbol{x}, \omega_o) = \int_{H^2} k_s f_{cook-torrance} L_i(\boldsymbol{x}, \omega_i)cos\theta_i d\omega_i \]

注意到这部分因为不仅仅跟 \(\omega_i\) 有关系，还跟 \(\omega_o\)，并且还跟 BRDF 中的 \(\alpha\)、\(F_0\)、\(\theta_h\) 有关系，因此对变量进行降维：

\[f_{cook-torrance}=\frac{DFG}{4(\boldsymbol{\omega_o} \cdot \boldsymbol{n})(\boldsymbol{\omega_i} \cdot \boldsymbol{n})} \]

Epic Games 提出了一种分离近似和的方法，将高光部分的积分拆分成两部分:

\[\begin{aligned} L_s(\boldsymbol{x}, \omega_o) &= \int_{H^2} k_s f_{cook-torrance} L_i(\boldsymbol{x}, \omega_i)cos\theta_i d\omega_i \\ &= \int_{H^2} f_r(p, \omega_i,\omega_o) L_i(\boldsymbol{x}, \omega_i)cos\theta_i d\omega_i \\ &\approx \int_{H^2} L_i(\boldsymbol{x}, \omega_i) d\omega_i \cdot \int_{H^2} f_r(p, \omega_i,\omega_o)cos\theta_i d\omega_i \end{aligned} \]

积分的第一部分称为 预滤波环境贴图，类似于辐照度图，也是预先计算的环境卷积贴图，但是这里需要考虑粗糙度，因此可以预先按照粗糙度级别，预计算不同粗糙度下的的贴图，后续根据真实的粗糙度来插值得到积分结果。

不同粗糙度下的，镜面反射的的波瓣（lobe），越粗糙，波瓣范围越大

这部分没弄太明白

最后只剩下最后一部分积分了，这里单独把菲涅尔方程提出来，然后展开菲涅尔项，主要是把 \(F_0\) 从积分项中拆出来

\[F_{Schilick}(\boldsymbol{h},\boldsymbol{v},F_0) = F_0 + (1 - F_0)(1 - cos \theta_h)^5 \]

这里会把半程角度 \(\theta_h\)，跟入射角度 \(\theta_i\) 都看成 \(\theta_i\)

\[\begin{aligned} \int_{H^2} f_r(p, \omega_i,\omega_o)cos\theta_i d\omega_i &= \int_{H^2} \frac{f_r(p, \omega_i,\omega_o)}{\textcolor{#008000}{\boxed{F(\omega_o, h)}}} \textcolor{#008000}{\boxed{F(\omega_o, h)}} cos\theta_i d\omega_i \\ &= \int_{H^2}\frac{f_r(p, \omega_i,\omega_o)}{\textcolor{#008000}{\boxed{F(\omega_o, h)}}} (F_0 + (1 - F_0)(1 - cos\theta_i)^5)cos\theta_i d\omega_i \\ &= F_0\int_{H^2}\frac{f_r(p, \omega_i,\omega_o)}{\textcolor{#008000}{\boxed{F(\omega_o, h)}}} (1 - (1 - cos\theta_i)^5) cos\theta_i d\omega_i + \int_{H^2}\frac{f_r(p, \omega_i,\omega_o)}{\textcolor{#008000}{\boxed{F(\omega_o, h)}}} (1 - cos\theta_i)^5cos\theta_i d\omega_i \\ &= F_0 \cdot A + B \\ &= F_0 * LUT.r + LUT.g \end{aligned} \]

其中：

\[\frac{f_r(p, \omega_i,\omega_o)}{\textcolor{#008000}{\boxed{F(\omega_o, h)}}}=\frac{DG}{4(\boldsymbol{\omega_o} \cdot \boldsymbol{n})(\boldsymbol{\omega_i} \cdot \boldsymbol{n})} \]