在上一篇文章《从拉普拉斯方程到球谐函数》中，我们详细的讲解了由拉普拉斯方程求解得到球谐函数的过程以及介绍了球谐函数的一些性质：

xiaomengge：【球谐光照】(一)从拉普拉斯方程到球谐函数?

zhuanlan.zhihu.com图标

在本文中将会介绍如何用球谐函数去近似游戏中的光源，也就是得到我们所说的球谐光照。本文主要分以下四个方面去介绍这一过程：

1.蒙特卡洛积分

2.正交基函数

3.光源投影

4.总结

一、蒙特卡洛积分

随机变数X会根据一个特定的函数来进行分布，我们记作 f(x) ~ p(x). 这里f(x)是关于随机变数X的函数，p(x)为随机变数X取特定值x对应的概率，也就是我们说的概率密度函数，概率密度函数在X的取值范围的积分为1。

int_{-infty}^{+infty}p(x)dx = 1，这里p(x) geq 0

所以随机变数X在[a,b]范围内发生的概率为概率密度函数在a到b上的积分：

P(xin[a,b])=int_{a}^{b}p(x)dx

每个关于随机变数X的函数f(x)都有一个期望值，它是你对f(x)做足够多的采样时将会趋近的一个值，对于关于连续随机变数的函数f(x)期望求解如下：

E[f(x)])=int_.f(x)p(x)dx

另一种求期望的方式是对f(x)进行大量的采样，然后求这些采样的均值，当采样的数量趋近无限时所得到的均值也就越趋近实际的期望值：

E[f(x)])approxfrac{1}{N}sum_{i=1}^{N}{f(x_{i})}

结合以上两种求期望的方法，我们就得到了一种在工程数学中对函数积分进行近似计算的巧妙方法——蒙特卡洛积分：

int_.f(x)dx=int_.frac{f(x)}{p(x)}p(x)dxapproxfrac{1}{N}sum_{i=1}^{N}frac{f(x_{i})}{p(x_{i})}

对于这个公式中，我们 frac{f(x_{i})}{p(x_{i})} 看为一个整体然后结合上面两种求期望的公式就很好理解了。我们令 w(x_{i})=frac{1}{p(x_{i})} 来表示每一个采样函数值对应的权重，这样可以得到蒙特卡洛积分的另一种写法：

int_.f(x)dxapproxfrac{1}{N}sum_{i=1}^{N}{f(x_{i})}{w(x_{i})}

二、正交基函数

基函数可以看做是一些很小的信号，通过对这些信号进行缩放然后组合就可以得到原函数的近似函数。计算每个基函数的缩放系数的过程称为投影，例如计算基函数B(x)的系数也就是计算原函数f(x)中包含了多少基函数B(x)，计算基函数系数只需对原函数和对应基函数的乘积进行积分即可：

然后我们用投影系数乘以对应的基函数：

然后把结果相加，得到原函数的近似函数：

上面是用多个基函数对原函数进行近似的过程，通常把多个具有正交性的基函数放到一个集合，我们称这个集合为正交多项式，正交多项式中任意两个多项式都满足：

int_{-1}^{1}F_{m}(x)F_{n}(x)dx=0,n e m\ int_{-1}^{1}F_{m}(x)F_{n}(x)dx=c,n= m

这里如果c=1，则称该多项式为标准正交多项式。

三、光源投影

在本节中，我们将会通过一个具体的例子来使用蒙特卡洛积分将一个光源投影为球谐函数。首先给出我们的光源函数-是由两个互成90度的单色光源组成，并且直接在球坐标系中定义：

light( heta,varphi) = max(0,5cos( heta)-4)+max(0,-4sin( heta-pi)*cos(varphi-2.5)-3)

下图中右边是光源在球坐标系下的函数图，左边是光源映射在单位球面上的效果。

首先我们知道在球面上对任意一个函数f(x)进行积分的方法为：

int_{0}^{2pi}int_{0}^{pi}f( heta,varphi)sin heta d heta dvarphi

通过上面利用积分求基函数的系数，球面坐标求光源light( heta,varphi) 对应自由球谐基函数(如果对球谐函数不是非常了解，请参照上一篇文章球谐光照【一】)的 y_{i}( heta,psi) 的系数应为：

c_{i}=int_{0}^{2pi}int_{0}^{pi}light( heta,varphi)y_{i}( heta,varphi)sin heta d heta dvarphi

函数中的每一组( heta,varphi) 表示的是球坐标系中单位球上的一个点 x_{j} ，带入上面式子同时对它进行蒙特卡洛积分可以得到以下式子：

c_{i}=frac{1}{N}sum_{j=1}^{N}{light(x_{j})y_{i}(x_{j})}w(x_{j})

这里令 f(x_{j})=light(x_{j})y_{i}(x_{j}) ，我们可以得到上面蒙特卡洛积分公式，因为在球面上进行的是无偏均匀的采样，单位球面积为4π，因此每个点采样概率为1/4π，所以这里的权重为w(x_{j}) = 4pi ，所以得到：

c_{i}=frac{4pi}{N}sum_{j=1}^{N}{light(x_{j})y_{i}(x_{j})}

计算球谐系数的代码如下：

typedef double (*SH_polar_fn)(double theta, double phi);

void SH_project_polar_function(SH_polar_fn fn, const SHSample samples[], double result[])
{
const double weight = 4.0*PI;
// for each sample
 for(int i=0; i<n_samples; ++i)
{
double theta = samples[i].sph.x;
double phi = samples[i].sph.y;
for(int n=0; n<n_coeff; ++n)
{
result[n] += fn(theta,phi) * samples[i].coeff[n];
}
}
// divide the result by weight and number of samples
 double factor = weight / n_samples;
for(i=0; i<n_coeff; ++i)
{
result[i] = result[i] * factor;
}
}

这里我们使用4阶的球谐函数，并且在球面上采样10000个点，这样我们得到每个球谐基函数 y_{i} 对应的 c_{i} 如下：

这样我们根据基函数计算光源的近似函数:

ilde{light(s)} = c_{1}y_{1}(s) + c_{2}y_{2}(s) + c_{3}y_{3}(s) +...=sum_{i=1}^{n^{2}} c_{i}y_{i}(s)

式子中 ilde{light(s)} 是光源light(s) 的低频近似光源，s表示球面上任意一点，我们根据 ilde{light(s)} 照亮球面：

用仅有的16个系数实现的效果不算差。注意在低频光源中有一些多余的小柱体出现(右图所示)，而在球面上的显示就是一些不可以预计的暗的部分被照亮(左图所示)。这是因为在原光源 light(s) 函数在0附近有高频信号分量，而我们用4阶的低频球谐函数重构光源导致高频部分丢失。如果用更高阶的球谐函数近似计算最终可以消除这种问题。

四、总结

在本文中我们用4阶的球谐函数对光源进行了近似计算，并且简单给出了进行近似光源的效果，这里我们涉及到了球谐函数、蒙特卡洛积分以及基函数对函数进行重建等相关内容。接下来我们会介绍使用球谐光源如何计算物体表面光照。

参考文献：

1.silviojemma.com/public/

2.基于物理的渲染-基于球面调和基的实时全局光照明 - UWA Blog

3.en.wikipedia.org/wiki/M

4.youtube.com/watch?

5.pdfs.semanticscholar.org


推荐阅读：
相关文章