王育坚，廉腾飞，吴明明，高 倩

(北京联合大学 信息学院，北京 100101)

基于八叉树与KD树索引的点云配准方法

王育坚，廉腾飞，吴明明，高 倩

(北京联合大学 信息学院，北京 100101)

针对点云配准算法中KD树多维查询效率较低的问题，提出一种基于八叉树和KD树多层索引结构的点云配准方法。首先为模型点云数据建立八叉树全局索引，然后在八叉树叶子结点构建局部数据的KD树索引。对传统的ICP点云配准算法进行改进，通过叶子结点的全局索引值快速定位局部点云数据块，利用局部KD树索引加快最近点的搜索，计算最近点时利用欧氏距离阈值、点对距离差值和法向量阈值剔除部分噪声点。实验表明，改进算法提高了点云配准的效率和精度。

点云配准；八叉树；KD树；ICP算法

三维重建在计算机视觉、虚拟现实、3D打印和逆向工程等方面有着广泛的应用，点云数据配准方法的优劣直接影响三维重建的效果。利用三维扫描仪多次扫描物体表面，得到不同视角下的三维点云数据，将扫描得到的多片点云数据进行配准，即可得到完整的模型[1]。点云配准包括粗配准和精配准，粗配准的目的是为精配准提供良好的初值，缩小相邻点云之间的旋转误差和平移误差。精配准的目的是找到最合适的旋转矩阵和平移矩阵，使得相邻点云配准的误差最小。

目前应用最广泛的点云配准算法是迭代最近点(Iterative Closest Point，ICP)算法[2]及其改进算法。ICP算法作为一种基于纯粹几何模型的配准方法，人们提出了多种改进方案[3]，改进的基本思路主要体现在点云数据模型设计和最近点搜索策略两个方面[4-5]。例如，在点云数据模型设计方面，利用Delaunay三角剖分[6]、KD树[7]等对ICP算法进行改进。在实际应用过程中，对于数据量很大的点云数据，ICP算法配准耗时较长，影响了算法效率。本文基于八叉树和KD 树的多层索引结构，提出一种点云精配准的改进方法。

1 点云数据结构

1.1 八叉树结构和KD树结构

八叉树结构是一种规则的数据结构，利用树形结构对模型进行递归，按X，Y，Z3个不同方向，将所要表示的三维空间实体分割成8个大小相等的子立方体。然后根据每个子立方体中所含的目标来决定是否对子立方体继续进行8等分的划分，一直划分到每个子立方体被一个目标所充满，或没有目标，或其大小已成为预先规定的体素为止。八叉树分解是将三维空间实体逐级分解，最终形成八叉树体素表示的结构。八叉树的主要优点是可以方便地实现物体的并、交、差等集合运算，适用于不同形状物体的建模。

KD树是一种把二叉查找树推广到多维数据的结构，实现多维空间数据的组织和存储。KD树利用超平面把一个空间划分成多个不相交的子空间，每一层都将所包含的空间分成两个子空间，顶层结点按一维划分，下一层结点按另一维划分，KD树所有维的属性在层间循环。任何一个非叶子结点的左右子树也是KD树，若结点的左子树不为空，则左子树上所有结点第d维的值均小于根结点第d维的值；若结点的右子树不为空，则右子树上所有结点第d维的值均大于等于根结点第d维的值。

KD树每一个结点划分结束的条件是结点中只包含一个数据或少于设定的上限为止。KD树可以用来建立多维空间数据集或数据块的索引。利用KD树进行数据查询时，每一步结点的条件判断只要比较其中的一个维。通过交替比较不同维的属性值，可以快速查找某个数据点的邻域，不需要知道数据之间的任何拓扑关系。

1.2 多层索引结构设计

