用于光线跟踪的高并行度表面积启发式（SAH）KD树构建

李建锋 谭耀华 廖胜辉

摘 要：提出一种用于光线跟踪的SAHKD树构建方法，解决当前KD树并行算法并行度不高且效率低的问题.算法首先对所有图元包围盒在三个维度按坐标轴左值排序，得到三维上有序的包围盒索引.然后使用层次遍历构建KD树，根据每个节点包围盒选择要划分的维度，并在当前层生成所有节点在该维度下的候选划分点序列.最后计算每个节点的空间树，在GPU中计算每个候选点的SAH值，选择每个节点的最小SAH值点进行划分.实验中采用4个常用场景进行测试算法性能，并同时比较了当前高效串行与并行算法，结果证明本文提出的算法在生成同等质量KD树的情况下达到对比串行方法4～6倍以及对比并行方法的1.3～1.5倍的计算速度，并且能在线程数成倍增加时达到相近倍数的加速比.

关键词：SAHKD树；空间树；并行计算

中图分类号：TP309.7 文献标志码：A

Abstract：This paper proposed a SAHKD tree construction method for ray tracing to solve the problem of low parallel degree and low efficiency in the existing algorithms. The algorithm first obtains the ordered index on the three latitudes by sorting the primitive bounding boxes according to the left value. Then， according to the node bounding box， we choose the dimension to be divided and generate the candidate partition points of each node under this dimension. Finally， the SAH value of each candidate point was calculated based on the spatial tree in GPU， and the minimum SAH value of each node was selected for partitioning. Four common scenarios were used for testing performance of the algorithm， and compared the efficient serial and parallel GPU algorithm. The results show that the proposed algorithm can achieve 4～6 and 1.3～1.5 times faster than the contrast method when generating the same quality KD tree. When the number of GPU cores fold increases， the speedup ratio reaches a near multiple.

Key words：SAHKDtree； space tree； parallel computing

K维树（KD树）是一类二进制分割树，作为加速结构广泛地应用于图形与查询领域，如光线跟踪中图元与光线的相交检测[1-2]和最近邻查询等.为了高效遍历KD树，通常在构建时使用表面积启发式划分法（SAH），通过最小化每个节点的期望遍历成本选择划分点.如何快速高效地构造SAHKD树一直是研究重點，主要研究方向集中在两种方法上.

一是串行构造方法，最简单的思路是通过全遍历计算SAHKD树所有划分点的SAH值，其计算复杂度为O（N2），由于计算复杂度较高，使得该方法无法应用于大场景.针对该问题，文献[3]和[4]提出将划分点进行预排序，通过迭代的方法依次求出每个划分点的SAH值，把每个节点的时间复杂度降低到O（NlogN），从而使得整棵树的构建复杂度降低至O（Nlog2N）.Wald和HaVRan在文献[3]和文献[5]的基础上，提出只在根节点排序，非根节点通过遍历其父节点序列得到有序序列，使得整棵树的构建复杂度降低至O（NlogN）[5]，此时已经达到了SAHKD树串行构建的理论下限.文献[6]提出一种自调整参数的通用方法，通过随机搜索来得到满足条件的划分点.该方法虽然可以提高构造速度，但无法得到稳定的结果.

二是并行化构建方法， 文献[7]和[8]提出在上层大节点通过一个高性能线程计算，直到每个节点计算量足够小，才将它们分配到多核系统中的构造.具体实践中，在4核CPU系统中能达到单核构造2.5倍的速度[9].Shevtsov等[10]通过并行方法计算中位图元取代SAH的计算来寻找上层的大节点的划分点，虽然在遍历性能较SAH方法有15%的退化，但该方法应用在Larrabee[11]的4路并行光线追踪器时达到3.9倍的构造与遍历的总性能.文献[12]先把图元进行层次分组，之后按组并行构建SAHKD树，在层次分明的场景中，12线程下性能最好时达到文献[4]中算法的130倍.该方法虽然平分了各个线程的计算量，但对于层次不明确的场景无法达到理想效果.自从CUDA与OPENCL等并行架构提出来后，GPU上SAHKD树的构建得到了更多的研究与改进，Kun Zhou[13]提出在大节点直接使用图元的中点划分，对于小节点则进行多路GPU构建.该方法在128线程上达到了单核构建的6～15倍性能，在252K个三角形的场景下达到了16线程的3～6倍.该方法虽然通过顶层退化提高了SAHKD树的构建效率，但在超大型场景遍历时有不同程度的性能损失.文献[14]在前者的基础上加入了CUDA自带的reduction和原子操作，在实际应用中构建效率得到了部分提升，但依然存在顶层退化的问题.

