MIT 6.837:Radiosity 辐射度全局光照

MIT 6.837:Radiosity 辐射度全局光照

在这之前,我们介绍的局部光照方法是对实际情况的一个近似。而在本文中,我们将基于光照的物理原理得到一种与现实的光照完全一致的方法,这种方法中采用了辐射度的概念与方法。

康奈尔盒(The Cornell Box)

69 - MIT 6.837:Radiosity 辐射度全局光照
70 - MIT 6.837:Radiosity 辐射度全局光照

上面两幅图片先后是光线追踪的结果与全局光照的结果。比较发现,全局光照的效果较亮,且能正确处理光线被多次反射后的效果,如玻璃球下的亮斑。

上面的结果说明,前面介绍的光线追踪算法并不能真正地模拟出真实世界中的光照。为了研究如何模拟真实光照,康奈尔大学光线测量实验室制作了这个康奈尔盒,它由现实世界中的一个盒子与建模得到的盒子场景构成,可以通过比较照片与渲染结果来验证渲染的正确性。

71 1024x457 - MIT 6.837:Radiosity 辐射度全局光照

正确实现的全局光照算法可以实现和现实光照完全一致的效果。

实现全局光照的方法可以分为两种,基于辐射度的方法和基于蒙特卡洛的方法。其中,辐射度方法只能计算散射材质,但结果与相机无关,而蒙特卡洛方法需要追踪很多间接光线,计算开销大。

渲染方程(The Rendering Equation)

光照方程是所有全局光照方法的基础,它以物理中辐射度的原理提出,表示如下

71 - MIT 6.837:Radiosity 辐射度全局光照

$$ L(x’, \omega’) = E(x’, \omega’) + \int \rho_{x’}(\omega, \omega’) L(x, \omega) G(x, x’) V(x, x’) \mathrm{d}A $$

这个方程的每一项含义如下:

  • $L(x’, \omega’)$:从点 $x’$ 处向立体角 $\omega’$ 射出的辐射率(radiance,$\mathrm{W \cdot m^{-2} \cdot sr^{-1}}$),表示点 $x’$ 处射出的光线
  • $E(x’, \omega’)$:点 $x’$ 处材质向立体角 $\omega’$ 方向的自发光辐射率
  • $\int \mathrm{d}A$:对其他物体的所有平面面积进行积分,从而统计从各个方向入射的辐射率
  • $L(x, \omega)$:从另一物体上的点 $x$ 向点 $x’$ 方向的立体角 $\omega$ 射出的辐射率,表示从另一物体射入的光线
  • $\rho_{x’}(\omega, \omega’)$:点 $x’$ 处材质将射入辐射转化为射出辐射的比例,又称为双向分布函数(BRDF)
  • $V(x, x’)$:用于表示点 $x$ 和 $x’$ 之间是否有遮挡,如遮挡为 0,表示辐射不能传播到点 $x’$,否则为 1
  • $G(x, x’)$:反映点 $x$ 和 $x’$ 所在面之间的关系

基于辐射度的全局光照方法

辐射度方法假设所有面都是理想散射材质(反射光线向所有方向均匀发出),而场景被划分为若干小的区域,认为区域内每一点辐射出的辐射率是相等的。

辐射度公式与辐射度矩阵

基于辐射度方法的描述,我们可以对渲染方程进行简化,将辐射率函数 L 换为每一区域中的辐射率 B,并将 BRDF 简化为常数,即

$$ B_{x’} = E_{x’} + \rho_{x’} \int B_x G(x, x’) V(x, x’) \mathrm{d}A $$

而我们将场景划分后,每一个区域之间的 G、V 关系是确定的,且辐射率为常数,因此可以得到离散化后的方程

$$ B_i = E_i + \rho_i \sum_{j=1}^n F_{ij} B_j $$

其中,$F_{ij}=\int G(x, x’) V(x, x’) \mathrm{d} A$ 因两区域的确定关系,可预先计算得到,称作视界因子(form factor)。

上面的方程是关于 $B_i$ 的线性方程组,因此可表达为矩阵形式