KD树通过左右孩子指针建立数据关系，索引指针数据占据了大量的内存空间。KD树采用一分为二的分割方式，由于点云数据量巨大，使得树的深度很大，增加了数据查找的时间，但影响KD树搜索效率的主要因素是回溯[8]。八叉树结构规则统一，树的深度大大降低，对于精确数据点查找，其性能较高[9]。但八叉树的动态性较差，在数据点集分布不均匀的情况下，树的平衡性不好。

考虑到八叉树数据存储结构特点和KD树搜索有效性，提出一种八叉树与KD树相结合的多层索引结构。设待配准点云为P，模型点云为Q。基于模型点云Q建立多层索引结构，在上层采用八叉树结构存储和管理全局的点云数据，在下层利用KD树组织和存储局部的点云数据，每个局部KD树索引的信息都保存在八叉树末端关联的叶结点中。全局八叉树和局部KD树的多层索引结构如图1所示，其中的虚线框内表示八叉树叶结点与KD树的对应关系。

图1 多层索引结构

根据这种多层索引结构，首先自上而下采用八叉树结构对模型点云进行空间分割。根据模型点云的最小和最大坐标，得到一个包围整个模型的立方体，该立方体作为八叉树的根结点。将该立方体分解成8个子立方体，作为8个孩子结点。根据配准精度要求确定八叉树的分割参数，得到每个叶子结点点云所包含最大的点数。当叶子结点内的点云数目小于最大点数时，结束八叉树的分解，不必像传统八叉树那样直到结点只有一个数据点或有限的数据点时才结束分解。

当利用八叉树分割点云模型生成叶子结点后，在八叉树非空的叶子结点设置指针值，指向关联的KD树。KD树将八叉树叶子结点中的点云数据一分为二，划分为两个子空间，然后再对每一个子空间进一步进行递归划分，最后得到一棵完整的KD树，即所有含有点云数据的叶子结点都建立关联的KD树。创建KD树采用一种改进的方法，即先对3个轴都进行采样[7]，然后选择合适的轴进行划分，以保证每次划分都能找到近似最优的分割位置。

全局八叉树是第一层索引，局部KD树是建立在八叉树叶子结点之下的第二层索引。通过分割参数设定八叉树叶子结点的数据量，如果设置的分割参数很大，则叶子结点关联的局部点云数据很多，KD树的深度较大，会降低KD树的搜索效率；反之，如果将结点划得很细，虽然KD树的深度降低，但增加了边界数据点的数量，势必影响配对的精度。因此，生成八叉树时要综合考虑模型的空间形状、点元数据规模和配准精度的要求，通过实验分析得出接近于最佳搜索的分割参数。

与传统的KD树相比，多层索引结构需要额外的空间存储八叉树，但由于八叉树没有分割到只有一个数据点才结束，八叉树的深度是可以控制的，需要的附加空间并不大。并且，结点之间的层次关系可以根据结点编码得到，配准时定位八叉树的叶子结点并不需要回溯，因此，八叉树每个结点可以不存储其父结点的指针，减少了模型的存储空间。

多层索引结构的八叉树结点用一个属性值表示叶子结点与KD树的关联，若八叉树叶子结点的属性值为“KD”，表示该叶子结点链接一个局部的KD树。通过叶子结点与其对应KD树的根结点的关联，形成局部KD树索引。搜索时根据待配准点的坐标和子立方体空间属性值快速完成最近点的自适应定位，最近点的搜索局限于数据点所在八叉树叶子结点的包围盒内，避免了回溯八叉树。

将索引代价较小的八叉树与KD树相关联，利用空间分块策略设计多层索引，这种索引结构在保证满足配准精度的前提下，能够大大优化最近点的搜索性能。建立全局八叉树和局部KD树多层索引结构的算法如下：

1)根据模型点云Q的大小选取一个包围盒立方体，将该立方体作为八叉树的根结点。

2)若当前结点包含的数据点数大于规定的点数，采用深度优先策略对该结点进行八叉树递归分解，将立方体分解成8个子立方体。