2 并行SAHKD树构造

并行化SAHKD树的构造，现有的并行方式都是基于单个节点，每个线程串行处理一个节点，并单独使用一个线程进行同步.由于上层节点的巨大计算量，通常会使用大小节点的方式分别计算：对于大节点，使用一种退化的方式划分，比如空间中点划分以及图元中点划分等，对于小结点，则通过上述高效的串行方法.理想的并行方法需要基于候选划分点，这样可以实现较高的并行度和每个线程极低的负载.

为了设计一棵基于候选划分点并行的SAHKD树，必须解决两个问题.第一是合理的同步机制.由于SAH值和空间樹的计算分别位于GPU和CPU中，因此计算前后的同步工作变得尤为重要.高效的同步方案应该减少空间分配以及内存与显存间互相拷贝的次数.第二是每个线程使用高效的线面相交检测的方法，由于GPU单核性能并不突出，如果直接使用暴力算法，在大场景中效率会十分低下，从而影响整体的构造性能.

本文提出的方法大致思想如下：首先生成一层每个节点在三个维度上图元包围盒的有序索引，然后根据每个节点包围盒选择一个维度进行划分，生成每个节点在该维度下候选划分点序列，最后一次性在GPU中并行计算一层每个划分点的SAH值，并选择每个节点最优的SAH值的点进行划分.

2.1 算法步骤

本文方法的算法步骤如下.

1）预处理：在SAHKD树的构建之前，对根节点图元进行预处理，其中包括在CPU中计算图元和根节点的AABB包围盒，之后在三个维度按包围盒左值的升序排列所有图元包围盒.

2）按层次递归构建KD树，对于该层的每个节点：

（a）数据准备：在CPU中根据节点包围盒选择划分的维度，生成该维度下的候选划分点序列和空间树，并传入GPU.

（b）计算SAH值：在GPU中利用空间树并行计算出每个划分点的SAH值并传回内存.

（c）图元指派：在CPU中选出最小SAH值对应的划分点，生成两个新节点，在三个维度将图元划入对应的节点中，得到新节点的有序图元包围盒.

2.2 计算候选划分点的SAH值

计算候选划分点的SAH值主要是为每个节点选出最优划分点.由于KD树每次只选择一个维度进行划分，所以第一步是找出要进行划分的维度.为了简便计算，通常对当前节点包围盒进行分析，选择包围盒最长的维度作为要划分的维度.

计算SAH值时需要统计划分点左右两侧图元数量，本文算法通过空间树来加速这一过程.空间树通过有序图元包围盒构造：有序序列可以看成一棵平衡二叉树，而空间树是对平衡二叉树加入一个副值来二次判断.以划分X轴为例，包围盒有两个坐标值Xleft和Xright，用来代表包围盒左右两端，此时加入一个辅助值T，其定义为：

1）如果该节点是叶节点，T定义为Xright；

2）如果该节点是非叶节点，T定义为两个子节点的T值与自己Xright的最大值.

很显然T值需要后序自底向上遍历，使用有序序列转换为平衡二叉树的算法实现T值的计算.为了节省函数调用的开销，本文使用非递归的方式实现.计算出空间树后，每个候选划分点都需要遍历相应节点的空间树求出其SAH值，主要复杂度在计算NL、NR和NP上，其中NL，NR，NP分别为位于候选划分点左，右以及被候选划分点穿过的图元数.过程中对于坐标为pos的候选划分点有：

1）如果当前节点的T值小于pos，则从left到当前节点都在划分点的左边，将NL加上当前区间图元的个数，并计算其右子节点.

2）如果当前节点的Xleft大于pos，则位于该节点右边的节点都在划分点的右边，将NR置为当前节点序号，并计算其左子节点.

3）否则，判断当前节点的图元包围盒与pos的位置关系，相应地修改NL和NP，并遍历其左右子节点.

4）得到NL，NR，NP后计算出SAH值.

上述空间树的遍历方法需要使用辅助栈，然而如果每个线程都使用辅助内存，在大场景下无法正常进行计算.在实际操作中可以使用二叉搜索树的思想进行改造，在构造空间树时得到节点遍历的顺序并设立一个状态标识，不但可以节省显存开销，同时使各个线程计算更加高效.

