这个系列是Physically Based Rendering: From Theory To Implementation的读书笔记，本节主要介绍光线求交的加速数据结构。
在光线追踪(ray tracing)的框架中最基本的运算是光线与场景中物体(三角网格)求交。如果直接使用循环来计算的话每条光线的计算复杂度是$O(n)$，这会导致整个渲染过程非常的低效。注意到在大多数情况下对任意的一条光线在计算相交时它最多只会和场景中的极小部分三角形相交，因此我们可以设计出更高效的数据结构来管理渲染的场景，使得在每次计算时可以忽略场景中大部分无关的三角形。这样的加速数据结构一般可以分为基于空间划分(spatial subdivision)的数据结构以及基于包围盒(bounding boxes)的数据结构。前者的代表是kd树(kd-tree)，而后者的代表是BVH(bounding volume hierarchies)。
Bounding Volume Hierarchies
在BVH中我们使用bounding box来表示场景中的每一个三角形，每个bounding box都是与世界坐标系对齐的长方体(axis-aligned bounding box, AABB)，因此可以高效地计算光线与AABB是否相交。显然如果光线不和AABB相交那更不可能和内部的三角形相交，因此在计算和三角形求交前需要判断光线和AABB的相交情况，如果光线和AABB相交则进一步计算是否和三角形相交。除此之外BVH还使用一棵二叉树来把场景中的三角形组织起来，每个三角形都存储在二叉树的叶节点中，这样理论上在遍历三角形时每次可以跳过一半数量的三角形，进一步提高光线求交的效率。
BVH Construction
BVH的构造过程可以分为3步：
- 初始化每个bounding box的信息，并将它们存储在一个数组中；
- 按照一定的方法将当前节点划分出2个子树；
- 最后将整棵树转换成一个更高效的表示。
在具体构造时首先对当前的所有bounding box取并集从而得到当前步的根节点，然后将这些bounding box划分成两部分再分别构造相应的树作为两个子节点。在构造子树前需要选择进行划分的轴，在一个好的轴上bounding box应该尽可能地分散，因此PBRT选择了bounding box中心范围最大的那个轴来进行划分：
The Surface Area Heuristic
PBRT中介绍了3种不同的划分方法，分别是选择当前范围的中点进行划分、将bounding box等量划分以及使用SAH(surface area heuristic)来进行划分。前两种方法比较容易实现，但在实践中发现这样划分往往会得到两个相互重叠的boudning box，使得在遍历树时的效率会降低。因此在PBRT中默认使用了SAH来进行划分。
SAH在每个节点上考虑把它作为一个叶节点以及对它继续进行划分两种可能性。直接构造叶节点的代价为：
\[\sum_{i=1}^N t_\text{isect} (i)\]而继续划分的代价为：
\[t_\text{trav} + p_A \sum_{i=1}^{N_A} t_\text{isect} (a_i) + p_B \sum_{i=1}^{N_B} t_\text{isect} (b_i)\]其中$t_\text{isect} (i)$为遍历单个三角形的代价；$t_\text{trav}$则是访问当前节点的代价；$p_A$和$p_B$分别是光线穿过左右节点的概率，在BVH中认为它们等于两个bounding box的表面积除以整个节点bounding box的表面积。
直接计算所有可行划分的代价还是太大，因此在PBRT中将所有的bounding box划分到nBuckets=12个bucket中，每个bucket对应一个大的bounding box。然后对这些bucket进行遍历并从中选择划分代价最小的那个来进行划分，对应的代价为对节点继续进行划分的代价。如果进行进行划分的代价小于把整体作为叶节点的代价就继续调用SAH划分左右两个子节点。
Compact BVH
建树的最后一步是把整个树转换成一个更紧凑的表示形式。在PBRT中整棵树的每个节点按照深度优先的顺序存储在内存中：每个内部节点的左节点直接放在它后面，同时每个内部节点还会存储右节点的偏移量。
int BVHAccel::flattenBVHTree(BVHBuildNode *node, int *offset) {
LinearBVHNode *linearNode = &nodes[*offset];
linearNode->bounds = node->bounds;
int myOffset = (*offset)++;
if (node->nPrimitives > 0) {
CHECK(!node->children[0] && !node->children[1]);
CHECK_LT(node->nPrimitives, 65536);
linearNode->primitivesOffset = node->firstPrimOffset;
linearNode->nPrimitives = node->nPrimitives;
} else {
// Create interior flattened BVH node
linearNode->axis = node->splitAxis;
linearNode->nPrimitives = 0;
flattenBVHTree(node->children[0], offset);
linearNode->secondChildOffset =
flattenBVHTree(node->children[1], offset);
}
return myOffset;
}
节点信息存储在结构体LinearBVHNode中，它被设计成32 bytes使得它可以被直接放入一个cache line中优化缓存性能。
struct LinearBVHNode {
Bounds3f bounds;
union {
int primitivesOffset; // leaf
int secondChildOffset; // interior
};
uint16_t nPrimitives; // 0 -> interior node
uint8_t axis; // interior node: xyz
uint8_t pad[1]; // ensure 32 byte total size
};
BVH Traversal
得到BVH树后最重要的应用是从上至下对树进行遍历。我们可以利用递归来实现遍历，但在PBRT中则是通过循环和栈来实现这个过程：
bool BVHAccel::Intersect(const Ray &ray, SurfaceInteraction *isect) const {
if (!nodes) return false;
ProfilePhase p(Prof::AccelIntersect);
bool hit = false;
Vector3f invDir(1 / ray.d.x, 1 / ray.d.y, 1 / ray.d.z);
int dirIsNeg[3] = {invDir.x < 0, invDir.y < 0, invDir.z < 0};
// Follow ray through BVH nodes to find primitive intersections
int toVisitOffset = 0, currentNodeIndex = 0;
int nodesToVisit[64];
while (true) {
const LinearBVHNode *node = &nodes[currentNodeIndex];
// Check ray against BVH node
if (node->bounds.IntersectP(ray, invDir, dirIsNeg)) {
if (node->nPrimitives > 0) {
// Intersect ray with primitives in leaf BVH node
for (int i = 0; i < node->nPrimitives; ++i)
if (primitives[node->primitivesOffset + i]->Intersect(ray, isect))
hit = true;
if (toVisitOffset == 0) break;
currentNodeIndex = nodesToVisit[--toVisitOffset];
} else {
// Put far BVH node on _nodesToVisit_ stack, advance to near node
if (dirIsNeg[node->axis]) {
nodesToVisit[toVisitOffset++] = currentNodeIndex + 1;
currentNodeIndex = node->secondChildOffset;
} else {
nodesToVisit[toVisitOffset++] = node->secondChildOffset;
currentNodeIndex = currentNodeIndex + 1;
}
}
} else {
if (toVisitOffset == 0) break;
currentNodeIndex = nodesToVisit[--toVisitOffset];
}
}
return hit;
}
其中currentNodeIndex记录了当前节点的编号；toVisitOffset记录了待访问节点的数量；而nodesToVisit[64]则记录了待访问节点的编号，它的作用相当于一个栈。在遍历时如果光线与当前节点的bounding box相交则继续向下遍历访问下一层的节点并将另一个节点推入栈中以便将来访问；如果光线不与bounding box相交则通过栈跳转到待访问的节点；当访问到叶节点时则按顺序访问叶节点内的三角形，让光线分别和它们计算是否相交同时更新光线的信息；如果访问节点时发现栈是空的则直接返回。
除此之外在遍历中还利用了光线的方向来选择访问子节点的顺序：如果光线方向都是正方向则优先访问左节点，否则优先访问右节点。
Kd-Tree Accelerator
除了BVH外还可以使用kd树来表示场景中的三角网格。不过在实践中发现使用kd树进行遍历的效率一般会低于BVH，因此目前的各种主流框架都基本不会使用kd树来管理场景。这里就先不介绍kd树相关的内容了，等将来有空的时候再补上。