点云数据三角化

第1章、

1.1空间曲面上散乱数据三角剖分的概念

目前有多种方法可以获得物理模型的形状信息.在制造工业中最常使用的是坐标测量机

(CoordinateMeasuringMachine,CMM)。坐标测量机能精确测量物体表面上点的位置，但其测

量速度较慢，当测量点数较多时，效率很低，一般用在对精度要求较高的场合，如检查零件的形

状精度、位置精度等。当需要大量获取零件表面的数据点时，一般使用激光扫描仪(Laser

Scanner)o激光扫描仪能在相对较短的时间内得到大量零件表面的数据点.另外一种在医学上常

用的测量设备是计算机断层扫描仪(ComputerizedTomography,CT).CT得到的是物体的轮廓线，

数据点呈层状分布，每一层代表物体的一个剖面.

这些测量设备得到的数据点形式各不相同，虽然在局部上某些数据点具有有组织的状态，如

激光扫描仪和CT所得的数据点呈现层状的特点，但在全局上基本均表现出散乱的特点.

所谓散乱数据的三角剖分就是给定一组散乱数据点，将各数据点之间以三角

形相互连接，形成一张三角网格。其实质是以三角网格反映数据点与其邻近点间的拓扑连

接关系。而正确的拓扑连接关系将有效揭示散乱数据集所蕴涵的原始物体表面的形状和拓扑结

构。

1.2空间曲面上散乱数据三角剖分的研究意义及应用范围

空间曲面上散乱数据的三角剖分是构造散乱数据插值曲面的必不可少的前置处理步骤，也

是最重要最关键的一步，基于散乱数据点三角剖分构造数乱数据插值曲面的过程如图1所示：

图1。1构造散乱数据插值曲面的步骤

空间曲面上散乱数据的三角剖分是在对测量数据点必要的处理之后进行的，它是构造散乱

数据插值曲面的前置处理步骤。平面域内的散乱数据的三角剖分研究己经经历了相当长的时间，

相关理论与算法己经相当成熟，特别是Delaunay三角剖分及其优化准则等研究成果使得平面内

的散乱数据点的三角剖分已经不再困难.但在把平面内的算法推广到空间曲面上时，由于空间曲

面散乱数据点之间拓扑关系的复杂性，对其直接剖分的理论和算法尚不完善。在处理空间曲面上

数据点时，一般的算法都是基于平面凸域的或者是已知各种约束条件的情况，对于与多值曲面对

应的空间曲面上散乱数据点的三角剖分算法研究的较少，且在算

法效率和剖分效果上还远远不能让人满意。

散乱数据曲面重构的难点在于如何在数据点集中快速自动得到邻近点间正确的拓扑连接关

系。目前的测量设备能够在短时间内得到数万乃至几十万数据点，所能体现的曲面形状信息越

来越精细和复杂，因此对曲面构造的效率和效果提出了较高的要求。研究散乱数据特别是大规

模散乱数据的三角剖分，对于迅速构建数据点之间的拓扑连接关系，进而构造插值曲面具有十

分重要的意义.

散乱数据的三角剖分不但是构造散乱数据插值曲面的重要基础，对三维数据场的可视化、快

速原型制造等新技术的研究也有巨大的推动作用。因此广泛应用于测量造型、计算机视觉、医

学、气像、勘探、环保等领域中。

本文所要研究三维散乱数据的直接剖分方法旨在探索解决与复杂曲面对应的散乱数据点的

三角化方法，快速准确地生成优化的三角网格，为构造散乱数据插值曲面做好准备。因此本文关

于空间曲面上散乱数据三角剖分方法的研究对散乱数据插值曲面的构造以及逆向工程曲面重构

方法的研究有着重要的现实意义，对三维数据场的可视化、基于CT图像的体数据的三维模型重

建等应用研究也有一定的借鉴与参考价值。

第2章2

2.1引言

空间曲面上散乱数据的三角剖分是科学计算与分析中的一种重要方法，在计算几何散乱数

据插值曲面构造和三维数据场可视化中得到了广泛应用.对空间曲面上散乱数据的三角剖分方

法的研究不论是在二维平面区域还是在三维空间区域上都己经有了很多成果，尤其是对二维平

面散乱数据三角剖分的研究，其理论和算法已经比较成熟，相对来说对空间曲面乱数据的三角剖

分，特别是曲面形状比较复杂、散乱数据点的数目很大时，目前的算法在适应性、执行效率等

方面还有待进一步提高。

2.2问题描述及分类

问题2.1：给定一组散乱数据点｛Vi｝(i=1,2,?,n),如何将各数据点之间以三角形相互连

