不规则体体积计算三维激光点云切片法

李 斌，魏俊博，马博超，王 璐,徐明霞

1. 长安大学地质工程与测绘学院，陕西 西安 710054； 2. 西部矿产资源与地质工程教育部重点实验室，陕西 西安 710054； 3. 陕西国防工业职业技术学院，陕西 西安 710300

体积是空间体对象形态分析的重要参数[1]，各种形体的体积计算，既是一个久远又新颖的命题[2-4],也是空间体对象形态分析的基本内容[5-6]。体积计算涉及规则几何体和不规则几何体两种类型。其中，规则几何体的体积计算有现成的公式甚至算法可用，操作简便、精度高且应用广泛[7]。相对而言，不规则几何体因形态各异，体积计算并无普适、统一的方法，甚至无章可循[8-13]，因此，仍是普遍面临和亟待解决的现实难题[14-17]。三维激光扫描技术的兴起和运用，客观上为不规则体的体积计算提供了新模式[18-19]。例如，基于三维激光扫描点云的不规则树冠体积[10]、商品包装的点云体积[20]、船舶排水量的点云体积[21]等的估算。由此可见，只要能用三维激光扫描方式完成空间体对象扫描并获取其整体外型轮廓的点云数据，不管其形态是否规则，都有可能求解其体积。另一方面，若三维激光扫描的空间体对象是规则的，则可用来检验所提计算方法的正确性和有效性。

探索激光扫描技术方法计算不规则体的点云体积，在很多领域已有应用且成效显著,但精度有待精化、适用性尚需提升[22-24]。目前，较为直接且简单的点云体积计算方法是“点云切段法[25]”(后简称“切段法”)。该方法的基本思想和流程可概括为“空间体对象激光扫描→点云数据切段→分段点云投影→点云片轮廓边界确定→点云片面积计算→点云体积计算”等6个步骤，而其中的“点云片轮廓边界确定”和“点云片面积计算”是该法求解三维激光扫描空间体对象点云体积的两个关键环节。由于存在可能的点云片外轮廓多边形“形态畸变”、固有的计算系统性“放大效应”和“三角形累加法”面积算法适用性不足等问题，造成该法不仅计算环节过多，而且计算结果偏差大，甚至可能错误等不良后果。因此，有必要提出适用三维激光扫描空间不规则体对象的点云体积计算得更为简洁、有效的方法。

本文从简化“切段法”计算环节、避免出现“放大效应”、防止产生点云片轮廓多边形“形态畸变”和选用更有效的面积算法4个方向入手，提出“点云切片法”(后简称“切片法”)。该法采用双向最近点搜索法取代可能引起平面点云外轮廓多边形“形态畸变”的扫描法，从算法层面消除了切片外轮廓多边形点云散点搜索出错的可能性；通过切片而非切段的技术路线，规避了“切段法”的投影环节，从而避免了由投影带来的系统性的三重“放大效应”，并有效减少了中间步骤；选用行列式模型法计算点云切片面积，从计算层面解决面积算法适用性的问题；用可靠的柱体体积计算方法计算点云分段体积，最终求和得到三维激光扫描空间体对象的整体体积，又确保了实体对象体积计算的精准。因此，从简化、改进、完善流程和算法等核心环节入手的“切片法”有望解决三维激光扫描前提下的不规则体的点云体积计算问题。

1 点云切片法

以水平方向切割点云为例，建立在算法与编程计算基础上的“切片法”的基本思想和流程可概括为“空间体对象激光扫描→点云数据切片→轮廓边界确定→切片面积计算→点云体积计算”5个步骤，具体为：

(1) 空间体对象激光扫描。用地面三维激光扫描仪扫描空间体对象，获得其点云数据D(如图1(a)所示)。

(2) 点云数据切片。在点云纵向最小值0与最大值H之间，用一组(设为n+1个)等间距(间距为h，如式(1)所示)的水平面自上而下顺序切割点云，依次得到系列水平点云切片Si(如图1(b)、式(2)所示)。

(3) 轮廓边界确定。使用双向最近点搜索法[26]取代扫描法对乱序的各平面点云Si进行排序，生成各点云切片散点外轮廓边界多边形Pi(i=0,1,…,m)。

(4) 切片面积计算。分别计算Pi围成的面积，即Si的面积Ai(i=0,1,…,n)(如式(4)所示)。

(5) 点云体积计算。累加Ai并乘以h得到点云体，也就是三维激光扫描物体的体积V(如式(3)所示)

(1)

i=0,1,…,n}

(2)

(3)

式中，n是切片数减1；x、y、z是点云点坐标。

在“切片法”中，切片平面点云外轮廓边界多边形正确确定与切片面积正确计算是“切片法”计算点云体积并确保正确的两个关键环节。