3)若当前结点包含的点数小于或等于规定的点数，停止分解，该结点为八叉树的叶子结点，生成结点数据域，设置结点的“KD”属性值。

4)针对八叉树叶子结点包含的局部点云数据创建KD树。根据3个轴的采样选择划分的坐标轴，将八叉树叶结点的点云数据一分为二，划分为两个平面；再对每个平面进一步递归划分，直到生成KD树的每一个结点。生成KD树后，建立八叉树叶子结点与KD树的索引。

5)将当前指针指向八叉树的下一个结点，重复步骤2)，直到处理完所有八叉树结点。

2 改进的点云配准算法

2.1 ICP算法的最近点搜索

ICP算法根据一定的准则确立对应点集P与Q，通过最小二乘法迭代计算最优的坐标变换，即旋转矩阵和平移矢量，将一个坐标系下的点数据变换到另一个坐标系，并使得误差函数最小。点云配准可以看作是求解变换矩阵的过程，决定ICP算法坐标变换的关键在于能否在模型点云Q中准确、快速地找到P中待配准点的最近点。因此，影响ICP算法效率和精度的主要因素是最近点的搜索方法[10]。

研究者基于最近点搜索提出多种ICP的改进方法。文献[11]提出基于重叠区域的近似KD树，在搜索重叠区域的子结点时不回溯，提高了搜索效率和精度。文献[12]基于点云单应性假设提出一种采用点和面对应的ICP改进算法，利用单应性剔除其余点对，算法具有较好的稳健性和收敛性。文献[13]在K近邻搜索中使用OKDT正交搜索树，利用主成分分析方法计算点集方差最大的正交轴方向，按照优化KD树划分方法划分子树，当点云数据量很大时该方法具有较好的搜索性能。此外，文献[14]根据模型点集的各维方差按坐标轴排序，能够快速搜索最近点搜索范围的边界。文献[15]采用点云欧氏距离阈值和方向矢量夹角阈值去噪。文献[16]提出曲率约束与对应点距离约束相结合的自适应噪点剔除策略。文献[5]首先采用中心重合法实现点云数据的粗配准，然后利用KD树快速搜索最近点对，完成点云数据的精配准。上述算法都在不同情况下提高了搜索的效率或精度。

2.2 基于KD树最近邻搜索算法的改进

通常情况下，计算P中每个点p的最近点，需要计算p与Q中所有点的欧氏距离。局部KD树是针对八叉树叶子结点的点云数据块建立的，通过模型点云的八叉树和KD树多层索引可以快速定位p的最近点所在的数据块，可以避免搜索所有点。最近点的搜索过程分为两步：第一步搜索全局索引层，通过搜索模型点云八叉树，定位最近点所在的点云数据块，即找到包含最近点的八叉树叶子结点；第二步，根据八叉树叶子结点中存储的信息找到对应的KD树根结点，即存储局部数据块的KD树，然后在局部点云数据块中搜索最近邻域点。

影响KD树最近邻搜索效率的主要因素是回溯，八叉树结点存储了相关的数据信息，不需要回溯就可以定位包含最近邻域点的局部点云块。八叉树中的每一个结点都对应一个空间包围盒，根据待配准点坐标pi(xi,yi,zi)与包围盒的空间位置和大小，确定需要继续搜索的八叉树叶子结点。设最近邻域点qj(xj,yj,zj)所在子立方体的空间索引值为(a,b,c)，与子立方体对应的结点的八进制编码为q=qn-1…qj…q1q0，qj(j=0, 1, … ,n-1 )表示叶结点到根结点的路径。根据式(1)可以得到八叉树中最近邻域点所在子立方体的全局索引值。

(1)

在局部KD树搜索最近邻域点时，对于位于分割子体边界上的点，搜索的结果有可能是错误的，可以将边界点作为噪声点处理，由于点云数据量巨大，这样处理提高了配准速度，而对精度没有太大的影响。研究表明，如果允许有少量错误的搜索结果，KD树的搜索效率会得到很大的提高[10]。

