BlueRose
BlueRose
文章97
标签28
分类7
回到首页
文章归档
关于我
© 2024 BlueRose Powered by Hexo & BlueRose
PBRT笔记(2)——BVH

PBRT笔记(2)——BVH

2018年12月31日 2.7k 字 大概 13 分钟
  1. 1. BVH
  2. 2. 分割
  3. 3. 紧凑BVH树

BVH

构建BVH树分三步：

  1. 计算每个图元的边界信息并且存储在数组中
  2. 使用指定的方法构建树
  3. 优化树，使得树更加紧凑 
    //BVH边界信息，存储了图元号，包围盒以及中心点
struct BVHPrimitiveInfo {
 BVHPrimitiveInfo() {}
 BVHPrimitiveInfo(size_t primitiveNumber, const Bounds3f &bounds)
     : primitiveNumber(primitiveNumber),
       bounds(bounds),
       centroid(.5f * bounds.pMin + .5f * bounds.pMax) {}
 size_t primitiveNumber;
 Bounds3f bounds;
 Point3f centroid;
};

    分割

    使用，子图元中质心距离最大的轴向作为分割方向。（另一种方法是尝试所有轴，之后再选择效果最好的那个轴作为分割方向。但是在实践中发现当前方案也有着不错的效果）
    int nPrimitives = end - start;
//只有一个图元，则生成叶子节点并且返回
if (nPrimitives == 1) {
 int firstPrimOffset = orderedPrims.size();
 for (int i = start; i < end; ++i) {
     int primNum = primitiveInfo[i].primitiveNumber;
     orderedPrims.push_back(primitives[primNum]);
 }
 node->InitLeaf(firstPrimOffset, nPrimitives, bounds);
 return node;
}

    int mid = (start + end) / 2;
//如果多个片元的组成的聚合体是0空间体，则生成叶子节点，并且返回（这是一种不寻常的现象）
if (centroidBounds.pMax[dim] == centroidBounds.pMin[dim]) {
 // Create leaf _BVHBuildNode_
 int firstPrimOffset = orderedPrims.size();
 for (int i = start; i < end; ++i) {
     int primNum = primitiveInfo[i].primitiveNumber;
     orderedPrims.push_back(primitives[primNum]);
 }
 node->InitLeaf(firstPrimOffset, nPrimitives, bounds);
 return node;
}

    int mid = (start + end) / 2;
if (centroidBounds.pMax[dim] == centroidBounds.pMin[dim]) {
<创建叶子节点>
} else {
 <使用对应方法切割图元>
 //使用递归构造BVH树
 node->InitInterior(dim,
 recursiveBuild(arena, primitiveInfo, start, mid,totalNodes,orderedPrims),
 recursiveBuild(arena, primitiveInfo, mid, end,totalNodes, orderedPrims));
}

    //中间对半切的方法
case SplitMethod::Middle: {
 Float pmid =
     (centroidBounds.pMin[dim] + centroidBounds.pMax[dim]) / 2;
 BVHPrimitiveInfo *midPtr = std::partition(
     &primitiveInfo[start], &primitiveInfo[end - 1] + 1,
     [dim, pmid](const BVHPrimitiveInfo &pi) {
         return pi.centroid[dim] < pmid;
     });
 mid = midPtr - &primitiveInfo[0];
 if (mid != start && mid != end) break;
}
    这里遇到指针相减的操作

    如果两个指针指向同一个数组,它们就可以相减,其结果为两个指针之间的元素数目。同理+1则表示内存偏移一个该类型空间。

[end - 1] + 1是会为了回避容器越界访问的问题，但是取了地址再偏移

end操作返回的迭代器指向容器的“末端元素的下一个”，指向了一个不存在的元素

但是如果有多个图元的边界盒处于与左右节点的边界盒重叠的状态，那分割很可能会失败。SplitMethod::EqualCounts排序速度会更快些

case SplitMethod::EqualCounts: {
    // Partition primitives into equally-sized subsets
    mid = (start + end) / 2;
    std::nth_element(&primitiveInfo[start], &primitiveInfo[mid],
                     &primitiveInfo[end - 1] + 1,
                     [dim](const BVHPrimitiveInfo &a,
                           const BVHPrimitiveInfo &b) {
                         return a.centroid[dim] < b.centroid[dim];
                     });
    break;
}