2.3 图元指派

2.4 复杂度分析

本文方法的特点一是按层级同时计算每一层节点的最优划分点.按层级计算可以使分配内存的次数以及GPU与CPU的交互次数只与层数有关，并且可以使用连续内存空间保存全部需要计算的图元，而不需要使用队列的辅助.在进行同步操作时，只需要在每个维度遍历一次上层节点就能得到下层节点的有序候选序列，避免了多次遍历.二是使用基于空间树的快速计数方法，使每个线程可以根据快速计算出划分点对应的NL、NR和NP，其计算复杂度和划分点穿过的图元数有关，考虑实际场景这个值总是远远小于总的图元数，使得计算速度大大提高.

本节对本文提出方法的复杂度从三个方面进行分析：

1）KD树构造前，需要先对图元包围盒在三个维度进行排序，其复杂度是O（NlogN）.

2）在KD树构造期间，对于每一个非根节点，CPU端首先需要构造出三个维度的有序图元包围盒和节点包围盒.节点包围盒可以根据最优划分点在O（1）的时间得出，将图元顺序指派至新的节点中，需要在三个维度上遍历所有图元包围盒序列一次，复杂度为O（N）.得出有序的图元包围盒后，层次递归每个节点的图元包围盒序列构造该节点的空间树，每一层总的复杂度为O（N）.在GPU端，根据空间树的特性，它需要搜索每个区间，如果区间内所有图元都在划分点的同一侧，则停止搜索这个区间.所以只有该区间包括了被划分点穿过的图元包围盒时，空间树才会在该区间的二分子区间中去搜索.又因为树的层高不会超过（logN），于是每次计算SAH值的复杂度最坏情况下不会超过O（NP+logN）.其中NP是该划分点穿过的图元数，正常场景中无论NP还是logN都是远远小于N的值.

3）整體分析：整体结合来看，每一层复杂度为O（N+NP+logN），而树的高度为（logN），因此该构造方法复杂度仍为O（NlogN）.由于CPU端的任务主要是指派图元与空间树构建，而不涉及任何浮点数计算，所以计算效率很高.处理中排序对算法性能影响很大，考虑到真实场景中一般都会有一个初始状态，我们可以基于该状态下提前离线进行预排序，从而节省初始的排序计算量.

3 实验结果及讨论

为了验证方法的有效性，我们使用统一计算设备架构（Compute Unified Device Architecture：CUDA）框架实现本文提出的算法.算法运行平台为：四核八线程INTEL酷睿E3处理器，GPU型号为GTX750Ti.测试数据为Toys、 Bunny 、Dragon 、buddha四个常用的场景，这些场景包含的图元的数从11 k到1 M不等.在最多64线程下本文算法与3种算法进行比较，分别是Wald提出的串行算法[2]、文献[14]算法，以及不使用空间树的本文算法（朴素GPU算法）做对比，算法时间（并行算法64线程）比较如表1所示.

表1结果表明，本文算法的构造时间在各个场景都优于其他算法，其主要原因是：本文算法是对图元包围盒进行排序，而不是候选点，使得需要排序的序列长度从2N减少到了N，在计算上又通过图元级并行处理，虽然算法均摊下来不如 Wald方法中O（1）简单，但从每个计算核心来看本文算法复杂度为O（NP+logN）优于Wald方法的O（N）.朴素GPU算法虽然不需要进行预排序，但每个线程都需要对当前所有图元进行相交测试，时间复杂度为O（N），复杂度也远高于本文算法，对比于Wald方法，二者复杂度虽然相近，但由于GPU单核心计算能力较低，所以速度较慢.文献[14]中算法虽然是退化的SAHKD树，但在实际场景中依然慢于本文算法，这归根于单节点运算的低并行度.此外，每次划分求出的SAH值是KD树性能的重要指标，但由于后三种方法都是通过寻找最优候选划分点，所以求出来的SAH值都差不多，并不存在遍历性能上的差别.

在具体的场景中，场景Toys的处理结果表明，在小场景下，本文算法与高效的串行算法差距并不大，主要原因是本文算法中指派图元是在CPU中完成，而少量的浮点数计算不会明显影响串行算法的计算时间.朴素GPU算法不需要最开始的预排序，且在小场景下空间树的效果并不明显，最终两者的结果相近.