接，形成一张三角网格，并使网格质量较优。该问题的解是散乱点集｛Vi｝的一个三角剖分T,其实

质是以三角网格M反映各个数据点与其邻近点之间的拓扑连接关系，从而揭示数据点之间的内

在本质联系。三角剖分所涉及的问题在实际应用中根据数据点位置的不同有三种情况：二维，

三维实体和空间曲面.根据这种情况，空间曲面上散乱数据的三角剖分可分为对空间曲面上散乱

数据投影域的剖分和在空间中直接剖分两种类型。空间曲面上散乱数据的投影域包括平面域和

球面域.直接三角剖分方法研究如何直接将空间曲面上散乱数据点在空间中连接成一个优化的

三角网格。

2.3.2三角剖分基本理论

2.3.2.1Voronoi图与Delaunay三角剖分

平面三角剖分的实质是以三角形反映数据点与其邻近点之间的拓扑连接关系。若能首先找

出平面上一点的所有邻近点，则问题2。1在平面上的情况就好办多了，为此需要解决下面的问

题：

问题2：给定平面中n个点的集合S=〈P1,P2,P3,P4…〉对于平面中比其它点更接近于P

的点(x,y)轨迹是什么？

图2。3

定义2。4:给定平面中n个点的集合s,V(pi)在平面S中比其它点更接近于的Pi点的

轨迹是n—1个半平面的交，它是一个不多于联1条边的凸多边形域，称V(pi)为关联于P的

Voronoi多边形或关联于P的Voronoi域.图2。4表示关联于pi的一个Voronoi多边形，它是

一个五边形，n=6o

图2.4

对于S中的每一个点都可以作一个Voronoi多边形，这样的n个Voronoi多边形组成的图成

为Voronoi图，记为Vor(S),如图2.5所示。该图中的顶点分别称为Voronoi顶点和Voronoi

边。Vor(S)的边是某点对的垂直平分线上的一条线段或者半直线，从而为该点对所在的两个多

边形域所共有。Vor(S)中有的多边形域是无界的。

在叙述Voronoi图的性质时作假设：原来集合S中没有四个点是共圆的。若此假设不成立,

则需要在证明中加入一些无关紧要的，但却是冗长的陈述。

在叙述Voronoi图的性质时作假设“原来集合S中没有四个点是共圆的。若此假设不成立,

则需要在证明中加入一些无关紧要的，但却是冗长的陈述。

图2。6

图2。6

定理2。1：对于S的Voronoi图每一个顶点v,圆C(v)不包含S的其它点。由于定理2。1：

Vor(S)又称为最近点意义下的Voronoi图。

定理2。2:在S中，Pi的每一个最近点确定Voronoi多边形v(pi)的一条边.

图2.7

到2。7Pl的每一个最近邻近点确定v(pi)

定义2.5:Voronoi图的直线对偶是由S的每个点对之间加入一个直线段以获得嵌入平面的

图，产生的图是原来n个点上的一个图。如图2.5中的虚线所示。对偶的重要性主要归于下面

的DeIsunay定理。

定理2.3:Voronoi图的直浅对偶是S的一个三角剖分。图2。5中的虚线即是S的一个三

角剖分。

Voronoi图的这些性质可以用来快速构造Voronoi图，且可应用它解决最邻近点问题和三角

剖分问题。

定义2.6：在一个三角剖分白，若每一个三角形的外接圆为空(即在外接圆中不含有其它点)，

则该三角剖分称为DeIaunay三南剖分。

对于给定的一组平面数据点集S,可以多种不同的三角剖分，其中Delaunay三角剖分是最

优的，二维Delaunay三角剖分由满足最小内角最大准则的三角形组成。下一节将介绍详细三角

剖分准则.

232.2三角剖分优化准则

在三角剖分过程中，往往用一种比较简单的方法构造散乱数据点的初始三角网格，然后对其

进行优化，以获得较优化的三角形网格。优化的方法和效果取决于所采用的优化准则。平面三

角剖分常采用的优化准则有Thiessen区域准则、最小内角最大准则、圆准则。

Sibson证明了这三个准则的等价性，并指出符合这三个准则的三角剖分只有一个，即

DeIaunay三角化。今Thiessen区域准则(Thiessenregioneriterion)Thiessen区域指Voronoi

分割后形成的多边形区域。如果两个区域具有非零长度的公共线段，则称这两个区域的生成点

为Thiessen强邻接点(strongThiessenneighbours);如果它们的公共部分仅是一个点,则称

为Thiessen弱邻接点(weakThiessenneighbours).一个严格凸的四边形至多有一对相对顶点

