Q132：PBRT-V3，BSSRDF（双向散射表面反射分布函数）（5.6.2章节、11.4章节）

一、BSSRDF的定义

BSSRDF描述的是次表面中光的传播。考虑次表面中光的传播的话，光线离开次表面时的点（Po）和光线进入次表面时的点（Pi）不再是同一位置。如下方图示：

BSSRDF的定义式：

BSSRDF对应的渲染方程的一般形式：

二、BSSRDF的建模与近似

在给定po、wo、pi、wi时，总不能根据BSSRDF的定义式来求解S(po, wo, pi, wi)吧。
最直接的方法是建立各种BSSRDF模型来近似模拟真实的物理情况。

书上提供了这个近似模型：

这个近似去除了位置（po, pi）和方向（wi，wo）之间的耦合。
其中各个分量的具体含义，原文截图如下：

三个分量，一个一个求解呗：

接下来，只剩下Sp分量啦，不过这个也是最为复杂的。

三、Sp(po, pi）

对Sp进行进一步近似和简化：

Sp（po, pi）表示Sp的值由po、pi的位置决定。
Sr（||po-pi||）表示Sr的值只由po、pi之间的距离决定。也就是说，到po相同距离的所有的pi对应的Sr值时相同的，这些所有的pi在一个以po为球心，||po-pi||为半径的球面上，“Sr”中的“r”表示的就是“radial”。称“Sr”为“半径分量”好了。用Sr近似Sp，简化了很多啊，有木有？

“半径分量Sr”的值在各种材料属性都确定的情况下，只和“半径r”有关。但是，不同的材料对应着不同的材料属性。所以，考虑“半径分量Sr”的一般形式：

尼玛，这个形式也忒复杂了吧，必须简化。

到目前为止，“半径分量Sr”已经简化到只和（rho, r）相关。接下来，怎么办呢？

书上讲到：

即：“半径分量Sr”的具体值存放在一个叫做“BSSRDFTable”的表中。
根据（rho, r）到“BSSRDFTable”的表中获取对应的“半径分量Sr”的值。

原来早就算好了所有的“半径分量Sr”的值，来来回回各种简化，最后原来是“查表”哈！

“profile”这个表里存放的是先前计算好的“半径分量Sr”的值，相当于是对Sr进行采样之后的采样点。

现在要用这个表里的数据，可不能只是简单“读取出来直接使用吧”。
应该怎么办？
“采样”对应“重构”嘛。
“profile”表中保存的是“半径分量Sr”的采样点，使用时当然得对采样点做相应的“重构”。

书上采样的重构方式是：双三次张量插值（bi-cubic tensor interpolation)。

什么意思？

先分别对rho_sample和r_sample进行双三次（样条）插值求得各自的权系数，然后求权系数的张量积，然后根据权系数的张量积对对应的“半径分量Sr”的值进行累加。

关于rho_sample和r_sample的“样条插值”，下图假设r=0.5, rho=0.2：

书上计算这个插值的函数是CatmullRomWeights()。

bool CatmullRomWeights(int size, const Float *nodes, Float x, int *offset,                       Float *weights) {    // Return _false_ if _x_ is out of bounds    if (!(x >= nodes[0] && x <= nodes[size - 1])) return false;    // Search for the interval _idx_ containing _x_    int idx = FindInterval(size, [&](int i) { return nodes[i] <= x; });    *offset = idx - 1;    Float x0 = nodes[idx], x1 = nodes[idx + 1];    // Compute the $t$ parameter and powers    Float t = (x - x0) / (x1 - x0), t2 = t * t, t3 = t2 * t;    // Compute initial node weights $w_1$ and $w_2$    weights[1] = 2 * t3 - 3 * t2 + 1;    weights[2] = -2 * t3 + 3 * t2;    // Compute first node weight $w_0$    if (idx > 0) {        Float w0 = (t3 - 2 * t2 + t) * (x1 - x0) / (x1 - nodes[idx - 1]);        weights[0] = -w0;        weights[2] += w0;    } else {        Float w0 = t3 - 2 * t2 + t;        weights[0] = 0;        weights[1] -= w0;        weights[2] += w0;    }    // Compute last node weight $w_3$    if (idx + 2 < size) {        Float w3 = (t3 - t2) * (x1 - x0) / (nodes[idx + 2] - x0);        weights[1] -= w3;        weights[3] = w3;    } else {        Float w3 = t3 - t2;        weights[1] -= w3;        weights[2] += w3;        weights[3] = 0;    }    return true;}

（原理在8.6.1章节，此处略去）

关于这行代码，再补充说明一下：

        // Cancel marginal PDF factor from tabulated BSSRDF profile        if (rOptical != 0) sr /= 2 * Pi * rOptical;

根据式子“11.9”，咱的目的是求解具体点pi对应的Sp，但是根据||po-pi||从Sr表格中查找出来的是“以po为圆心，||po-pi||为半径的po处surface的切平面上一个圆周”对应的Sr。所以，针对单个点pi的Sr，需要在查找结果的基础上除以“周长”。

归根到底是：简化，近似，查表

流水账式的记录，差不多就是这样~