1.1 轮廓边界确定

点云切片外轮廓边界多边形的正确生成是“切片法”切片面积能否正确计算的基础。相对而言，“切段法”中采用的射线360°扫描法虽然思路简单、计算复杂度低，且在特定的平面点云外轮廓多边形中可得到极佳的形态效果(如图2(a)所示)，但当平面点云轮廓为极端凹多边形时，“切段法”就难以规避“形态畸变”的致命缺陷，即当平面点云外接矩形的重心落在平面点云外面时，搜索得到的平面点云外轮廓边界不能保证不走样、不变形(如图2(b)所示)，因此，这就有可能引起平面点云面积计算错误，进而导致体积计算畸变。

由此可见，为避免“切段法”点云片的“形态畸变”出现，“切片法”的双向最近点搜索法就成为正确生成点云切片外轮廓边界多边形的关键举措,其核心思想如下：

(1) 选择边界起点。在切片点中任选一点(如y坐标值最小)PS作为起点，以其最近点PE为终点生成多边线。

(2) 确定最近点。取多边线一端PS(PE)，在剩余点中找到距PS(PE)点最近的点P。

(3) 生成边界边。分别计算PS和PE到P的距离dS和dE，并比较dS和dE的大小，若dS

(4) 搜索并判断进程。判断切片上的点是否搜索完成。若是，结束搜索，生成点云切片外轮廓多边形边界；否则转到步骤(2)。

当切片这种点云片外接矩形重心落在切片外时，取代射线360°扫描法的双向最近点搜索法就能确保对平面点云散点数据进行正确的排序(如图2(c)所示)。

1.2 切片面积计算

作为激光扫描对象点云体积计算的早期算法，“放大效应”是“切段法”的又一痼疾。在“切段法”中，因投影而生的点云环带的外轮廓边界就成为相应点云段块投影后的最大外延，这样，有赖“点云分段投影”的“切段法”就无形之中放大了平面点云外轮廓的边界，从而造成“放大轮廓边界”的第一重放大效应；由于面状物体是以其轮廓边界构成的多边形来表示的[1]，因此，平面点云轮廓边界多边形放大的实质，就是投影面面积的放大，由此导致“切片面积放大”的第二重放大效应；而事实上，点云切段的柱体体积由切片面积Ai和段距h两个因素决定，Ai放大而h一定，则必然引起“点云体积放大”的第三重放大效应。由此可见，连锁出现的三重“放大效应”，就使激光点云对象的体积呈现单边放大的系统性后果。

若规定m个顶点p0，p1，…，pm=p0，按逆时针回路首尾相接构成多边形并设为Pi，则Pi的面积Ai可由行列式(如式(4)所示)计算得到[1]

(4)

式中，xj、yj为切片平面点云外轮廓多边形Pi(i=0,1，…,n)的顶点pj(j=0,1,…,m)的坐标；j为点云切片外轮廓边界多边形的顶点编号；i为点云切片的编号；n为点云切片个数减1。

“切片法”采用的双向最近点搜索法，聚焦点云切片排序问题，从算法层面有效规避了外轮廓边界多边形生成时出错的可能性，为点云切片面积正确计算创造了有利条件；而在计算点云切片面积时，“切段法”采用的“三角形累加法”是以投影切片外接矩形的重心为起点，依次与投影面上的两个相邻点构造三角形，并将三角形面积累加得到投影切片的面积；同样，由于外接矩形重心位置的不确定，该法也不能适用于极端凹多边形时的面积计算。为解决面积计算方法适用性的问题，本文采用了上述行列式的面积计算模型用于点云切片面积的计算，不仅克服了“三角形累加法”适应性不强、局限性大的缺点，且计算流程清晰、易于编程、结果可靠，从计算层面保证了点云体积的计算建立在面积计算准确无忧的基础之上。

1.3 要点与优势

对不同物体而言，点云切片的外轮廓边界形态虽然各异，但因“双向最近点搜索”算法的普适性，绕开了与中心位置关系密切的扫描法的排序陷阱，从而既确保了外轮廓边界多边形的计算与点云切片的凹或凸的形态无关，也保证了点云切片面积计算的正确性，从而为后续体积计算结果的准确无误奠定了基础。事实上，“双向最近点搜索”算法也可用在传统的“切段法”中，以解决由“射线360°扫描”算法引起的“形态畸变”的可能，使“切段法”有所改进——不至于因点云片外轮廓多边形的计算错误而导致后续面积计算乃至体积计算出错的严重后果。

2 点云体积计算

