Delaunay三角剖分及matlab实例

Delaunay三⾓剖分及matlab实例

鉴于Delaunay三⾓剖分在点云拟合⽅⾯有不错的应⽤，现对该算法的原理进⾏简单的汇总~

----------------------------原理部分------------------------

1、三⾓剖分与Delaunay剖分的定义

如何把⼀个离散⼏何剖分成不均匀的三⾓形⽹格，这就是离散点的三⾓剖分问题，散点集的三⾓剖分，对数值分析以及图形学来说，都是极为重要的⼀项处理技术。该问题图⽰如下：

1.1 三⾓剖分定义

【定义】三⾓剖分：假设V是⼆维实数域上的有限点集，边e是由点集中的点作为端点构成的封闭线段，E为e的集合。那么该点集V的⼀个三⾓剖分T=(V,E)是⼀个平⾯图G,该平⾯图满⾜条件：1、除了端点，平⾯图中的边不包含点集中的任何点。2、没有相交边。//边和边没有交叉点

3、平⾯图中所有的⾯都是三⾓⾯，且所有三⾓⾯的合集是散点集V的凸包。

//：⽤不严谨的话来讲，给定⼆维平⾯上的点集，凸包就是将最外层的点连接起来构成的凸多边型，它能包含点集中所有的点。

1.2 Delaunay三⾓剖分的定义

在实际中运⽤的最多的三⾓剖分是Delaunay三⾓剖分，它是⼀种特殊的三⾓剖分。先从Delaunay边说起：【定义】Delaunay边：假设E中的⼀条边e(两个端点为a,b)，e若满⾜下列条件，则称之为Delaunay边：存在⼀个圆经过a,b亮点，圆内（注意是园内，圆上最多三点共圆）不含点集V中任何其他的点，这⼀特性⼜称空圆特性。【定义】Delaunay三⾓剖分：如果点集V的⼀个三⾓剖分T只包含Delaunay边，那么该三⾓剖分称为Delaunay三⾓剖分。

1.3 Delaunay三⾓剖分的准则

要满⾜Delaunay三⾓剖分的定义，必须符合两个重要的准则：1、空圆特性：Delaunay三⾓⽹是唯⼀的（任意四点不能共圆），在Delaunay三⾓形⽹中任⼀三⾓形的外接圆范围内不会有其它点存在。如下图所⽰：

2、最⼤化最⼩⾓特性：在散点集可能形成的三⾓剖分中，Delaunay三⾓剖分所形成的三⾓形的最⼩⾓最⼤。从这个意义上讲，Delaunay三⾓⽹是“最接近于规则化的”三⾓⽹。具体的说是在两个相邻的三⾓形构成凸四边形的对⾓线，在相互交换后，两个内⾓的最⼩⾓不再增⼤。如下图所⽰：

1.4 Delaunay三⾓剖分的特性

以下是Delaunay剖分所具备的优异特性：1、最接近：以最接近的 三点形成三⾓形，且各线段（三⾓⾏的边）皆不相交。2、唯⼀性：不论从区域何处开始构建，最终都将得到⼀致的结果。3、最优性：任意两个相邻三⾓形构成的凸四边形的对⾓线如果可以互换，那么两个三⾓形六个内⾓中最⼩⾓度不会变化。4、最规则：如果将三⾓⽹中的每个三⾓形的最⼩⾓进⾏升序排列，则Delaunay三⾓⽹的排列得到的数值最⼤。5、区域性：新增、删除、移动某⼀个顶点只会影响邻近的三⾓形。6、具有凸边形的外壳：三⾓⽹最外层的边界形成⼀个凸多边形的外壳。

1.5局部最优化处理

理论上为了构造Delaunay三⾓⽹，Lawson提出的局部优化过程LOP(Local Optimization Procedure),⼀般三⾓⽹经过LOP处理，即可确保为Delaunay三⾓⽹，其基本做法如下所⽰：1、将两个具有共同边的三⾓形合成⼀个多边形。2、以最⼤空圆准则作检查，看其第四个顶点是否在三⾓形的外接圆内。3、如果在，修正对⾓线即将对⾓线对调，即完成局部优化过程的处理。LOP处理过程如下图所⽰：