case SplitMethod::SAH:
default: {
    //图元数量较少时没必要使用SAH
    if (nPrimitives <= 2) {
        // Partition primitives into equally-sized subsets
        mid = (start + end) / 2;
        std::nth_element(&primitiveInfo[start], &primitiveInfo[mid],
                         &primitiveInfo[end - 1] + 1,
                         [dim](const BVHPrimitiveInfo &a,
                               const BVHPrimitiveInfo &b) {
                             return a.centroid[dim] <
                                    b.centroid[dim];
                         });
    } else {
        // Allocate _BucketInfo_ for SAH partition buckets
        PBRT_CONSTEXPR int nBuckets = 12;
        BucketInfo buckets[nBuckets];

        //计算出当前图元的质心处于第几个桶，centroidBounds为当前节点中所有图元的边界盒
        //之后进行BucketInfo统计，并且调整对应桶的边界盒
        for (int i = start; i < end; ++i) {
            int b = nBuckets *
                    centroidBounds.Offset(
                        primitiveInfo[i].centroid)[dim];
            if (b == nBuckets) b = nBuckets - 1;
            CHECK_GE(b, 0);
            CHECK_LT(b, nBuckets);
            buckets[b].count++;
            buckets[b].bounds =
                Union(buckets[b].bounds, primitiveInfo[i].bounds);
        }

        Float cost[nBuckets - 1];
        for (int i = 0; i < nBuckets - 1; ++i) {
            Bounds3f b0, b1;
            int count0 = 0, count1 = 0;
            //第一个桶到[i]桶，所有的图元数与边界盒
            for (int j = 0; j <= i; ++j) {
                b0 = Union(b0, buckets[j].bounds);
                count0 += buckets[j].count;
            }
            //[i+1]桶到最后一个桶，所有的图元数与边界盒
            for (int j = i + 1; j < nBuckets; ++j) {
                b1 = Union(b1, buckets[j].bounds);
                count1 += buckets[j].count;
            }
            //bounds变量为所有图元的边界盒
            //通过面积模型计算，相交开销设为1
            cost[i] = 1 +
                      (count0 * b0.SurfaceArea() +
                       count1 * b1.SurfaceArea()) /
                          bounds.SurfaceArea();
        }

        //计算出最小消耗的切割位置与消耗的量
        Float minCost = cost[0];
        int minCostSplitBucket = 0;
        for (int i = 1; i < nBuckets - 1; ++i) {
            if (cost[i] < minCost) {
                minCost = cost[i];
                minCostSplitBucket = i;
            }
        }

        Float leafCost = nPrimitives;
        //图元数超过最大值（255）或者最小消耗小于所有图元都创建叶子节点的消耗（切割具有良好效果的情况）
        //使用partition算法，按使用SAH找到最小消耗的切割位置，对桶进行排序，之后进行下一次遍历
        //不然就直接生成节点，结束遍历
        if (nPrimitives > maxPrimsInNode || minCost < leafCost) {
            BVHPrimitiveInfo *pmid = std::partition(
                &primitiveInfo[start], &primitiveInfo[end - 1] + 1,
                [=](const BVHPrimitiveInfo &pi) {
                    int b = nBuckets *
                            centroidBounds.Offset(pi.centroid)[dim];
                    if (b == nBuckets) b = nBuckets - 1;
                    CHECK_GE(b, 0);
                    CHECK_LT(b, nBuckets);
                    return b <= minCostSplitBucket;
                });
            mid = pmid - &primitiveInfo[0];
        } else {
            // Create leaf _BVHBuildNode_
            int firstPrimOffset = orderedPrims.size();
            for (int i = start; i < end; ++i) {
                int primNum = primitiveInfo[i].primitiveNumber;
                orderedPrims.push_back(primitives[primNum]);
            }
            node->InitLeaf(firstPrimOffset, nPrimitives, bounds);
            return node;
        }
    }
    break;
}

紧凑BVH树

使用SAH有两个缺点：

  1. 花费过多时间在使用SAH构建树上。
  2. 自上而下的BVH树的构建很难并行化。有一个解决方法就是构建若干个独立子树，不过这反过来限制了并行的可伸缩性。这也是GPU渲染所要面对的问题。

Linear bounding volume hierarchies (LBVHs) 就是为了解决这个问题而开发的。
LinearBVHNode大小为32字节，以满足32位内存对齐要求。
LBVH是基于莫顿码，讲原本的多维数据转换为排序过的一维的数据。
也就是将树中节点的相对位置按照规律排序：
如果用uint8来表示一个二叉树的2维状态：

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
};

