基于三维激光扫描数据的三维实景构建

赫春晓，李 浩

（河海大学 地球科学与工程学院，江苏 南京 210098）

地面三维激光扫描作为一种可以快速获取空间数据的新型技术，越来越多地被应用于城市建模、工程建设等领域。三维激光扫描技术又称“实景复制技术”［1］，不同于传统测量技术的单点定位方式，而是利用激光测量单元进行从左到右、从上到下的全自动高精度步进测量，进而得到完整的、全面的、连续的、关联的全景点的三维坐标［2］。本文以Trimble地面三维激光扫描仪为例，以地物作为采集对象，进行三维构建。

1 基于网格的数据压缩

三维激光扫描获得的点云数据，其数据量比较大，在建立三角网时，会导致三角形过小、过多，导致计算慢、效率不高；而在本文三维场景的构建中，不需要过于细密的三角网，所以要对点云数据进行压缩。

本文采用了基于网格的压缩方法，其基本原理为将原始数据置于一个长方体中，长方体的长、宽、高分别为：L＝Xmax－Xmin，W＝Ymax－Ymin，H＝Zmax－Zmin，然后将该长方体分成l×m×n个小长方体。统计第（i，j，k）号长方体内的点个数，然后对这些点计算特征量，并保留特征量最大点［3］。

在该原理的基础上，本文做了一定改进，采用了自适应边长，该边长由公式l＝得到，然后将点云的各个点依次放入相应小立方体，由式（1）可得某点隶属的立方体：

再根据每个小立方体内的点云点计算得到其均值作为中间值，然后保留到该中间值距离最小的点，删除剩余点。

2 三角网数字地面模型的建立

2.1 Lawson逐点插入算法

Lawson逐点插入算法［4］的基本思想是：在一已存在的三角网中插入一个点，将该点与包含它的三角形的3个顶点相连，形成新的3个三角形，然后用对角线交换法来优化新形成的三角形，从而保证所建三角网为Delaunay三角网。

其基本步骤为：

1）建立包含所有数据点的初始多边形，该多边形可为1个正三角形或2个直角三角形组成的矩形。

2）从数据域中取出任一点P，做如下工作：

找出包括点P的三角形T，设T的3个顶点为V1，V2，V3；

P与T3个顶点相连，形成3个新的三角形T1，T2，T3；

对所有新形成的三角形，用LOP算法进行优化。

3）裁剪数据范围外的多余三角形。

2.2 点在三角形中的判断

算法中要确定插入点所在的三角形，根据三角网拓扑关系和三角形面积坐标可以对该点是否在三角形内作出判断。

设三角形3个顶点为V1（x1，y1），V2（x2，y2），V3（x3，y3），任给点为P（x，y），按有限元理论，P在三角形△V1V2V3内的面积坐标L1，L2，L3定义为

式 中：A为△V1V2V3面 积；A1，A2，A3分 别 为△V2V3P，△V3V1P，△V1V2P的面积。

当L1，L2，L3都大于0时，该点在该三角形内；否则则在该三角形外。

2.3 三角化过程三角形拓扑关系的创建

在当前三角网中插入一点P时，利用2.2方法判断插入点的位置（当插入点落在三角形边上时，三角形的拓扑关系维护比较麻烦。为维护程序的稳定性和通用性，在保证计算精度和数据范围的条件下，本文对共线点位置施加一微小的改正数，即将共线点转化为点在三角形内［5］），确定其所属三角形，同时修改新增加三角形数拓扑组成。由于每插入一点三角形数目增加2个，故拓扑的修改只要对这3个三角形即可，其伪代码为：for i＝0to n

tempTr//存储当前的三角形

trNumID［0］＝curTrNumIndex//当前的三角形号

trNumID［1］＝trNum＋＋//trNum为当前三角形个数

trNumID［2］＝trNum＋＋

newTr1.pt［0］＝tempTr.pt［0］

newTr1.pt［1］＝tempTr.pt［1］

newTr1.pt［2］＝i//插入点的点号

newTr1.tr［0］＝tempTr.tr［0］；

newTr1.tr［1］＝trNumID［1］

newTr1.tr［2］＝trNumID［2］

newTr2.pt［0］＝tempTr.pt［1］

newTr2.pt［1］＝tempTr.pt［2］

newTr2.pt［2］＝i//插入点的点号

newTr2.tr［0］＝tempTr.tr［1］；

newTr2.tr［1］＝trNumID［2］

newTr2.tr［2］＝trNumID［0］

newTr3.pt［0］＝tempTr.pt［2］

newTr3.pt［1］＝tempTr.pt［0］

newTr3.pt［2］＝i//插入点的点号

newTr3.tr［0］＝tempTr.tr［2］；

newTr3.tr［1］＝trNumID［0］

newTr3.tr［2］＝trNumID［1］

End

依照2.1节步骤遍历点云中的所有点构建三角网，完成地物的建模过程。

3 纹理映射

本文在纹理映射过程中，不同于计算机图形学一般常用的cutmall、二步法等纹理映射方法，创新性将摄影测量空间后方交会及共线方程引进来并进行映射。

3.1 三维激光扫描数据的二维球面映射

利用Trimble三维激光扫描仪获取点云的空间坐标（X，Y，Z）和激光强度（R）。由三维激光扫描仪测量原理及其激光点坐标计算公式［6］，本文将三维激光数据映射到球面即α，θ面，其中α为横向扫描角度观测值，取值范围为（0，2π），θ为纵向扫描角度观测值，取值范围为（－π/2，π/2），具体公式为