1.6三维空间情形

Delaunay三⾓⽹格的性质可以推⼴到⾼维。⼀个三维点集的三⾓剖分是由四⾯体构成。下图显⽰了⼀个由两个四⾯体构成的简单三维Delaunay三⾓剖分。⼀个四⾯体的外接球来说明空接球准则。

2.Delaunay剖分的算法

Delaunay剖分是⼀种三⾓剖分的标准，实现它有多种算法。

2.1.Lawson算法

逐点插⼊的Lawson算法是Lawson在1977年提出的，该算法思路简单，易于编程实现。基本原理为：⾸先建⽴⼀个⼤的三⾓形或多边形，把所有数据点包围起来，向其中插⼊⼀点，该点与包含它的三⾓形三个顶点相连，形成三个新的三⾓形，然后逐个对它们进⾏空外接圆检测，同时⽤Lawson设计的局部优化过程LOP进⾏优化，即通过交换对⾓线的⽅法来保证所形成的三⾓⽹为Delaunay三⾓⽹。上述基于散点的构⽹算法理论严密、唯⼀性好，⽹格满⾜空圆特性，较为理想。由其逐点插⼊的构⽹过程可知，遇到⾮Delaunay边时，通过删除调整，可以构造形成新的Delaunay边。在完成构⽹后，增加新点时，⽆需对所有的点进⾏重新构⽹，只需对新点的影响三⾓形范围进⾏局部联⽹，且局部联⽹的⽅法简单易⾏。同样，点的删除、移动也可快速动态地进⾏。但在实际应⽤当中，这种构⽹算法当点集较⼤时构⽹速度也较慢，如果点集范围是⾮凸区域或者存在内环，则会产⽣⾮法三⾓形。

2.2.Bowyer-Watson算法 (推荐)

⽬前采⽤逐点插⼊⽅式⽣成的Delaunay三⾓⽹的算法主要基于Bowyer-Watson算法，Bowyer-Watson算法的主要步骤如下：1、构造⼀个超级三⾓形，包含所有散点，放⼊三⾓形链表。2、将点集中的散点依次插⼊，在三⾓形链表中找出其外接圆包含插⼊点的三⾓形（称为该点的影响三⾓形），删除影响三⾓形的公共边，将插⼊点同影响三⾓形的全部顶点连接起来，从⽽完成⼀个点在Delaunay三⾓形链表中的插⼊。3、根据优化准则对局部新形成的三⾓形进⾏优化。将形成的三⾓形放⼊Delaunay三⾓形链表。4、循环执⾏上述第2步，直到所有散点插⼊完毕。这⼀算法的关键的第2步图⽰如下：

2.3伪代码表述

https://github.com/ironwallaby/delaunay对算法进⾏了优化，算法伪代码：

1.

input: 顶点列表(vertices) //vertices为外部⽣成的随机或乱序顶点列表

2. output:已确定的三⾓形列表(triangles)

3.

初始化顶点列表4.

创建索引列表(indices = new Array(vertices.length)) //indices数组中的值为0,1,2,3,......,vertices.length-15.

基于vertices中的顶点x坐标对indices进⾏sort //sort后的indices值顺序为顶点坐标x从⼩到⼤排序（也可对y坐标，本例中针对x坐标）6.

确定超级三⾓形7.

将超级三⾓形保存⾄未确定三⾓形列表（temp triangles）8.

将超级三⾓形push到triangles列表9.

遍历基于indices顺序的vertices中每⼀个点 //基于indices后，则顶点则是由x从⼩到⼤出现10.

初始化边缓存数组（edge buffer）11.

遍历temp triangles中的每⼀个三⾓形12.

计算该三⾓形的圆⼼和半径13.

如果该点在外接圆的右侧14.

则该三⾓形为Delaunay三⾓形，保存到triangles15.

并在temp⾥去除掉16.

跳过17.

如果该点在外接圆外（即也不是外接圆右侧）18.