hierarchical linear bounding volume hierarchy (HLBVH)
讲质心坐标转化为莫顿码，之后统计对应莫顿码（桶中）图元数量，之后以莫顿码作为index，将图元放入容器的对应位置，从而完成排序（使用了基数排序）。

ParallelFor([&](int i) {
    // Initialize _mortonPrims[i]_ for _i_th primitive  
    //mortonScale=1024;
    PBRT_CONSTEXPR int mortonBits = 10;
    PBRT_CONSTEXPR int mortonScale = 1 << mortonBits;
    mortonPrims[i].primitiveIndex = primitiveInfo[i].primitiveNumber;
    //因为bounds.Offset返回的是[0,1]区间的百分比，所以为了转换成莫顿码需要再乘以一个比例系数
    //PBRT使用int32来存储莫顿码，因为需要存储x、y、z三个维度，所以每个维度占用10个bit，所以比例系数为10，对于二进制莫顿码来说乘以10等于左移10个位，也就是2^10=1024
    Vector3f centroidOffset = bounds.Offset(primitiveInfo[i].centroid);
    //以边界盒的min为零点建立坐标系，使用质心位置*莫顿码缩放值，计算莫顿码
    //EncodeMorton3，使用LeftShift3()分别计算x、y、z，之后再将y、z分别往左偏移1与2位
    //LeftShift3()参看书中的图就可以明白了
    mortonPrims[i].mortonCode = EncodeMorton3(centroidOffset * mortonScale);
}, primitiveInfo.size(), 512);

之后以索引对莫顿码进行排序（为了追求效率而没有选择std::sort）,这里如果不明白基数排序，就很难看懂。

for (int pass = 0; pass < nPasses; ++pass) {
//pass如果各个位都为1，着in为tempVector的引用，否则则为v
//每次循环out与in都进行互换，将上一次排序结果接着进行排序（第一次直接用外部未排序的容器，之后将每个pass都排序一遍，所有莫顿码即排序完成）
std::vector<MortonPrimitive> &in = (pass & 1) ? tempVector : *v;
std::vector<MortonPrimitive> &out = (pass & 1) ? *v : tempVector;
}

//一共64种可能的莫顿码
PBRT_CONSTEXPR int nBuckets = 1 << bitsPerPass;
int bucketCount[nBuckets] = {0};
PBRT_CONSTEXPR int bitMask = (1 << bitsPerPass) - 1;
for (const MortonPrimitive &mp : in) {
    int bucket = (mp.mortonCode >> lowBit) & bitMask;
    CHECK_GE(bucket, 0);
    CHECK_LT(bucket, nBuckets);
    //取得当前的pass的莫顿码，并且进行统计
    ++bucketCount[bucket];
}
//计算每个桶到第一个桶的莫顿码总量，也就是index的偏移量，毕竟当前位置的index+该位置桶内元素数量，就等于下一个桶到第一个桶的偏移量。
int outIndex[nBuckets];
outIndex[0] = 0;
for (int i = 1; i < nBuckets; ++i)
    outIndex[i] = outIndex[i - 1] + bucketCount[i - 1];
//通过莫顿码讲图元节点放到对应的位置，从而完成排序。++是为了偏移放置下一个该桶元素的位置
for (const MortonPrimitive &mp : in) {
    int bucket = (mp.mortonCode >> lowBit) & bitMask;
    out[outIndex[bucket]++] = mp;
}

    // 寻找每个小数中图元的间隔
    std::vector<LBVHTreelet> treeletsToBuild;
    for (int start = 0, end = 1; end <= (int)mortonPrims.size(); ++end) {
#ifdef PBRT_HAVE_BINARY_CONSTANTS
      uint32_t mask = 0b00111111111111000000000000000000;
#else
      uint32_t mask = 0x3ffc0000;
#endif
      //遍历所有莫顿码高位不同的图元
      if (end == (int)mortonPrims.size() ||
            ((mortonPrims[start].mortonCode & mask) !=
             (mortonPrims[end].mortonCode & mask))) {
            // Add entry to _treeletsToBuild_ for this treelet
            int nPrimitives = end - start;
            int maxBVHNodes = 2 * nPrimitives;
            BVHBuildNode *nodes = arena.Alloc<BVHBuildNode>(maxBVHNodes, false);
            treeletsToBuild.push_back({start, nPrimitives, nodes});
            start = end;
        }
    }