再将映射到球面的点云以Δα，Δθ等角划分，形成二维矩阵，建立点云各个点对应的矩阵位置坐标，其计算公式为

其中：M为二维矩阵的行数，N为二维矩阵的列数。

从已知的点云数据中，无法直接得到Δα，Δθ的值，本文提出一种自适应方法通过点云空间坐标间接得到Δα，Δθ的值：由点云空间坐标计算得到α，θ的最值，得总面积（αmax－αmin）（θmax－θmin），每个点的面积为Δα×Δθ，从而得到Δα，Δθ的值为

式中n为点云总数。

该方法对于点云数据球面映射后分布越均匀、扫描时水平间距和垂直间距越接近，得到的结果精度越高。

通过上述方法得到三维激光数据的二维矩阵形式，再将该矩阵对应栅格中的点云强度计算该位置的特征量，将该特征量作为最后显示图像的局部像素值，这里特征量以栅格内各点云强度的均值作为最终值。将得到的位图作为内外方位元素的计算时同名控制点的选取所用的图像。

二维激光强度图像矩阵为

3.2 影像内外方位元素的计算

单像空间后方交会中，计算内外方位其中一种方法是共线条件方程解法，即以像点的像空间坐标和相应物点的工程坐标系坐标为已知条件，将共线方程线性化，列出条件方程，经过最小二乘平差以及迭代计算，求出像片内外方位元素［7］。本文将激光数据按照3.1节所述的拓扑结构映射到二维平面以图像形式显示，再同实际相片比较，选择同名点由共线方程经过平差计算得到内外方位元素。

计算内外方位元素时，要保证三维激光坐标需同近景摄影测量坐标系方向大致相同，才能进行计算。本文三维激光扫描仪的坐标系中，X轴、Y轴都在横向扫描面内，Z轴垂直横向扫面［8］；而近景摄影测量的像空间辅助坐标系是以摄影中心S为坐标原点，X轴、Y轴同像平面坐标系X，Y轴平行，Z轴与主光轴重合。而主光轴的方向基本同三维激光坐标的横向扫描面一致，因此三维激光坐标系的Z轴和近景摄影坐标系z轴基本上是垂直的，无法直接解算相机的内外方位元素，需要进行坐标变换。转换公式为

其中：α＝（maxα＋minα）/2－π/2，α为横向扫描角度观测值。

根据求得的内外方位元素，可以计算得到目标点云在相应摄像机图像处的位置坐标，再分别除以摄像机图像原尺寸的宽度和高度，得到相应的纹理坐标，最后利用OpenGL相关函数完成纹理映射，纹理映射步骤如下：

1）将外部图像数据读入内存，用auxDIBImage－Load（）函数可以实现；

2）调用glGenTextures（）函数创建一个纹理，该函数参数为第一步返回的指针；

3）调用glBindTexture（）函数绑定纹理；并进一步通过glTexImage2D（）函数生成纹理，同时还必须指定纹理数据的指针、纹理的大小、纹理类型等；

4）调用glTexParameter（）函数控制纹理的映射方式，包括：纹理滤波、重复与缩放；

5）将纹理坐标作为参数调用glTexCoord（）函数。

4 实验与结论

试验数据由Trimble三维激光扫描仪采集得到，以Visual Stdio6.0作为开发平台，用VC＋＋语言开发完成了整个三维场景的显示。图1为试验结果图。

原始点云总数331 883，压缩后点云为8 152，压缩率97.6%；构建三角网用时0.063s，共15 420个三角形，构网效率较高。

原始点云总数99 239，压缩后点云为4 111，压缩率95.9%；构建三角网用时0.031s，共7 463个三角形，构网效率较高。

从试验结果可以看出，本文中压缩算法压缩率高，比较适用于一般精度要求不高、效率要求较高的情形；三角网构建用时较短、效率较高，构建的三角网符合Delaunay三角网的要求，在实验二构建的三角网中可以看到部分三角网突起，这是由于扫描仪分辨率的设置过大导致该处点云的部分丢失，出现了孤立的点云，而本文没有对其进行粗差剔除；从最终纹理的映射结果可知，本文提出的利用空间后方交会及共线方程求纹理坐标的方法是可行的，能够满足一定的精度要求，同时也可以看到有些部位的纹理放大后不清晰，这是由于该处点较少，导致三角网中的三角形过大，纹理不丰富，放大后会出现模糊，造成纹理拉花现象；点较多的部分纹理丰富，不会出现模糊现象。

［1］范海英，杨伦，邢志辉，等.Cyra三维激光扫描系统工程应用研究［J］.矿山测量，2004（3）：16－19.

［2］毛方儒，王磊.三维激光扫描测量技术［J］.宇航计测技术，2005，25（2）：1－6.

［3］吴杭彬，刘春.三维激光点云数据的空间压缩［J］.遥感信息，2006（2）：22－24.

［4］LAWSON C L.Software for C1Surface Interpolation［M］.Mathematical SoftwareⅢ.J.Rice，Ed.New－York：Academic Press，1977.161－194.

［5］刘学军，王永君，车伟涛，等.三角剖分中拓扑关系的动态创建与维护［J］.湖南科技大学学报：自然科学版，2008，23（2）：79－83.

［6］宋宏.地面三维激光扫描测量技术及其应用分析［J］.测绘技术装备，2008，10（2）：40－43.

［7］李德仁，周月琴，金为铣.摄影测量与遥感概论［M］.北京：测绘出版社，2001.

［8］吴静，靳奉祥，王健.基于三维激光扫描数据的建筑物三维建模［J］.测绘工程，2007，16（5）：57－60.