是Thiessen强邻接点。Thiessen区域准则指对一个严格凸的四边形三角化时，将Thiessen强

邻接点相连，若两对顶点都是Thiessen弱邻接点，则任选一对相连，这样构造的三角化是最优的,

见图2.8:

图2.8Thiessen区域,准则

最小内角最大准则

在平面内，对一个严格凸的四边形进行三角化时，有两种选择，最小内角最大准则就是要

保证对角线两侧两个三角形中的最小内角最大。如图2。9所示：

图2.9最小内角最大准则

圆准则

严格凸的四边形中的三个点确定一个圆，如果第四个顶点落在圆内，则将第四个顶点与其相

对的顶点相连，，否则将另两个相对顶点相连，如图210所示：

图2.10圆准则

平面三角剖分优化准则理论为平面内散乱点集的快速三角剖分提供了判断标准，也为三维

空间散乱点集的三角剖分优化标准提供了借鉴依据。

2.3.3经典三角化算法

三角化算法虽然很多，但大多算法生成的是Delaunay三角网格，即以空外接圆(球)准则为

优化准则。根据实现方法的不同，Delaunay三角化方法。

主要有三类：

换边法(Swappingedges):首先构造非优化的初始三角化，然后对2个共边三•角形形成的

凸四边形迭代换边优化。以Lawson为代表的对角线交换算法属于换边法，换边法适用于二维

Delaunay三角化，对于三维Delaunay三角化，则需要对共面四面体进行换面优化。

加点法(Addingpoints):从一个三角形开始，每次加一个点，保证每一步得到的当前三角

化是局部优化的.以英国Bath大学数学分校Bowyer,Green,Sibson为代表的计算Dirichlet

图的方法属于加点法，是较早成名的算法之一；以澳大利亚悉尼大学地学系Watson为代表的空

外接球法亦属于加点法。加点法算法简明，是目前应用最多的算法。分割占有法(Devideand

conquer):将数据域递归细分为子块，每块实现局部优化三能化，然后合并。对应于上述三类算

法各有一些著名的算法.

2.3.3.1Bowyer算法

该算法基于Dirichlet图的构造，适用于n维空间中的离散点集的网格剖分，是英国Bath

大学教学分校的A.Bowyer在总结该校RobinSibson教授、PeterJ.Green博士七十年代所做工

作的基础上，于1981年提出的。

算法思路

Bowyer算法是针对Voronoi图的实现，首先开始于一个由n+1个数据点形成的DeIaunay单

纯体，这样将得到一个包含一个真实顶点、而其它顶点为无穷远点0的Voronoi图：如图277

所示，往己有数据形成的Voronoi图中加入新点Q,从包含新点的Voronoi多边形开始，利用

Voronoi多边形相邻关系，找到最近点，构造新加入点的Voronoi多边形(粗虚线).

图2-17在已有Vororwi图中如入新点Q生成其Voronoi多边形的过程

算法分析

该算法整体效率为0(n(k+1/k)),k为空间的维数。基于Voronoi图的Bowyer算法计算

复杂，空间消耗大，算法时空效率不高。

233.2Watson算法

该算法是澳大利亚悉尼大学GeoIogy与Geophyics系D.FWatson于1981年提出的，是用

计算机构造晶体模型的研究成果。

算法思想

首先给出由一个或多个外接球不包含其它数据点的单纯体组成的初始网格，然后往其中加

入一个数据点，考察外接球的包含情况，去除包含新点的n维单纯体，用n+2个点能组合的、外

接圆不包含其它点的单纯体取代。如图278所示是Watson算法的具体思路。具体实现时，可

一次性全部找出并删除所有包含新加入点的单纯体，得到一个包含新点的空洞，空洞的边界面与

新点相连，得到新点加入后的Delaunay网格，这样可避免进行新生成单元的外接圆是否包含old

points的计算，具体加入一点的算法流程叙述如下。

算法流程

(1)加入新点，搜索单纯体链表，查找外接球包含新点的单纯体，所有包含新点的单纯体组

成一个多面体。

(2)包含新点的单纯体的各个面加入一个临时链表，若一个面加入两次，说明是两个单纯

体共享的面，必然位于多面体的内部，从链表中删除，若出现新点位于外接球上的退化情形，则抛

弃链表和新点，改用其它方法处理。

(3)若未出现退化情形，则将新点与多面体的各个面相连，得到新的单纯体。新点加入过

程结束。

E2-13Wauon算法流修””

算法分析

Watson算法的思路非常简明，易于编程实现。但当出现k+2(k为空间维数)点位于同一圆