之后对16个区块进行构建子树
emitLBVH

CHECK_GT(nPrimitives, 0);
//如果图元小于一定数量或者 则创建叶子节点
if (bitIndex == -1 || nPrimitives < maxPrimsInNode) {
    (*totalNodes)++;
    BVHBuildNode *node = buildNodes++;
    Bounds3f bounds;
    int firstPrimOffset = orderedPrimsOffset->fetch_add(nPrimitives);
    for (int i = 0; i < nPrimitives; ++i) {
        int primitiveIndex = mortonPrims[i].primitiveIndex;
        orderedPrims[firstPrimOffset + i] = primitives[primitiveIndex];
        bounds = Union(bounds, primitiveInfo[primitiveIndex].bounds);
    }
    node->InitLeaf(firstPrimOffset, nPrimitives, bounds);
    return node;
} else {
// 从高位开始分割 （29-12）
    int mask = 1 << bitIndex;
    //如果所有图元都集中在分割的一边，则开始分割下一个子树
    if ((mortonPrims[0].mortonCode & mask) ==
        (mortonPrims[nPrimitives - 1].mortonCode & mask))
        return emitLBVH(buildNodes, primitiveInfo, mortonPrims, nPrimitives,
                        totalNodes, orderedPrims, orderedPrimsOffset,
                        bitIndex - 1);

    // Find LBVH split point for this dimension
    int searchStart = 0, searchEnd = nPrimitives - 1;
    while (searchStart + 1 != searchEnd) {
        CHECK_NE(searchStart, searchEnd);
        int mid = (searchStart + searchEnd) / 2;
        //用二分法在这个轴线寻找分割点
        if ((mortonPrims[searchStart].mortonCode & mask) ==
            (mortonPrims[mid].mortonCode & mask))
            searchStart = mid;
        else {
            CHECK_EQ(mortonPrims[mid].mortonCode & mask,
                     mortonPrims[searchEnd].mortonCode & mask);
            searchEnd = mid;
        }
    }
    int splitOffset = searchEnd;
    CHECK_LE(splitOffset, nPrimitives - 1);
    CHECK_NE(mortonPrims[splitOffset - 1].mortonCode & mask,
             mortonPrims[splitOffset].mortonCode & mask);

    //因为已经用莫顿码排序过，所以这里直接递归分割
    (*totalNodes)++;
    BVHBuildNode *node = buildNodes++;
    BVHBuildNode *lbvh[2] = {
        emitLBVH(buildNodes, primitiveInfo, mortonPrims, splitOffset,
                 totalNodes, orderedPrims, orderedPrimsOffset,
                 bitIndex - 1),
        emitLBVH(buildNodes, primitiveInfo, &mortonPrims[splitOffset],
                 nPrimitives - splitOffset, totalNodes, orderedPrims,
                 orderedPrimsOffset, bitIndex - 1)};
    int axis = bitIndex % 3;
    node->InitInterior(axis, lbvh[0], lbvh[1]);
    return node;
}

所有子树都创建后，buildUpperSAH将会构建所有子树的BVH。这里的操作和SAH差不多。之后就是flattenBVHTree，以深度优先顺序，对所有节点进行排序。
以下是便利过程：

 while (true) {
 const LinearBVHNode *node = &nodes[currentNodeIndex];
    //是否与当前节点的边界盒相交
    if (node->bounds.IntersectP(ray, invDir, dirIsNeg)) {
        if (node->nPrimitives > 0) {
            //如果是叶子节点，就对节点下的所有图元进行相交测试
            for (int i = 0; i < node->nPrimitives; ++i)
                if (primitives[node->primitivesOffset + i]->Intersect(
                        ray, isect))
                    hit = true;
            //如果已经遍历了另外所需遍历的节点，则停止循环，不然设置下一个循环NodeIndex
            if (toVisitOffset == 0) break;
            currentNodeIndex = nodesToVisit[--toVisitOffset];
        } else {
            //如果是子树节点，先判断方向正负，是正就访问左节点，并且将右节点放入需要遍历的数组中，待之后的循环进行相交测试
            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];
    }
}
图形学
本文作者：BlueRose
本文链接：http://www.blueroses.top/2018/12/31/pbrt-bi-ji-2-bvh/
版权声明：本文采用 CC BY-NC-SA 3.0 CN 协议进行许可