分而治之
在上一章 软光栅下体积云理论与实现(一) 中我们聊到了怎么去定义云这种材质，通过使用体素对于问题进行分化，把非标准的云划分成了一个个标准的体素，这对于接下来的渲染降低了难度。
这回我们聊聊怎么设计光线在云中的传播，从而完成渲染任务。我们一直说“云这种材质”，有点拗口，“云”作为参与到了光线传播的媒介之一，我们后面在通用的一些地方就叫它： participating media.
首先回忆一下最简单一种光线传播方式：光线在bling-phong模型下的传播
可以分为四步：
- 构建一条射线作为光路
- 看场景中是否存在交点
- 获得直接光照的能量
- 根据材质选择一条光路，计算间接光照的能量
伪代码大致是
获得radiance(一条射线){
构建光线
判断场景中的交点
if (存在交点) {
直接光照的能量 = 根据材质、光照位置计算得到
在半球分布上构建下一条光路
间接光照 = 获得radiance（下一条光路）
} else {
return 无能量
}
根据RR、bounce次数等条件 判断停止
}
这样就构成了一幅图像，blingphong的物体会展示出来，并把后面的物体挡住。
对于光线在participating media中的传播，和普通的表面材质非常相似，只有一个地方改变了
“在半球分布上构建下一条光路 ” -> "在整个球分布上构建下一条光路 "。
这里的光路就变成了：
为了方便理解，这里定义两个假设（后面会打开这两个限制）
- 光线每次就只前进一个体素的距离，就一定要bounce一次
- 光在传播过程中就像在真空一样，没有任何折损（先不考虑任何光的散射效果）
就像blingphong模型下的表面反射的能量 由直接光照和间接光照组成，在participating media中的散射也可以按照“直接”+“间接”的思路来理解。
光线每次在体素中bounce一次，都可以看做接受了一次间接光照的能量贡献（上图虚线部分），最后光线从云中出射到背景的物体，也会作为“间接”光照的能量参与到整体能量的共享中。
伪代码可以是
获得radiance(一条射线){
构建光线
if (在 “云” 材质中) {
判断材质中的 体素 的交点
} else {
判断场景中的交点
}
if (存在交点) {
直接光照的能量 = 根据材质、光照位置计算得到
if （在 “云” 材质中） {
在整个球分布上构建下一条光路
间接光照 = 获得radiance（下一条光路）
} else {
在半球分布上构建下一条光路
间接光照 = 获得radiance（下一条光路）
}
} else {
return 无能量
}
根据RR、bounce次数等条件 判断停止
}
根据上面提到的两点限制，我们在伪代码中增加了两个判断，
- 在 “云” 材质中，由于我们要求每次只能前进一个体素，所以对于场景物体的求交改为了体素中的求交
- 在 “云”材质中，散射不局限于半球范围，而是整个球
思路上这样就没问题了，不过读起来着实有点吃力，if套if不是好的工程实践，我们稍微改一下伪代码：
获得radiance(一条射线){
构建光线
for (bounce次数 < 最大bounce次数) {
//tips：间接光照 （例如第2次bounce）时 光线方向已经改变
if (在 “云” 材质中) {
判断材质中的 体素 交点
} else {
判断场景中的 表面 交点
}
if (存在体素交点) {
计算光照能量
在整个球分布上构建下一条光路
更新光线方向
} else if （存在表面交点） {
计算光照能量
在半球分布上构建下一条光路
更新光线方向
} else {
return 无能量
}
if （命中RR） {
break;
}
}
return 光照能量
}
这样问题就解决了…吗？
表面与次表面
并没有。
这里存在一个问题。
如果按照上面这样设计光路，那不就和一块玻璃一样了吗？
就像做Whitted style ray tracing 也会根据fresnel公式来计算整个球的反射方向，难道云就是多bounce了几次的玻璃？
云材质效果
玻璃材质效果
从感性的认知来看肯定是不对的，云和玻璃显然表现是非常不一样的。差在哪里呢？
上面的光线传播提到了两个假设
- 光线每次就只前进一个体素的距离，就一定要bounce一次
- 光在传播过程中就像在真空一样（先不考虑任何光的散射效果）
这里要打破第二个假设。
在玻璃中，玻璃内部的光线不会受到任何影响，其方向、光的能量都不会有改变。所以如果也用体素的思路去构建一块立方体玻璃也是没有问题的。
光线在玻璃内部同等于在真空传播
但光并不能在云中像在真空一样无所影响的传播，上一节聊到了云的材质其中有无数的水滴或冰晶，即使我们已经把云划分成了一个个的相同小体素，体素内部也会对光线造成影响。
这个影响分为两个部分，一个是光的能量，一个是光的方向。
散射效果对能量的影响
根据上一章 软光栅下体积云理论与实现(一) 已经了解到了光在体素中会有4种散射效果，接下来就来看他们是如何影响到光的能量，是会增加还是减少，具体又改变了多少。
只要获得了光在每个体素中是如何变化，上一节所说的光的传播线路就可以被计算出来。
我们可以朴素的建立一套伪代码：
改变的光能量 = Absorption的能量 + Emission的能量 + in-scattering能量+ outscattering能量
- Absorption 效果
每一个体素中有若干个分子，这些分子有的会参与到散射效果，有的不会，这种不确定的问题，照旧用我们的老办法：概率统计。
对于每个体素可以继续拆分为更小的体素，可以假设：
在大体素中每走过的每个小体素，都有确定影响，例如：光的能量减少1%，那么走过的小体素越多，能量少的就越多，最终可以得到大体素的absorption效果
可以得到一套伪代码：
Absorption的能量 = 入射的光能量 * 每个小体素的减少比例 * 光路的距离
- Emission 效果
同样使用概率统计的思路，伪代码为
emission的能量 = 入射的光能量 * 每个小体素的释放比例 * 光路的距离
接下来的工作就都是相同的了，in-scattering和out-scattering的效果也是如此，这样可以得到新的伪代码，我们把所有小体素的效果比例称之为“系数”就可以得到：
改变的光能量 = 四种效果的系数和 * 光路距离 * 入射的光能量
其中，in-scattering和out-scattering会根据方向所调整，从而合并得到一个scattering系数。
而absorption作为负向影响，可以与正向的scattering再次合并，从而得到一个attenuation efficient记作-σ_t；
而类似渲染方程会把自发光项提出来一样，对于participating media中的散射效果的影响也会如此处理。
最终可得：
改变的光能量 = 自发光 + attenuation_efficient * 光路距离 * 入射的光能量
tips：
这里为了方便理解，有许多隐式假设，将会在后面打开
从数学上来说，上面的式子是没问题的，但是从工程上来说，还可以再优化一下，显而易见有几个要求：
- 改变的能量需要小于总能量且大于零
- 光路的距离，也就是participating media的“厚度”thickness ,与光的能量变化成正相关
这两个要求翻译一下，对于工程的限制就是：
限制一：改变的能量 / 入射能量 的值在 [0,1] 的范围内
限制二：thickness即使无限增加，也应该满足限制一
结合上面的限制，可以设计一个用于描述光能量衰减的比例 beam transmittance 记为TR
以e为底的指数函数就可以满足上面的要求，e^{-x} 的图像如下：
其中的连续性与可导性非常适合我们的工程要求，并且指数衰减也符合电磁波的物理规律。
结合工程限制以及考虑散射效果，现在就可以计算出光在体素内一次bounce会改变多少radiance了：
L_o = L_i * Tr
单个体素的计算完成，就可以计算整个participating media中的能量情况了，对于整条光路，亦是单个体素的叠加：
所以结合各个体素情况可得：
Tr_0n = Tr_01 * Tr_12 * Tr_2n // Tr_01 代表 从p0到p1 的TR，下同
L_o = Li * Tr_0n
最终就可以得到出射光线的radiance。
散射效果对方向的影响
文章开头聊到了“ ‘在半球分布上构建下一条光路’改变为 ‘在整个球分布上构建下一条光路’ ”。
散射对于光的影响不仅是上一节提到的能量衰减，同时还是散射方向的改变。
可以进行一下对比：
反射仅在法线方向上进行的分布；散射会在整个球的范围进行分布。
反射有BRDF描述分布情况；散射也必然有对应的模型来进行描述。
这就是phase function。
Phase function
Phase function在职责上与BRDF相比，不能说毫无区别，也可以称为一模一样。都负责两件事：
- 确定出射方向
- 确定强度分布
不过相对于BRDF一般由严谨的物理性质推导而来，phase function更多是经验模型和概率统计而来；
类比与BRDF有很多模型：phong、microfact，同样的phase function也有很多模型，这里介绍广为使用的“The Henyey-Greenstein Phase Function”。
对于出射方向与强度分布，最简单的一种情况是在整个球上随机出射、整个球上均匀分布，这种情况下，伪代码为：
w_i = randomSoildAngle();
f = 1/4 * INV_PI;
也有一种可能是phasefunction也像BRDF一样，不同的材质有不同的lobe，类似于lobe会有自己的倾向，不一定是均匀的
diffuse module 的lobe 分布
Glossy specular 的lobe 分布
这种情况下能量并不会均匀的分布在整个球上，而是有明显的方向倾向性，可能会沿着入射方向散射，也有可能向着相反方向散射，于是我们就可以使用“The Henyey-Greenstein Phase Function” 来描述这种倾向性：
这里引入了参数 g （asymmetry parameter） ，是一个根据材质预先设计好的常数，表示光的前向散射倾向性：
图1 “The Henyey-Greenstein Phase Function” 公式
g=-0.7 光线向入射方向的反方向散射
g=0.7 光线向入射方向的同方向散射
g越大的情况下，散射方向就会越靠近入射方向。根据入射与出射的cos夹角，可以得出对应方向的强度分布。
对于出射的方向，我们当然可以继续均匀的进行全球采样，phase function自然会给出对应强度情况。
只是这样会不会有什么问题？
这个后面再讨论，我们先把最小流程跑通再说。
综上可以得出伪代码：
w_i = randomSoildAngle();
f = PHG(Dot(w_i,W_o),g);// PHG即为图1公式
次表面散射光路
这下终于完成了所有前置条件，我们可以把散射对于能量和方向的影响都加入到path tracing的逻辑中了，伪代码为：
获得radiance(一条射线){
// 构建光线
ray = Ray();
L = Vector(0.f);// radiance
f = 1.f; // 散射率
for (bounce次数 < 最大bounce次数) {
//tips：间接光照 （例如第2次bounce）时 光线方向已经改变
if (在 participating media 中) {
判断材质中的 体素 交点
} else {
判断场景中的 表面 交点
}
if (存在体素交点) {
// 计算当前点的“直接光照”radiance
L = 光照到达体素交点的radiance();
// 根据 beam transmittance 计算达到能量的损失情况
L = L * f * tr;
// 在整个球分布上构建下一条光路
w_i = randomSoildAngle();
// 根据participating media 预设的g 计算能量比例
f = PHG();
// 根据出射方向更新光线方向
ray = swapRay(w_i);
} else if （存在表面交点） {
计算光照能量
在半球分布上构建下一条光路
更新光线方向
} else {
return 无能量
}
if （命中RR） {
break;
}
}
return 光照能量
}
tips: 上面的模型都存在一个前提，participating media是各向同性的，一旦是各向异性的，不仅Tr的计算需要考虑入射方向，PHG的计算也要考虑，将在后续下文中迭代处理。
这样一来在participating media下的次表面散射效果的基本框架就建立了出来，逻辑上也很通顺。
但
真的会如此顺利吗？
均匀的散射分布到底会不会有问题？
数学上没问题，工程上有大问题到底是什么挑战？
单个的散射组件怎么集成到复杂的piccolo系统中？
如果你觉得看伪代码不方便，也可以直接看楼主写的项目代码：Github
后续越来越复杂啦，争取一周内更新。
连载不易，原创费时，各位认可的求求给个免费的赞再走吧~
如有错误之处，恳请大家批评指正，学习交流。
创作等级
会员等级