球面时，则三角化结果不唯一，这种情形称为退化情形,如图279,为二维空间的退化情形.

除了如图2-20所示的规则数据，实际应用中散乱数据点集很少出现退化情形。但由于计算机的

计算精度是有限的，当新点与外接球球面之间的距离小于预设计算精度，则认为新点位于球面

上；这种新点与球面距离的计算误差可能引起拓扑关系的不一致。

2-19四点共Bl的退化情形S2-20烷则数塞点生成的网格

2.333四叉树、八叉树方法

四叉树、八叉树本身是数据结构，当用于空间编码时，可进行实体造型，适当修改即成为网

格剖分算法。MAYerry.MSShephard于1983年、1984年发表13了四叉树、八叉树在二维、

三维网格剖分中的应用，有些文献中将他们的算法称为Shephard-Yerry算法.Shephard-Yerry

算法只适用于域剖分。网格单元可以是四边形/六面体，也可以是三角形/四面体。这种算法的特

点是：与实体造型相结合，自动化程度高，网格密度可调整，剖分速度快，内部单元形状比较

好，但边界单元形状很难保证。在二维空间，以矩形网格为例，边界单元共有16种被裁情形；

三维六面体网格单元被边界面截切后共有4096种情形，被切后的单元与原有网格单元拓扑关系

可能不一致，需要进行拓扑变换.此外，这种方法不具备几何不变性，即剖分对像旋转后，剖分

结果发生变化.

韩国中央大学YHJung和KLee于1993年提出了一种直接基于四面体的八叉树空间编码法。

这种算法一步到位，内部单元为四面体，边界单元易于处理为四面体，不需要拓扑变换。总之，

Shephard—Yerry算法由于与实体造型技术相结合，是一种很有前途的方法，一旦边界单元得到

很好处理，则必在实体网格剖分方面占据主要地位。

2.3.3.4换边》奂面

1977年，美国加利福利亚工学院喷气推进实验室的CharlesL。Lawson提出了基于边交换

的二维DeIaunay三角化。力口拿大埃德蒙顿Alberta大学计算机科学系BarryJoe分别于1989

年、1991年给出了基于局部换面的三维Delaunay三角化的算法和证明。

算法分析

换边换面法适用于离散点集剖分和域剖分，算法过程简明，容易实现.但是三维三角化算法

过程中需要检测的面很多，有效控制变换的范围，合理的数据结构和快速的查询法是提高速度

的关键。换边换面法的最大的困难在于如何处理不可变换的情形，这是本文研究的关键问题。只

有解决不可变换的情形，才能得到delaunay三角剖分网格；否则，结果是非delaunay的。

2.3.3.5网格的前百生成法。

网格的前沿生成法从提出到现在发展最快，现在有了很多的算法。最早法国学者SHL。于

1985年在文献中提出了网格前沿法的雏形，只将网格前沿法作为节点连元的方法，没有与节点生

成同时考虑。英国学者J.Peraire,M。Vabdati,K.Morga.等于1987年实现了按方向细化的网

格前沿法。此后研究的各种网格前沿法最大的特征是能够生成复杂形状的非结构网格，其按方

向细化的特点特别适合于三维可压缩流的优化笄法。

算法思路

以剖分域的边界为网格的初始前沿，按-设定的网格单元的形状，尺度等要求向域内生成节

点，连成单元，同时更新网格前沿，如此逐层向剖分域内推进，直到所有的空间都被剖分。

算法分析

网格前沿法能够处理比较复杂的对像，其主要特点是提供了控制网格密度和质量的手段。

但是网格前沿法中存在大量的查询操作以及网格前沿面的相交检测，这些操作是很浪费时间的。

很多算法在这两点上作了快速的处理。下图是网格前沿法生成的三角网格实体。

第3章空间曲面上散乱数据点的快速三角剖分算法实现

3.1采用的数据结构

点类，点节点类，点链表类，边类，边节点类，边链表类（多边形类），三角形类，三角形节

点类，三角形链表类。

以下对这几种结构做简单的介绍（只介绍主要的）：

点类中的数据成员

cIassPoint_T

（

doubIem_X;//空间点的坐标

doubIemY:

doubIem_Z;

〃下面是一些方法

）

点节点类中的数据成员.

cIassPNode_T

（

Pont_Tm_Point;

PNode_T

PNode_T*m_Right;

〃下面是一些方法

）

点链表类中的数据成员。

cIassPNList_T

（

PNodeT*m_Head;

PNode_T*m_TaiI；

//下面是一些方法略

）

边类中的数据成员。

cIassEdge_T

（