$$
\begin{bmatrix}
1-\rho_1F_{11} & -\rho_1F_{12} & \cdots & -\rho_1F_{1n} \\
-\rho_2F_{21} & 1-\rho_2F_{22} & \cdots & -\rho_2F_{2n} \\
\vdots & \vdots & \ddots & \vdots \\
-\rho_nF_{n1} & -\rho_nF_{n2} & \cdots & 1-\rho_nF_{nn}
\end{bmatrix}
\begin{bmatrix} B_1 \\ B_2 \\ \vdots \\ B_n \end{bmatrix}
=
\begin{bmatrix} E_1 \\ E_2 \\ \vdots \\ E_n \end{bmatrix}
$$

这个方程的解可以通过 Gauss-Seidel 迭代方法获得,这是一种求出线性方程组的近似解的方法,利用上一次的解迭代得到本次的解,若干轮迭代后解将收敛。

$$
\begin{bmatrix} B_1 \\ B_2 \\ \vdots \\ \textcolor{red}{B_i} \\ \vdots \\ B_n \end{bmatrix}
=
\begin{bmatrix} E_1 \\ E_2 \\ \vdots \\ \textcolor{red}{E_i} \\ \vdots \\ E_n \end{bmatrix}
+
\begin{bmatrix} \\ \\ \\ \textcolor{red}{\rho_iF_{i1}} & \textcolor{red}{\rho_iF_{i2}} & \textcolor{red}{\cdots} & \textcolor{red}{\rho_iF_{in}} \\ \\ \end{bmatrix}
\begin{bmatrix} \textcolor{red}{B_1} \\ \textcolor{red}{B_2} \\ \textcolor{red}{\vdots} \\ \textcolor{red}{B_i} \\ \textcolor{red}{\vdots} \\ \textcolor{red}{B_n} \end{bmatrix}
$$

但渲染时,我们实际上需要的是顶点的颜色值,而区域的辐射度表示的是区域内部的颜色值,我们可以通过对顶点周围的区域颜色取平均值获得顶点颜色,并在区域内插值得到所有像素的颜色值。

求解视界因子

在辐射度矩阵中出现的形状因数 $F_{ij}$,在辐射度意义下,表示从区域 j 发出的辐射中,被区域 i 吸收的部分,这一关系可以根据图示中的两来刻画

72 - MIT 6.837:Radiosity 辐射度全局光照

$$ F_{ij} = \frac{1}{A_i} \int_{A_i} \int_{A_j} \frac{\cos \theta_i \cos \theta_j}{\pi r^2} V_{ij} \mathrm{d}A_j \mathrm{d}A_i$$

上式中的关系中,$\cos \theta$ 表示因斜角度导致的功率损失,而 $1/r^2$ 代表距离上的损失,$1/\pi$ 为理想散射材质的 BRDF。

理想散射材质的 BRDF 推导如下:

对于理想散射材质,假设其入射光辐照度 $L_i$(radiance,$\mathrm{W \cdot m^{-2}}$),BRDF 为常数 $\rho$(单位 $\mathrm{sr^{-1}}$),出射光辐射度 $L_o$。出射光辐射度可以通过对立体角积分得到为

$$ \begin{aligned}
L_o &= \int_\Omega \rho L_i \cos \theta \mathrm{d}\omega \\
&= \int_0^{2\pi} \int_0^{\pi/2} \rho L_i \cos\theta \sin\theta \mathrm{d}\theta \mathrm{d}\varphi \\
&= 2\pi \rho L_i \int_0^{\pi/2} \cos\theta \sin\theta \mathrm{d}\theta \\
&= 2\pi \rho L_i \cdot \frac12 = \pi \rho L_i
\end{aligned} $$

由于能量守恒,材质无自发光与吸收能量,入射辐射应等于出射辐射,故有 $\pi \rho=1 \Rightarrow \rho=\dfrac{1}{\pi}$。

73 - MIT 6.837:Radiosity 辐射度全局光照

一种计算视界因子的方式是 Nusselt analog,考察面 $A_j$ 对于 $A_i$ 的视界因子,只需要将 $A_j$ 投影到以 $A_i$ 一点为球形的单位球面上,并将球面区域再投影至平面上,求投影占单位圆面积的比例。投影到球面上等效为求立体角(即考虑距离因素),第二次投影等效为考虑斜角度,求比例等效为应用 BRDF。