本文选用规则的圆锥体模型(如图3(a)—3(c)所示)和不规则的墩台柱体实物(如图4(a)—4(c)所示)两种对象，使用三维激光扫描仪获取研究对象的点云数据，并对其进行拼接、去噪[27-29]等处理后(其流程详见图5)，通过圆锥模型计算、圆锥墩台点云计算及其结果对比，验证“切片法”点云体积计算算法的可行性、正确性、高效性和适用性。其中，点云体积算法可行性、正确性的判断，源于对圆锥体的计算与分析。圆锥体是规则几何形体，其模型尺寸可以量算得到(如表1、表2及图6所示)，体积可直接根据模型尺寸用圆锥体体积公式(如式5所示)计算获得(详见表3)，并可作为似真值使用

(5)

式中，r是圆锥模型底半径；h是圆锥模型高度；V是圆锥模型体积。

图1 点云数据及其切片示意图Fig.1 Diagram of the point cloud data and its slicing

图2 平面点云外轮廓多边形形态示意图Fig.2 The shape of the contour polygon of the slice point cloud

图3 圆锥模型、尺寸及点云数据Fig.3 Conical model, size and its point cloud data

图4 石狮实物图及其整体与墩台点云数据Fig.4 Physical map and point cloud data of stone lion, and the point cloud data of the stone lion base

因为模型尺寸可以精确量算，因此，由模型尺寸和式(5)计算所得的圆锥体模型体积值真实、可靠，可作为圆锥的体积真值看待和使用。而通过与圆锥体体积近似值，即点云方式体积计算值的比较，则既能对形态规则的圆锥体点云体积计算方法可行性进行判断，也能对体积计算算法正确性进行验证。墩台是不规则的几何形体，其精确几何尺寸未知，但当墩台切割间距为最小分辨率时，用验证过并证明为正确的“切片法”算法计算求得的墩台点云体积值就可以作为近似真值看待；当然，用“切段法”也能计算墩台点云的体积，并可与“切片法”进行比较以进一步判断和验证两种点云体积计算算法的效率优劣和适用情况。

图5 三维激光扫描对象点云数据处理流程Fig.5 Flow charts of 3D laser scanning point cloud data processing

表1 圆锥模型直径d两种方法的量测结果

表2 圆锥体模型高度h的量算结果

图6 圆锥模型尺寸与体积Fig.6 Geometric size measurement of conical model

表3 圆锥模型尺寸与体积

2.1 圆锥体体积计算

2.1.1 圆锥体点云数据体积计算

以正轴、横轴和斜轴(如图7(a)—7(c)所示)3种方式检验点云“切片法”和“切段法”的正确性和效能。首先，在点云纵向最小值0与最大值H之间，进行对象点云切片，依次得到系列水平点云切片Si；其次，使用双向最近点搜索法确定点云外轮廓边界多边形Pi；再次，分别计算Pi围成的面积，即点云切片Si的面积Ai；最后，累加Ai并乘以h得到三维激光扫描的圆锥的体积。

图7 圆锥体三轴体积计算示意图 Fig.7 Diagram of conical volume calculation in three directions

2.1.2 圆锥体体积计算结果分析

分别用“切片法”与“切段法”计算正轴、横轴与斜轴3种工况下的圆锥体点云的体积(详见表4—6)。

表4 正轴圆锥体切片(切段)计算结果对比分析

由表4数据可知：①两种方法的两种误差与切距大小均成正比；②计算用时切片法与切距成反比，切段法与切距成正比，而在最小切距时，两者用时相当。

表5 横轴圆锥体切片(切段)计算结果对比分析

由表5数据可知：①两种方法都保持间距小、精度高的基本趋势；②虽然当切距最小时，切片法计算用时最长，切段法计算用时最短，但任何情况下切片法效率都优于切段法；③对规则圆锥模型而言，误差较正轴时均大幅减少，说明横轴是更好的切割姿态。

表6 斜轴(α=50°)圆锥体切片(切段)计算结果对比分析

由表6数据可知：①两种方法都保持了间距小、精度高的特点；②间距越小，切片法计算用时越长，切段法计算用时反而越短，任何情况下切片法效率都高于切段法；③误差大小与横轴时大致相当，说明斜轴(50°)是可以接受的切割姿态。

对比分析可知：因圆锥模型体积似真值(即公式计算值，如表3所示)与点云体积计算值(如表4—6所示)足够接近，再结合算法分析与精度效率试验结果(详见表4—6,图8)，可以认定“切片法”算法可行、结果正确且计算高效。

2.2 墩台柱体体积计算

2.2.1 墩台柱体点云数据体积计算

用激光扫描仪获取墩台柱体三维点云坐标数据；分别用“切片法”与“切段法”计算正轴、横轴与斜轴(如图9(a)—9(c)所示)3种工况下的墩台柱体点云体积。

2.2.2 墩台柱体体积计算结果分析

以正轴、横轴和斜轴(见表7—9，图10)3种方式检验点云“切片法”和“切段法”的计算效率与适用性(表中墩台皆以“切片法”正轴1 mm切距时的点云体积计算结果473 142 270 mm3为似真值进行计算和比对分析)。