Pont_Tm_Point1;

Pont_Tm_Point2;

PontTmPoint3;〃所在三角舫第三个点

booIm_Used;〃表示这条边是活边还是死边

//下面是一些方法略

}

边节点类中的数据成员.

cIassENode_T

（

Edge_Tm_Edge;

ENode_T*m_Right;

ENode_T*m_Left;

〃下面是一^些方法略

）

多边形类中的数据成员。

cIassPoIygone_T

（

ENode_T*m_Head;

ENode_T*m_TaiI;

Intm_LivingEdgeNum;//边界前沿中的火边数量。

//下面是一些方法略

）

三角形类中的数据成员。

cIassTriangle_T

（

Pont_Tm_Point1;

PontTmPoint2;

Pont_Tm_Point3;//所在三角影第三个点

//下面是一些方法略

)

三角形节点类中的数据成员

cIassTriNode_T

(

Triangle_Tm_Triangle;

TriNode_T*m_Right;

TriNode_T*m_Left;

〃下面是一^些方法略

)

三角形链表类中的数据成员

cIassTriNList_T

(

TriNode_T*m_Head;

TriNode_T*m_TaiI;

//下面是一些方法略

}

点链表：存储点云数据的一个单向链表。

三角形链表：生成的三角形存储在单链表中.

前沿边界(多边形)：多边形采用边的双循环链表的形式存储。

3.2三角剖分的过程

定义:

内边：三角剖分结果中被两个三角形公用的边。

外边：三角剖分的结果中只被一个三角形拥有的边C

内点：如果一个点的相连边都是内边这个点就是内点。

活边：没有经历找点生成新边过程的边.这里要强调，活边经历找点生成新边过的程，可能

找到了匹配点，可能没找到匹配点，可能会产生0条新边，1条新边或2条新边.

死边：经历一次找点生成新边过程的边就是死边。死边可以是边界边，也可以是外边。

外边框：由外边首尾相连构成的空间多边形就是外边框。

边界点：在外边框中的顶点.

外点：没有被选择过的点.

3.2.1算法简单描述：

(1):点云的预处理。

(2):构造一个相对饱满的初始三角形作为种子三角形。

(3):开始循环：如果外边界框中的活边的数量不是零，就继续下面的步骤。否则算法结束。

(4)：对外边框做循环：ENMov=：外边框的头节点.

(5):如果ENMov是活边：选择匹配的点，如果选择到了匹配点就更新点链表，更新三角形

链表，更新外边框链表；如果没选择到匹配点，把ENMov标志成死边used=1。如果ENMov是死

边：什么也不做。

(6):ENMov=：外边框更新前ENMov的下一个节点.如果ENMov是空的，那么返回(3),否则

返回(5)。以下对•算法中的各个步骤做详细的解释.

3.2.2数据点的预处理

(1):把点云数据文件中的点读到点单链表中.

(2):PNMov=：点链表的头节点.

(3):在链表中删掉PNMov节点后面到PNMov距离小干DIST的所有的节点.其中DIST是一

个可以指定的数，其数值越大剩下的点越少。

(4)：PNMov=：点链表中PNMov的下一个节点。

(5):PNMov如果不是链表的尾节点：返回(3)。否则结束.处理后的点相对少得多，也保证

处理后的点云要相对的均匀，最主要的是去掉了曲率过大的点。

3.2.3初始三角形的建立

(1):先得到点云链表中的第一个点作为firstPoint.

(2)：然后得到在firstPoint后面，距离firstPoint最近点作为第二个点secondPointo

(3):这两个点构成了一条边.在secondPoint后面的点中寻找距离第一条边的两个端点的

距离和最小的点作为第三个点thirdPointo

(4):如果这个三角形的最小的内角小于30度。那么选择点云链表firstPoint节点的下一

个节点为firstPoint,返回(2).否则，结束。

3.2.4活边寻找最佳匹配点的算法

为了加快运算的速度我们采用包围盒算法。一个点要成为当前边的候选点要具备的必要条

件。

i:这个点和当前边构成的三角形于当前边所在的三角形钩成的面角要大于给定的阀值

MINFACEANGLEo

下面介绍空间中两个面的二面角的计算方法。

P2

p4

p3

三角形p1p2p3和三角形p1p2p4所成的二面角的大小等于向量p5p4与向量p6p3所成的角。

P1p5是p1p4在p1p2上面的投影向量，p1p6是向量p1p3在p1p2上的投影向量.

Ii:候选点和当前边构成的三角形的最小内角不能小于MININANGLE.

Iii:候选点必须是包围盒中的点，包围盒是这样定义