Delaunay三角网生成算法的研究与实现

摘 要 Delaunay三角网作为一种主要的数字地形模型表示法，经过二十多年来的研究，它的生成算法已趋于成熟。本文在简单回顾和评价了分割—归并法、逐点插入法、三角网生长法等三类主流算法的基础上，介绍并实现了一个融以上算法优点于一体，兼顾空间与时间性能的合成算法。
关键字 数字地层模型； 三棱柱； Delaunay；三角网； 生成算法

0 引言
计算机图形学是利用计算机研究图形的表示、生成、处理、显示的学科。经过30多年的发展，科学可视化已成为计算机图形学中最活跃的分支之一，并得到了广泛的应用。在地质领域,由于大量珍贵的地层钻探数据需要用有效的方式进行直观地表达，因而致使可视化技术成为地质研究和工程勘查领域必不可少的手段。
在建模中，2.5维的分析处理由DTM(数字地形模型)模型进行。DTM主要由栅格与TIN（不规则三角网）两种数据格式来表示[1，2]，而以后者更为重要。TIN的生成算法中，最终有三种为普遍接受和采用，它们是分割—归并法、逐点插人法和逐步生长法。本文在简要分析了上述算法所有缺点的基础上，实现了一种合成算法。
1 Delaunay三角网生成算法回顾
Tsaj根据实现过程，把生成Delaunay三角网的各种算法分为三类：分治算法；逐点插入法；三角网生长法。Tsai为比较算法性能，给出了一张各种算法的时间复杂度对照表，如表1所示。
表中，N为数据点数。0(f(N))表示算法的时间复杂度，它以算法中频度最大的语句频度f(N)来度量。
上述三类算法中，三角网生长法在80年代中期以后就很少用到，较常见的是分治算法和逐点插入法，而这两类算法又各有其长处和短处。逐点插入法虽然实现过程相对简单，所需内存较小，但它的时间复杂度高。所以从时间复杂度方面看，分治算法最好。但由于算法中存在递归，它需要较大内存空间。在普通的计算机平台上，运行速度慢和占用较大内存都是应该尽量避免的。本次设计中，我们引入并实现了一种合成算法，将逐点插入法植入到了分治算法中，互相取长补短，从而达到了较好的时空性能，也很好地体现了两者的优势。

表1 几种Delaunay三角网生成算法的时间复杂度[3]
算法 一般情况 最坏情况
分 Lewis和Robinson(1987) O(NlogN) O(N2)
治 Lee和Schachlter(1980) O(NlogN) O(NlogN)
算 Dwyer(1987) O(NloglogN) O(NlogN)
法 Chew(1989) O(NlogN) O(NlogN)
逐 Lawson(1977) O(N4/3) O(N2)
点 Lee和Schachlter(1980) O(N3/2) O(N2)
插 Bowyer(1981) O(N3/2) O(N2)
入 Watson(1981) O(N3/2) O(N2)
法 Sloan（1987） O(N5/4) O(N2)
三 Green和Sibson(1987) O(N3/2) O(N2)
角 Brassel和Reif(1979) O(N3/2) O(N2)
网 MaCullagh和Ross(1980) O(N3/2) O(N2)
生 Mirante和Weigarten(1982)O(N3/2) O(N2)
成
法

注：表中N表示数据量
2 合成算法的研究与实现
2.1 算法基本思想和步骤
把分治算法与逐点插入法结合起来的具体做法是,以分治算法为主体,当递归分割数据集的过程进行到子集中的数据量小于一个预定值——分割阈值时终止,然后用逐点插入法在子集中生成子三角网。我们把这一新的算法称为合成算法。
合成算法的基本步骤是【4】：
Begin
把点集V以横坐标为主，纵坐标为辅按升序排序，然后递归地执行以下步骤：
if V中数据量大于一给定值，把V分为近似相等的两个子集VL和VR
在VL和VR中用合成算法生成三角网；
找出连接VL和VR中两个凸壳的底线和顶线；
由底线至顶线合并VL和VR中两个三角网；
else
生成基本三角网；
End
2.2 算法的实现
2.2.1主要模块
基本三角网的生成模块：
(1)生成凸壳过程。（基于“数据点集”之上进行处理）
(2)生成初始三角网过程。（基于“凸壳”之上进行处理）
(3)点插入过程。（基于“初始三角网”与“点集”之上进行处理）
(4)局部优化过程。（基于“三角网”之上进行处理）
三角网归并模块：
（1）寻找要归并的两个相邻子三角网的上下顶线与底线。（基于凸壳）
（2）两个子三角网归并过程。（基于三角网）
涉及到三个数据对象：初始点集，凸壳，三角网，作为程序中的三个主要对象。
2.2.2数据类型定义
在程序设计中,数据类型定义的好坏对程序设计的影响是很大的,良好合理的数据类型定义,将极大地方便程序中数据的处理,本文中,根据ＴＩＮ在数据文件中的存储方式定义其相应的数据结构。
对于离散数据点,结合点号与单项链表结构定义如下：

typedef struct vertex
{
//顶点
int x;
int y;
int z;
}Vertex;

typedef struct node
{
//链表中的一个节点
Vertex data;
struct node *next;
}Node;
其中 ,ｘ,ｙ为数据点的位置坐标,ｚ为数据点高程值,next为指向下一个点数据的指针,点和点之间以此指针相连形成单向链表。
对于三角形,结合三角形号与单向链表定义三角形Ｃ++数据类型如下：
typedef struct triangle{//三角形
Vertex n;
struct triangle *adj[3];//三个邻接三角形
struct triangle *next;
}Triangle;
其中,ｎ为三角形号,adj为指向三角形3个顶点数据的指针,next为指向下一个三角形数据的指针,它使该数据结构形成单向链表。这种数据结构既便于数据查找,又便于内存分配。
3 结论及展望
本文讨论了合成算法的优势和实现办法，此算法的时间复杂度为O (N m)[4]。在数据量比较大时, 合成算法所耗费的时间与点数几乎成线性的关系。
笔者利用VC++6.0的MFC和OpenGL在微机上实现了对地层模型的三维显示与剖分算法，实验表明，这种算法在实际应用中具有较好的效果。
参考文献
〔1〕李青元. 三维矢量结构GIS拓扑关系研究,博士论文,北京: 中国矿业大学，1996
〔2〕柯正谊，建邦，天河编著。数字地面模型[M]. 北京: 中国科学技术出版社，1993
〔3〕V. J. D. Tsai. Delaunay Triangulations in TIN Creation：an Overview and a Linear-time Algorithm[J]. Int .J. of GIS. 1993,7(6):501-524
〔4〕武晓波,王世新,肖春生. Delaunay三角网的生成算法研究. 测绘学报[J]. 北京:中科院遥感应用研究所
〔5〕B. J. Larkin. An ANSI C program to determine in expected linear time the vertices of the convex hull of a set of planar points [J]. Computers & Geosciences, 1997, 17(3) : 431－443
〔6〕Sedgewick, R. Algorithms [M]. Addison2Wesley, New York, 1998
〔7〕贾志刚. 精通OpenGL[M]. 北京: 电子工业出版社,1998
〔8〕徐青. 地形三维可视化技术[M]. 北京: 测绘出版社,2000
〔9〕李志林,朱庆. 数字地面模型[M]. 武汉: 武汉测绘科技大学出版社,2000