改进算法利用欧氏距离阈值剔除噪声点[15]。通过模型点云八叉树叶子结点关联的KD树，搜索与点pi欧式距离最近的3个点s1，s2和s3；若pi与s1，s2，s3构成的平面的距离超出阈值E，则剔除该对应点；否则，以这3个点中距离pi最近的点sj作为对应点。阈值为E=c*d，其中c为控制系数，d为点云中相邻点间的平均距离。此外，边界噪声点也有可能形成多个对应点对，可以比较两组对应点对p1和q1，p2和q2，若出现式(2)情况，即它们之间的差超过F，则视为噪声点，也予以剔除。

|dist(p1-p2)-dist(q1-q2)|≥F.

(2)

利用欧氏距离阈值虽然剔除了大量的噪声点，但对应点对仍然存在噪声点。考虑到待配准点云和模型点云虽然处于不同的坐标系，但其空间拓扑关系应该一致，点云之间除了有平移量，还有旋转量。因此，对匹配点对的法向量夹角设置一个阈值，以进一步剔除错误的点对。

(3)

应用最小二乘法，可以得到以下3*3矩阵A。

(4)

可以证明，A的最小特征值对应的特征向量即可作为法向量ni的近似值。

通过以上方法得到两个对应点集中各点的法向量。对于任意对应点对pi和qj，它们的法向量分别为ni和nj。两者的法向量差别越大，夹角的余弦值越小，即ni·nj就越小。因此，根据式(5)对向量夹角余弦设置阈值，即将法向量的乘积小于G的点对视为噪声点，予以剔除。

ni·nj≥G.

(5)

2.3 算法步骤

改进的点云配准算法的主要步骤如下：

1)针对模型点云Q建立八叉树与KD树多层索引结构。

2)迭代初始化：选择初始目标点集P0=P，设定最大迭代次数Kmax，给定法向阈值V。

3)根据待配准点云Pk中的每个点pi的坐标pi(xi,yi,zi)和八叉树子立方体包围盒的空间位置及大小(x,y,z,l)，在八叉树中定位最近点qi所在的叶子结点。

4)通过八叉树的叶子结点找到包含局部点云的KD树，基于KD树搜寻数据点集Pk中每一个点pi的最近点，得到对应点集Qk。

5)利用欧氏距离阈值和点对距离差值剔除部分噪声点。

6)根据法向量阈值剔除错误的匹配点对。

7)利用四元数法对式(6)进行最小化，求出旋转矩阵Rk和平移矢量Tk。

(6)

8)根据旋转矩阵Rk和平移矢量Tk得到新的数据点集：Pk+1=RkPk+Tk。

9)R=RkR,T=RkT+Tk，重复进行步骤3)～8)，直至前一次最近点之间的距离与后一次最近点之间的距离满足条件：dk-dk+1

10)利用变换矩阵参数R和T将初始目标点云数据变换到参考点云所在的坐标系，完成点云数据的配准。

设模型点云Q有N个数据点，待配准点云P有M个数据点，模型点集中数据点的总数N与八叉树一个叶子结点包含的数据点之比为K(分割参数)，模型点云八叉树深度为h。可以算出，八叉树每个叶子结点包含的点数为N/K，叶子结点数为K，将八叉树简化为满八叉树处理，可以计算出八叉树深度为

h≥|log8K|+1.

(7)

对于有M个数据点的待配准点云P，对应八叉树的搜索时间为O(Mlog8K)；在局部模型点云块KD树中搜索最近点，时间为O(Mlog2(N/K))。因此，搜索最近邻的总时间为