由表7数据可知：①两种方法的两种误差与切距大小都成正比；②计算用时切片法与切距成反比，切段法与切距成正比，而在最小切距时，用时两者相当。

由表8数据可知：①间距变化与切片(切段)数变化相反；②间距变化与绝对误差和相对误差变化紊乱；③间距与切片法用时成反比，与切段法用时成正比。

图8 两种方法圆锥体积计算精度与效率走势图Fig.8 Calculation accuracy and efficiency charts of cone volume in two methods

图9 墩台三轴体积计算示意图Fig.9 Diagram of stone lion base volume calculation in three directions

表7 正轴墩台切片(切段)计算结果对比分析

表8 横轴墩台切片(切段)计算结果对比分析

表9 斜轴(α=50°)墩台切片(切段)计算结果对比分析

由表9中数据可见：①间距变化与绝对误差和相对误差变化一致；②间距变化与切片法用时成反比，与切段法用时成正比；③误差与轴向的关系：误差均匀，各轴向均适宜。

由表7—9和图10可知：“切片法”在精度和效率上均优于“切段法”，且随着切距的增大，优势持续扩大。

3 讨 论

3.1 试验综述

圆锥、墩台两种对象，3种工况下的点云体积算例计算结果(如表4—表9、图8、图10所示)表明：

(1) “切片法”的计算结果更接近真值，算法正确性和可靠性不仅与前述理论分析一致，而且得到了试验数据的有效验证；

(2) 从时间上看，相同切距条件下，“切片法”的计算用时更短，说明计算流程简约的“切片法”算法效率更高、更有效；

(3) 点云切距大小可控，又与精度有关的事实说明，点云体积计算的精度可由切距大小控制，即切距越小精度越高，切距越大精度越低。

因此，计算正确、流程简洁、结果可靠、算法高效、精度可控的“切片法”，不仅全面优于“切段法”，而且适用性更强，能够解决不规则体体积计算的难题，当然前提是该不规则体对象可用三维激光扫描仪进行整体扫描，且在点云数据处理时应顾及并选择多个方向切割、计算并相互印证为好。

3.2 误差分析

“切片法”在切片间距之间出现点云形态极端变化时，体积计算可能因补偿性(即某种不确定的随机性差异互补性)或其他不确定性，而造成切距小、误差反倒大的反常情况出现(如表5“切片法”中1 mm切距的误差反倒比2 mm的误差大，表6中也有类似反常现象，且更为明显)，因而误差小即精度高不足以说明准确度高，但在没有补偿或补偿等不确定性因素影响不明显情况下间距小准确度高则基本是可以肯定的，即无论切段间距还是切片间距越小，体积计算精度越高，结果也越可靠。

对规则圆锥体而言，正轴时“放大效应”这种系统性偏差会更大(1 mm间距都能差到5%以上)，产生此情形的原因是圆锥体出现了圆柱化现象，且可能在圆锥体正轴时达到极点，这说明切割方向很重要。

另外，墩台体点云体积横轴切片或切段计算可能出现误差的原因是由于沿横轴切割会出现一个切片(切段)上有两个或以上多边形切片(或段块)单元即所谓“多环[30]”的情况，从而导致前述切片(或段块)面积偏离实际甚至错误的结果(如图1(b)中S7、图10(c)及图11所示)。

4 结 论

本文提出的“切片法”用双向最近点搜索的方式，消除了平面点云“形态畸变”的隐患；用切片计算体积的技术路线，根除了三重“放大效应”的弊端；用行列式模型法计算点云切片面积的办法，解决了面积算法适用性不足的问题；用柱体体积计算方式计算点云体积的方式，确保了对象体积计算的正确。因此，从流程算法入手的“切片法”，不仅使高效、精确计算不规则体的体积成为可能，而且能够解决特定条件下不规则体的点云体积计算问题。

理论分析和算例结果一致表明，点云体积计算是不规则体体积计算的有效方式，“切片法”易于编程实现、推广应用前景广阔，整体优于“切段法”。但所遗留的多环问题尚有待继续探索，可以预见的基本思路是：首先，进行聚类分析[1]以区分多环；其次，待截面多环个数及各环切片边界依次正确确定后，逐一进行截面各环切片面积的计算；最后，依据切割顺序逐环理清各环切片纵向对应关系，再设法解决三维激光扫描对象点云体积先分部、后整体的计算问题。此过程较为复杂，既需要前述工作成果作为基础，也需要破解许多新问题，这些都有待跟进并可望解决。

点云体积计算方法从“切段法”演进到“切片法”，标志着基于三维激光扫描的不规则体对象体积计算的方式方法从有到好、走向成熟。