则该三⾓形为不确定 //后⾯会在问题中讨论19.

跳过20.

如果该点在外接圆内21.

则该三⾓形不为Delaunay三⾓形22.

将三边保存⾄edge buffer23.

在temp中去除掉该三⾓形24.

对edge buffer进⾏去重25.

将edge buffer中的边与当前的点进⾏组合成若⼲三⾓形并保存⾄temp triangles中26.

将triangles与temp triangles进⾏合并27.

除去与超级三⾓形有关的三⾓形28.

end

⽤图来解释伪代码过程，先⽤三点来做实例：

如图，随机的三个点

根据离散点的最⼤分布来求得随机⼀个超级三⾓形（超级三⾓形意味着该三⾓形包含了点集中所有的点），⽅法是根据相似三⾓形定理求得与矩形⼀半的⼩矩形的对⾓三⾓形，扩⼤⼀倍后则

扩⼤后的直⾓三⾓形斜边经过点(Xmax,Ymin)，但是为了将所有的点包含在超级三⾓形内，在右下⾓对该三⾓形的顶点进⾏了横和⾼的扩展，并要保证这个扩展三⾓形底⼤于⾼，才能实现包

含。这样求得的超级三⾓形不会特别⼤使得计算复杂，⽽且过程也简单，并将超级三⾓形放⼊temp triangles中。

未确定三⾓形存⼊temp triangles，Delaunay三⾓形保存到triangles，边存⼊edge buffer。

接下来就像是伪代码中描述的那样，对temp triangle中的的三⾓形遍历画外接圆，这时先对左边的第⼀个点进⾏判断，其在圆内，所以该三⾓形不为Delaunay三⾓形，将其三边保存⾄edgebuffer中，temp triangle中删除该三⾓形。（第⼀个判断的三⾓形是超级三⾓形）

将该点与edge buffer中的每⼀个边相连（红⾊线），组成三个三⾓形，加⼊到temp triangles中，

再将重复对temp triangles的遍历并画外接圆，这时使⽤的是第⼆个点来进⾏判断

1. 该点在三⾓形1外接圆右侧，则表⽰左侧三⾓形为Delaunay三⾓形，将该三⾓形保存⾄triangles中

2. 该点在三⾓形2外接圆外侧，为不确定三⾓形，所以跳过（后⾯会讲到为什么要跳过该三⾓形），但并不在temp triangles中删除

3. 该点在三⾓形3外接圆内侧，则这时向清空后的edge buffer加⼊该三⾓形的三条边，并⽤该点与edge buffer中的三⾓边进⾏组合，组合成了三个三⾓形并加⼊到temp triangles中

再次对temp triangles进⾏遍历，这⾥该数组⾥则含有四个三⾓形，⼀个是上次检查跳过的含有第⼀个点的三⾓形和新根据第⼆个点⽣成的三个三⾓形

1. 该点在三⾓形1外接圆右侧，则该三⾓形为Delaunay三⾓形，保存⾄triangles中，并在temp triangles中删除

2. 该点在三⾓形2外接圆外侧，跳过

3. 该点在三⾓形3外接圆内侧，将该三边保存⾄temp buffer中，并在temp triangles中删除

4. 该点在三⾓形4外接圆内侧，将该三边保存⾄temp buffer中，并在temp triangles中删除

这时，temp buffer 中有六条边，triangles中有两个三⾓形，temp triangles中有1个三⾓形，对temp buffer中的六条边进⾏去重，得到五条边，将该点与这五条边组合成五个三⾓形并加⼊到temp triagnles 中，这时temp triangles中有6个三⾓形

由于三个点已经遍历结束，到了不会再对第三个点形成的三⾓形做外接圆，这时则将triangles与temp trianlges合并，合并后的数组表⽰包含已经确定的Delaunay三⾓形和剩下的三⾓形。这时除去合并后数组中的和超级三⾓形三个点有关的所有三⾓形，即进⾏数组坐标的限定，则得到了最后的结果：

这是⽤最少的三个点来做讲解，点数越多的话计算量会越⼤，但是都是在上⾯步骤下进⾏的。