其他三个场景中随着图元数量的增加，本文算法达到4～6倍于串行算法的时间效率，且随着图元数变大稳步增加.串行算法在每层构建时通常需要O（N）级别的浮点数计算，当数据量大时计算量远超过条件判断并且越来越明显.空间树的效果也在大场景中体现出来，对比于朴素GPU算法，由于条件判断较多消耗了大量线程的计算能力，而空间树可以大大降低这一过程，代价仅需在CPU端对有序的图元进行一次遍历.

构建SAHKD树过程中，随着层次的加深，每个节点所包含的图元数越少，但遍历时需要更多的相交测试才能达到叶节点，为了平衡这一现象，在实际应用中都会制定相应的早停策略.图1列出buddha场景下SAHKD树计算1～7层时，通过公式（7）计算得到的算法时间效率（文献[14]的算法在上层使用中位图元划分，不适用于这类分析）.

4 总 结

本文提出了一种使用CPU和GPU并行构造一棵基于层次的SAHKD树，解决了当前SAHKD树并行度较低和需要区分大小节点的问题.由于KD树的构建是一个静态过程，本文所有测试都是在静态场景中进行.但是真实场景许多是动态的，以往研究者通常是对静态场景进行重复计算来模拟动态场景.我们接下来的研究方向是如何在动态场景中快速生成SAHKD树，基于本文算法在预排序阶段进行相关的优化，期望达到实时光线跟踪的效果.

参考文献

[1] HUSSAIN S， GRAHN H. Fast kdtree construction for 3Drendering algorithms like ray tracing[C]//Lecture Notes in Computer Science. Berlin： Springer， 2004，4842：681-690.

[2] WALD I. Realtime ray tracing and interactive global illumination[J]. Information Technology，2004，48 （4）：242-245.

[3] MACDONALD J D， BOOTH K S. Heuristics for raytracing using space subdivision[J]. Visual Computer， 1990，6 （3）：153-165.

[4] SZECSI L. An effective implementation of the kdtree[M]. Rockland： Charles River Media Inc，2003：315-326.

[5] WALD I， HAVRAN V. On building fast kdtrees for ray tracing and ondoing that in O（NlogN）[C]//Symposium on Interactive Ray Tracing 2006（RT）.Salt Lake City：IEEE，2006 ：61-69.

[6] TILLMANN M， PFAFFE P， CKAAG W F.Onlineauto tuning of parallel SAH kdtrees[C]//International Parallel and Distributed Processing Symposium（2016）.Chicago：IEEE，2016：628-637.

[7] LEE J， SHIN Y.Realtime ray tracing on coarsegrained reconfigurable processor[J].International Conference on Fieldprogrammable Technology，2014， 149 （3）：192-197.

[8] POPOV S， GUNTHER J. HP seidel.experiences with streaming construction of SAH kdtrees[C]//Symposium on Interactive Ray Tracing （2006）.Salt Lake City：IEEE，2006：89-94.

[9] HUNT W，MARK W R，STOLL G.Fast kdtree construction with an adaptive errorbounded heuristic[C]//Symposium on Interactive Ray Tracing （2006）.Salt Lake City：IEEE，2006：81-88.

[10]SHEVTSOV M， SOUPIKOV A， KAPUSTIAN A. Highly parallel fast kdtree construction for interactive ray tracing of dynamic scenes[J].Computer Graphics Forum， 2007，26 （3）：395-404.

[11]SPJUT J， BOULOS S， KOPTA D. A multithreaded architecture for realtime ray tracing[J].Symposium on Application Specific Processors，2008，28 （12） ：108-114.

[12]KANG Y S， NAH J H， PARK W C. gkDtree： A groupbased parallel update kdtree for interactive ray tracing[J].Journal of Systems Architecture the Euromicro Journal ， 2013，59 （3）：166-175.

[13]ZHOU K， HOU Q， WANG R.Realtime kdtree construction on graphics hardware[J].Acm Transactions on Graphics，2008，27 （5）：1-11.

[14]CHANGH B，SEOH W，IHM I. On the efficient Implementation of a realtime kdtree construction algorithm[M].Singapore： Springer，2015：207-219.

[15]MACDONALD J D， BOOTH K S.Heuristics for raytracing using space subdivision[J].Visual Computer，1990，6 （3）：153-166.

[16]SANTAL L A.Integral geometry and geometric probability[M]. London： AddisonWesley Pub Co，1976 ：91-113.