O(Mlog8K+O(Mlog2(N/K)).

(8)

研究表明[14]，在实例随机分布的情况下，KD树最近邻的深度优先搜索的时间为O(Mlog2N)，在最差回溯情况下的时间为O(3MN2/3)。因此，在保证配准精度的前提下，选择合适的K值，满足式(9)就能保证改进后的算法优于随机条件下的传统KD树最近邻搜索算法。

Mlog8K+Mlog2N/K

(9)

即：

log8K

(10)

显然很容易满足上述算式，分割参数K越大，即点云模型分割越细，时间效率越高，当然，前提是需要保证配准的精度。

3 实验与分析

为了验证改进算法的正确性和有效性，分别选择了三组点云模型进行配准实验。实验的系统环境为内存4G、Window7操作系统，软件为MATLAB R2014a。第一组实验选择经典的Bunny模型，图2(a)为配准前的点云(灰色部分表示模型点云，黑色部分表示待配准点云)，待配准点云包含10 753个点。图2(b)为传统ICP算法的配准结果，平均配准时间为26.588 s。图2(c)为本文改进算法的配准结果，平均配准时间为9.695 s。实验结果显示，采用本文提出的改进算法，配准速度大大提高，同时配准精度得到改善。

图2 Bunny模型配准

第2组实验采用David 3D三维扫描仪，从两个不同角度分别对一个瓶子模型进行扫描，对获得的两个点云模型先粗配准再精配准。图3(a)为精配准前的点云，待配准点云包含2 473个点。图3(b)为传统ICP算法的配准结果，平均配准时间为3.785 s。图3(c)为本文改进算法的配准结果，平均配准时间为2.327s。第3组实验采用Cat模型，图4(a)为配准前的点云，待配准点云包含21 530个点。图4(b)为传统ICP算法的配准结果，平均配准时间为87. 829 s。图4(c)为本文改进算法的配准结果，平均配准时间为22.451 s。

对于不同规模的点云模型分别采用传统ICP算法、文献[15]提出的算法和本文提出的改进算法进行配准，不同规模点云情况下3种算法的配准时间如表1所示，3种算法的配准精度如表2所示。可以看出，本文首先利用八叉树对空间进行分割，建立局部KD树索引，有效减少最近邻搜索时间，使算法的效率有较大的提高，并且在精度上有较好的改善。在数据量非常大的海量点云情况下，降低配准时间的效果更加明显。

图3 杯子模型配准

图4 Cat模型配准

表1 不同数据规模点云配准的时间比较

表2 不同数据规模点云配准的精度比较

为了分析分割参数K对最近邻搜索效率和配准精度的影响，分别采用有5 456、10 753、21 530个点的点云模型进行实验。图5所示是不同分割参数K对最近邻搜索时间的影响，图6所示是不同分割参数K对配准精度的影响。显然，K越大，搜索效率越高，但误差越大。对于不同形状的三维点云模型，通过反复实验和分析，可以得到接近于最佳效率和精度的分割参数K值。

图5 分割参数K对搜索时间的影响

图6 分割参数K对配准精度的影响

4 结束语

本文对ICP及其改进的点云配准算法进行深入研究，针对大规模点云数据配准KD树查询效率较低的问题，提出一种基于八叉树与KD树索引的点云精配准方法。利用八叉树空间结构特点和KD树搜索特性，采用欧氏距离阈值、点对距离差值和法向量阈值剔除错误点对。将索引代价较小的全局八叉树与高效的局部KD树相关联，在保证满足配准精度的前提下，能够不同程度提高配准算法的时间效率，特别适合于较大规模的点云模型的配准。通过算法分析和实验可以看到，在实际应用中，设置合适的模型点云分割参数，是影响改进算法配准效率和精度的关键。

[1] XIE J, HSU Y F, FERIS R S, et al. Fine registration of 3D point clouds fusing structural and photometric information using an RGB-D camera[J]. Journal of Visual Communication & Image Representation, 2015, 32:194-204.

[2] BESL P J, MCKAY N D. A method for registration of 3-d shapes[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 1992, 14(2): 239-256.

[3] 刘丰华. 复杂模型三维点云自动配准技术的研究[D]. 天津：天津大学，2013.

[4] LIU Z X，AN J，JING Y. A simple and robust feature point matching algorithm based on restricted spatial order constraints for aerial image registration[J]. IEEE Transactions on Geoscience and Remote Sensing, 2012, 50(2): 514-527.

[5] 刘江，张旭，朱继文. 一种基于K-D树优化的ICP三维点云配准方法[J]. 测绘工程, 2016, 25(6): 15-18.

[6] MULCHRONE K F. Application of delaunay triangulation to the near_est neighbor method of strain analysis[J]. Journal of Structural Geology, 2003, 25(5): 689-702.

[7] 何婧，吴跃，杨帆，等. 基于KD树和R树的多维云数据索引[J]. 计算机应用，2014，34(11):3218-3221.

[8] 杨建思. 一种四叉树与KD树结合的海量机载LiDAR数据组织管理方法[J]. 武汉大学学报(信息科学版)，2014，39(8): 918-922.

[9] WANG Yujian, TAN Shaowei, DONG Weiwei, et al. Research on 3D modeling method based on hybrid octree structure[J]. The Open Electrical & Electronic Engineering Journal, 2014, 8: 323-329.

[10] ARYA S, MOUNT D M, NETANYAHU N S, et al. An optimal algorithm for approximate nearest neighbor searching fixed dimensions[J]. Journal of the ACM, 1998, 45(6): 891-923.

[11] 郑明玲, 许柯, 刘衡竹, 等. 基于重叠区域的高性能近似KD树算法[J]. 计算机辅助设计与图形学学报, 2015, 27(6): 1053-1059.

[12] 韦盛斌, 王少卿, 周常河, 等. 用于三维重建的点云单应性迭代最近点配准算法[J]. 光学学报, 2015, 35(5): 31-37.

[13] LIAW Y C, LEOU M L, WU C M. Fast exact k nearest neighbors search using an orthogonal search tree[J]. Pattern Recognition, 2010, 43(6): 2351-2358.

[14] 祝继华, 尹俊, 邗汶锌, 等. 面向低维点集配准的高效最近邻搜索法[J]. 模式识别与人工智能, 2014, 27(12): 1071-1077.

[15] 钟莹, 张蒙. 基于改进ICP算法的点云自动配准技术[J]. 控制工程, 2014, 21(1):37-40.

[16] 李聪波，肖卫洪，杜彦斌，等. 基于改进ICP算法的损伤零部件精确配准方法[J]. 计算机集成制造系统, 2016, 22(4):1021-1028.

[责任编辑：张德福]

Point cloud registration based on octree and KD-tree index

WANG Yujian, LIAN Tengfei，WU Mingming,GAO Qian

(School of Information, Beijing Union University, Beijing 100101, China)

A multilayer index structure based on octree and KD tree is reported for low query efficiency problems in multi-dimensional queries of KD tree. First, the octree global index for model point cloud is established. Then, the local data KD-tree indexes are built in the octree leave nodes. To improve the traditional Iterative Closest Point algorithm, the local point cloud data is quickly located based on the global index of leave nodes. Using the local KD tree indexes, the searching speed of closest point is sped up. Part of the noise points are removed by the euclidean distance threshold, the difference of interval of point pair and the normals threshold. Experimental result indicates that the proposed method can improve the efficiency and accuracy of registration.

point cloud registration; octree; KD-tree; ICP algorithm

2016-12-20

国家自然科学基金资助项目(61271369)

王育坚(1963-)，男，教授.

廉腾飞(1991-)，女，硕士研究生.

著录：王育坚，廉腾飞，吴明明，等.基于八叉树与KD树索引的点云配准方法[J].测绘工程，2017，26(8)：35-40.

10.19349/j.cnki.issn1006-7949.2017.08.008

TP391

A

1006-7949(2017)08-0035-06