74 - MIT 6.837:Radiosity 辐射度全局光照

在实际计算中,常用半立方体方法(Hemicube Algorithm)来处理视界因子。对于每一个区域,在其中心周围构造半立方体,半立方体的面被划分为类似于像素的小区域,其他区域在立方体上用类似光栅化与 Z 缓冲区的方式计算每个像素的视界因子。

75 - MIT 6.837:Radiosity 辐射度全局光照

另一种计算的方法是基于光线追踪的方法,在两区域之间追踪多条光线(4~32 条),用于得到其可见性,并近似得到视界因子的积分数值。由于计算量大,可以在计算时同时保存 $F_{ij}$ 和 $F_{ji}$ 以减少冗余计算。

渐进式完善

渲染完整个场景的计算量较大,呈现给用户之前需要有较长的准备时间,如果有一种方法可以在较短时间内得到一个完成度较高的结果,并在后续实时地完善,就能够提升用户体验。渐进式完善的的辐射度全局光照基于一个思路,每一次迭代并非只更新一个区域的辐射度,而是更新所有区域(即整个场景)的辐射度。

76 - MIT 6.837:Radiosity 辐射度全局光照

考虑将辐射的发射与接收方翻转,每一次迭代计算一个区域对其他区域的辐射度,这种方法也被称为 Southwell 迭代,矩阵表示如下。

$$
\begin{bmatrix} \textcolor{red}{B_1} \\ \textcolor{red}{B_2} \\ \textcolor{red}{\vdots} \\ \\ \textcolor{red}{\vdots} \\ \textcolor{red}{B_n} \end{bmatrix}
=
\begin{bmatrix} \textcolor{red}{B_1} \\ \textcolor{red}{B_2} \\ \textcolor{red}{\vdots} \\ \\ \textcolor{red}{\vdots} \\ \textcolor{red}{B_n} \end{bmatrix}
+
\begin{bmatrix}
\cdots & \textcolor{red}{\rho_1F_{1i}} & & \cdots \\
\cdots & \textcolor{red}{\rho_2F_{2i}} & & \cdots \\
\\ \\ \\
\cdots & \textcolor{red}{\rho_nF_{ni}} & & \cdots
\end{bmatrix}
\begin{bmatrix} \\ \\ \vdots \\ \textcolor{red}{B_i} \\ \vdots \\ \end{bmatrix}
$$

辐射度全局光照的高级话题

自适应划分

通常阴影边缘的辐射度变化较快,而阴影、物体内部的变化较缓慢。因此,应该在细节较多的地方划分为更小的区域以提升渲染画面的品质。

78 1024x674 - MIT 6.837:Radiosity 辐射度全局光照

不连续网格

78 - MIT 6.837:Radiosity 辐射度全局光照

如上图所示,非点光源会在遮挡物下形成本影和半影,这些阴影互相重叠,形成了多个不连续区域。不连续网格将这些区域的边界计算出,并分区域地计算其辐射度。这种方法能够解决阴影效果的锯齿,并能得到品质高的阴影边界,但需要引入新的几何,从而增加计算的复杂度。下图中蓝色线条就是多个非点光源光照被椅子及桌子遮挡形成的阴影中的不连续处。

79 - MIT 6.837:Radiosity 辐射度全局光照

层次化方法

可以将光照品质要求不高的区域合并,以减少总的区域数目。但合并的策略需要选择,并且实际实现时并不简单,可能需要耗费内存。

80 - MIT 6.837:Radiosity 辐射度全局光照

总结

辐射度方法是一种实用方法,它即可用于离线领域,也可用于实时领域,例如物理仿真、游戏光照图的预处理等。但它的应用收到许多限制,例如要求材质必须是理想散射材质,并且区域网格的划分、视界因子的计算等也存在局限。

另一种全局光照方法被称作蒙特卡洛方法,这种方法能处理的情况更加多,通用性更好,且是研究热点之一。这一方法将在后面的内容中研究。

参考资料



发表评论

您的电子邮箱地址不会被公开。

此站点使用Akismet来减少垃圾评论。了解我们如何处理您的评论数据