光线求交加速算法：边界体积层次结构(Bounding Volume Hierarchies)2-表面积启发式法(The Surface Area Heuristic)

光线求交加速算法：边界体积层次结构(Bounding Volume Hierarchies)2

上篇的两种图元分区方法(Middle,EqualCounts)对于某些图元分布可以很好地工作，但是在实践中它们常常选择性能较差的区，从而导致光线访问树的更多节点，因此在渲染时不必要地降低了光线-图元相交的计算效率。当前用于构建射线追踪加速结构的最佳算法中，大多数都是基于“表面积启发式”(SAH)的。

 

表面积启发式法(The Surface Area Heuristic)

该算法提供了扎实的成本模型来估计了哪个划分位置会导致射线和图元相交最为便宜。SAH模型估计执行射线相交测试的计算成本，包括遍历树的节点所花费的时间，以及针对特定图元分区的射线与图元相交测试所花费的时间。然后，用于构建加速结构的算法以将总成本降至最低为目标。通常，使用一种贪婪算法，该算法可最大程度地减少单独构建的层次结构中每个节点的成本。

SAH成本模型背后的思想很简单，它有两种方法计算。第一种：在构建自适应加速结构(图元细分或空间细分)的任何时候，我们都可以为当前区域和几何图形创建叶节点，即直接将当前结点下所有图元归为一个叶子结点内计算求交花费。在这种情况下，穿过该区域的任何射线都将针对所有重叠的图元进行测试，并且将产生以下花费：

其中N是图元的数量。是计算光线和第i个图元相交的时间。

第二种是区域划分：

其中是遍历内部节点并确定射线穿过哪个子节点所花费的时间，和是射线穿过每个子节点的概率，和是两个子节点中图元的索引，和分别是与两个子节点的区域重叠的图元数。图元如何划分的选择会影响两个概率的值以及拆分每一侧的图元集。

我们将简化假设，即对于所有基元都是相同的。这个假设可能与现实相去不远，它引入的任何错误似乎都不会对加速器的性能产生太大影响。当然我们也可以向Primitive添加一个方法，该方法返回交集测试所需的CPU周期数的估计值(但没什么必要)。

概率和可以使用根据几何概率得出的想法来计算。可以看出，对于包含在另一个凸体积B中的一个凸体积A，穿过B的均匀分布的随机射线也将穿过A的条件概率是它们的表面积之比和：

当我们的splitMethod的值为SplitMethod :: SAH时，我们选择SAH方法作为分区的算法来构建BVH。我们通过计算多个可能分区，选择它们中最小SAH估计成本的拆分方案。值得注意的是当递归执行到图元个数非常少时，在使用此方法会造成比较大的浪费，我们选择直接利用等量划分的方法(EqualCounts)。我们列出SAH划分代码：

case SplitMethod::SAH:
default: {
    if (nPrimitives <= 4) {
    //+当图元个数小于4时，用等量划分法。(+)表示代码展开
    } else {
        //+为SAH分区存储桶分配桶信息BucketInfo，并初始化
        //+计算每个桶的分区花费(划分点在桶位置后)
        //+寻找最小的SAH花费
        //+在选定的SAH存储桶中创建叶子结点或拆分图元集
    }
}

这里的实现方式不是沿轴详尽地考虑所有可能的分区。而是沿轴将范围划分为少量相等范围的存储区，且只考虑位于存储桶边界的分区。我们沿轴计算每个桶边界的SAH以选择最佳分区。 这种方法比考虑所有分区更有效，而通常仍会产生几乎一样有效的分区。 下图说明了这个想法：

 如图，将图元边界框的质心投影到所选的分割轴上。 每个图元根据其边界的质心沿轴放置在存储桶中。 然后，该实现会估算沿每个桶边界(蓝色实线)拆分图元，并估算成本。最后选择以表面积启发法给出的最小成本。

我们写出初始化桶的代码：

// +为SAH分区存储桶分配桶信息BucketInfo，并初始化
constexpr int nBuckets = 12;
struct BucketInfo {
    int count = 0;
    Bounds3f bounds;
};
BucketInfo buckets[nBuckets];

for (int i = start; i < end; ++i) {
    int b = nBuckets * centroidBounds.Offset(primitiveInfo[i].centroid)[dim];    //b为哈希桶索引
    if (b == nBuckets) b = nBuckets - 1;
    buckets[b].count++;
    buckets[b].bounds = Union(buckets[b].bounds, primitiveInfo[i].bounds);
}

其中centroidBounds.Offset(const Point3<T> &p)函数返回点p相对于边界框centroidBounds的位置(范围为0到1)。因此哈希桶索引被约束在[0，12)。每次有图元分配到桶内，桶内元素个数加1，并合并它们的边界框。

接下来我们循环遍历所有存储桶，并初始化cost [i]数组以存储估计的SAH成本，以便在第i个存储桶之后进行拆分(最后一个桶之后不考虑)。我们代码如下：

// +计算每个桶的分区花费(划分点在桶位置后)
Float cost[nBuckets - 1];
for (int i = 0; i < nBuckets - 1; ++i) {
    Bounds3f b0, b1;
    int count0 = 0, count1 = 0;
    for (int j = 0; j <= i; ++j) {
        b0 = Union(b0, buckets[j].bounds);
        count0 += buckets[j].count;
    }
    for (int j = i+1; j < nBuckets; ++j) {
        b1 = Union(b1, buckets[j].bounds);
        count1 += buckets[j].count;
    }
    cost[i] = .125f + (count0 * b0.SurfaceArea() +
                       count1 * b1.SurfaceArea()) / bounds.SurfaceArea();
}

我们将估计的交叉点成本()设置为1，然后将估计的遍历成本()设置为1/8。它们其中任意一个都可以设置为1，因为确定其影响的是估计的遍历和相交成本的相对值。经验算，相交成本约为遍历成本8倍。

要获取最低成本值，只需要遍历一次数组就能获得：

// +寻找最小的SAH花费
Float minCost = cost[0];
int minCostSplitBucket = 0;
for (int i = 1; i < nBuckets - 1; ++i) {
    if (cost[i] < minCost) {
        minCost = cost[i];
        minCostSplitBucket = i;
    }
}

如果选择的用于划分分区的存储区边界的估计成本(第二种)低于直接使用叶结点估计的成本(第一种)，或者如果存在的节点数量超出了节点中允许的最大基元数量，则可以使用std :: partition()函数来完成工作对原始信息数组中的节点进行重新排序:

// +在选定的SAH存储桶中创建叶子结点或拆分图元集
Float leafCost = nPrimitives;//第一种方法计算的花费，求交花费为1所以累加是nPrimitives
if (nPrimitives > maxPrimsInNode || minCost < leafCost) {
    BVHPrimitiveInfo *pmid = std::partition(
        &primitiveInfo[start], &primitiveInfo[end - 1] + 1,
        [=](const BVHPrimitiveInfo &pi) {//[=]捕获了nBuckets，centroidBounds，dim
            int b = nBuckets * centroidBounds.Offset(pi.centroid)[dim];
            if (b == nBuckets) b = nBuckets - 1;
            return b <= minCostSplitBucket;//找到第一个b大于minCostSplitBucket的指针
        });
    mid = pmid - &primitiveInfo[0];
} else {
    // +创建叶子结点
}

这里创建叶子结点的代码和之前的依旧是一样的。