【Aegisub】不规则三角网(TIN)、德劳内(delaunay)三角网的生成算法

三角网当然就是一堆点连起来形成的三角形网络
在生成不规则三角网TIN(即Triangulated Irregular Network)的时候,非常常见的就是delaunay三角网了
那么什么样的三角网是德劳内三角网呢?就是如果有一堆离散的点,比如其中有个点A,你找到这堆离散点里和A最近的两个点,连接起来,就构成了一个delaunay三角形,所以这一堆点都和与它最近的两个点相连,连完以后,就成了一张狄洛尼(德劳内)三角网了
这样形成的德劳内三角网有一些特性,比如“空圆特性”:
即 对于delaunay三角网中任意两个有公共边的三角形,它们中任意一个三角形的外接圆中都不会包含有另一个三角形的顶点。如图1.01

在图1.01中,ABD和BCD就是两个德劳内三角形,因为C点不会落在ABD的外接圆内,同时A点也没有落在BCD的外接圆内,这个就满足“空圆特性”
知道了delaunay三角网是啥以后,接下来就是用什么算法实现了。实际在各种算法中,都会利用到这个“空圆特性”,比如你求出某个三角形的外接圆,然后看这个外接圆里面有没有其它的点在里面,如果有,那这个三角形就不是delaunay三角形,因为刚刚说了,德劳内三角形的外接圆里面不会有其它的点落在里面。
本文接下来要介绍的算法是插入法(有兴趣的也可以自己去试试别的方法,比如分治法)
插入法的大概步骤是:
1、随机生成一些离散点(当然注意随机生成的点不能重合、两点不能坐标都一样了)
2、将点集的点按照x坐标升序排列
3、计算出一个能包含住点集中所有点的大大的三角形(所有的点都在这个特大三角形内)
4、将排好序的点逐一插入、并利用“空圆特性”连线
OK,一个个的解决啊,首先是生成一堆离散点:
先说最简单的,在一个矩形范围内生成一堆点。那假设你用math.random的话,随机产生的数就可能会一样,如果直接随机,那么产生了这些点以后就要去除重复、把坐标相同重复多余的点去掉。有没有什么办法可以让产生点集的时候不产生重复的点呢,想想啊,不重复意味着坐标不同,就是说x坐标或者y坐标但凡有一个不一样就可以了,所以只要在随机的时候控制比如x的坐标,让x坐标不一样即可,比如图1.02

图1.02,假设矩形的左上角坐标是x0,y0、矩形宽w高h,现在想产生n个点。建立一个表coor用来装随机点。因为要设定每个点的x不同,所以让每一次循环时x的大小都绝对大于上一次。一共要产生n个点,所以把矩形的横宽w分为n份(n个区间),每一份宽就是w/n,所以每一次循环时x增加最大不过w/n即可。那么这里当然得用到math.random了,math.random一般生成整数很友好,但是现在不是要生成整数的时候啊,毕竟w/n多半会是小数,所以此时用math.random()来产生小数,math.random在不填任何参数时,会产生0到1的随机小数(不含端点、不会随机出0或1)
那么区间长w/n,w/n*math.random()就是x在这个区间范围内随机了。好,x坐标保证每次循环结果不一样以后,y现在就可以随便了,y直接在y0到y0+h随机无所谓
然后,还可以顺带计算一个特大三角形(超级三角形),因为反正都要算超级三角形嘛。刚刚提到,超级三角形只要能让所有点都在该三角形内部即可,所以怎么样的三角形“都可以”,只要够大就可以了,那还不容易吗?比如图1.03这一堆点

你可以用如图1.04这样的三角形包住它们,也可以用如图1.05这样的超级三角形包住他们


总之确保特大三角形包含点集中的所有点即可。那么我这里写的代码就是:首先这些点集的范围是知道的,因为这些点都在我们设定的矩形框内,所以超级三角形只要框得住这矩形即可,所以让超级三角形左下的点x坐标在x0减掉10倍矩形宽w、y坐标在y0+h的基础上加上10倍的矩形高h即可,然后超级三角形的其它点的坐标同理计算一下就可以了,反正三角形够大就可以了,你也可以弄比如图1.06这样的等边三角形

图1.06就是所有点都在一个矩形框内,根据这矩形算出一个等边三角形够大。
甚至你说你生成的点集所在范围就不大(反正视频屏幕就那么大),你直接设定特大三角形的左下x是-9999999,y是9999999这种只要保证能包含住所有离散点就都可以
然后因为现在产生的随机点刚好是x升序的,所以就不用排序了(如果要排序的点集,后面会讲怎么排序)
那么此时就到了“最后一步”,逐点插入了。
先来具体解释一下逐点插入是怎么个方法,然后再慢慢地讲算法
假设现在一共有3个离散点如图1.07,那么连成delaunay三角网就得到了图1.08对吧


问题就是这是怎么个连法的呢,用什么方法连接这些离散点的呢:
(1).计算出一个能包含所有点的超级三角形ABC(图1.09)

当然此时ABC中还没有开始插入点的,只是为了方便看到特大三角形包含住了所有点,所以图1.09的三角形ABC里才有3个点,实际上还没有开始插呢
(2).插入第一个点P1(注意之前说了点集中点的顺序是需要按x升序的),如图1.10

(3).检查P1是否在ABC的外接圆内,如图1.11

如果在,刚刚说了,则ABC不是德劳内三角形,则不需要三角形ABC,所以要删除ABC这个三角形,然后重新连线,怎么连,因为你发现P1在△ABC外接圆内,所以将P1点和A、B、C连、形成三个三角形,如图1.12

如图1.12,因为P1在△ABC外接圆内部,所以现在不要△ABC了,而是重新连出了3个△,分别是△AP1B、△P1BC和△AP1C,注意,现在已经没有△ABC了!!△ABC已经删除了
(4).插入第二个点P2,如图1.13

(5).现在有的三角形是△AP1B、△P1BC和△AP1C,遍历每个三角形并做其外接圆,然后检查P2是否在其外接圆内部,如图1.14

这里呢,P2在△AP1B外接圆的右侧,那么直接判定△AP1B是德劳内三角形,所以我们要把它保存起来,那么当然是存在一个表里,比如表起名叫u6h2吧,而且当然这个表用来储存我们最后要获得到的delaunay三角形,所以这个表在你写代码的时候肯定要一开始就建立好,你最后函数返回的就是这个装有delaunay三角形的表。记得一开始就要建立这个表哟,用来放我们找到的delaunay三角形的!
OK,因为这里说了离散点都是按x升序排好的,所以当遍历(即插入)下一个点时,若该点在外接圆的右侧,则表示以后所有的点都在该外接圆的右侧,则保证了Delaunay三角形的空圆特性,所以说,P2在△AP1B外接圆的右侧,就直接判定△AP1B是合法的德劳内三角形,所以添加进一开始就建好的名叫u6h2的表里面。然后注意,既然已经把△AP1B添加进u6h2表里了,那么现在就要删除掉△AP1B了,就像刚刚删除掉△ABC一样,只是现在不需要重新连线了
然后再看,发现P2在△AP1C的外接圆外边,此时这种情况就不做任何改变、不删除、不连线、也不把△AP1C添加进u6h2里。就是说当你发现这个点它既不在这三角形外接圆内也不确定是在外接圆右侧的,那么就判定不做任何改变!然后继续检查下一个三角形△P1BC
现在发现P2在△P1BC的外接圆内部,所以说,△P1BC一定不是德劳内三角形,所以删除△P1BC,并重新连线,就像刚刚一样,将P2和P1、B、C连出三个三角形,那么现在的三角形就是如图1.15这样的了

现在就只有△AP1C、△P1BP2、△BP2C和△P1CP2了
捋一下啊,之前的△AP1B因为是德劳内三角形,所以人家已经去u6h2了,这里就删除掉它了,而△AP1C保持不变,然后原本的△P1BC被分成了3个三角形,△P1BC也被删除了,所以当点P2插入以后就经历了这些变化,三角形就变成只有△AP1C、△P1BP2、△BP2C和△P1CP2这几个了
(6).插入第三个离散点P3,如图1.16

(7).当然继续遍历现在所有的三角形,检查P3是否在其外接圆的内部,如图1.17

P3不在△P1BP2的外接圆内部、并且在△P1BP2的外接圆右侧,所以△P1BP2是合法的德劳内三角形,所以将△P1BP2放进表u6h2里,并删除掉这个三角形。P3在△BP2C的外接圆外面,但是不在内部,所以和刚刚一样,就不做任何的改变、不做任何改动。继续看,P3在△AP1C的外接圆内部、同时也在△P1CP2的外接圆内部,所以△AP1C和△P1CP2不是合法的德劳内三角形,要删除掉△AP1C和△P1CP2,并重新连线,也就是P3和A、P1、P2、C连成4个三角形
△AP1P3、△AP3C、△P1P2P3、△P2P3C,然后这样的话,现在剩下的三角形就有: △BP2C、△AP1P3、△AP3C、△P1P2P3和△P2P3C,如图1.18

(8).因为现在已经插入完所有的离散点了,所以在最后就别忘了,把剩下的这几个三角形放进u6h2这个装德劳内三角形的表里
这样一来,u6h2里面就有原来的△AP1B、△P1BP2和最后加进去的△BP2C、△AP1P3、△AP3C、△P1P2P3、△P2P3C,这些三角形都是合法的德劳内三角形,因为满足“空圆特性”、因为不满足的都删除掉且重新连线了。
那么当然最后要返回的是一开始的离散点连接成的三角网,所以,还要再把u6h2这个表里和超级三角形有关的三角形都删除掉
现在u6h2里有△AP1B、△P1BP2、△BP2C、△AP1P3、△AP3C、△P1P2P3、△P2P3C,其中△AP1B的两个顶点和超级三角形的顶点一样、△P1BP2中的点B是超级三角形ABC的顶点、△BP2C有和超级三角形有关系的顶点B和C、△AP1P3有和特大三角形有关的点A、△AP3C有点A和点C、△P2P3C有和特大三角形有关系的点C,所以△AP1B、△P1BP2、△BP2C、△AP1P3、△AP3C、△P2P3C这些三角形都要从u6h2这个表里删除掉,于是乎,表u6h2里就剩下△P1P2P3了,这个就是我们要的离散点连接成的德劳内三角网,也就是图1.08
以上就是插入法的大概方法
那么知道了怎么用插入法来连出delaunay三角网以后,具体的算法要怎么规划呢?
1、一个三角形用什么表示:一个三角形就是3个点,它们的顺序怎么换,都能连接成三角形,如ABC、BCA、CAB等等等等都是同一个三角形,所以直接把3个点的坐标装进一个table,这个table就表示了一个三角形,如{{0,0},{23,66},{-5,2.5871}}
2、如何判断插入的点是否在某三角形外接圆内部或右侧或外部:根据三角形的3点坐标求得外接圆的中心坐标,以及外接圆的半径,那么插入的一点P(x,y)就能判断在不在圆内了。假如说点P(x,y)到你求出的外接圆圆心的距离比你外接圆的半径小,那么当然就说明P在圆内;假设点P的横坐标x比你求出的外接圆圆心的横坐标大,并且该数值比外接圆半径还大的话,那当然说明P在外接圆的右侧,很简单,如图1.19是个圆

那么一个点如果在这个圆右侧,是怎样的,当然就是图1.20这样的了

所以如果点P的横坐标比外接圆圆心的横坐标大,且大了半径r还多,即判断P在圆右侧。若外接圆圆心横坐标为cx,而P的横坐标为x的话,那当cx+r<x的时候,就说明P在右侧
然后当然,如果P到圆心的距离大于半径了,就说明P在外部,这个就能判断出P是否在外接圆外部了。所以总共3种需要判断的位置关系的方法就有了(①P在外接圆内部 ②P在右侧 ③P在外接圆外)
3、还需要准备哪些东西?
①需要有临时存储三角形的表,就像前文中会有三角形的删除和重连操作,那这些三角形就需要临时的表来装。
②需要有临时储存三角形的边的表,用来连线,因为连线的时候你发现,如果你插入的点它同时在两个三角形的外接圆里的话,结果就是重新连成4个三角形(而不是6个),因为如果你插入的点只在一个三角形的外接圆内部的话,那直接重新连成3个三角形即可,但是如果插入的点同时在两个三角形的外接圆内的话,直接每个三角形重新连出3个新的就不行、因为那样会连出总共6个三角形,所以这里需要去掉重复的,而且是完全去掉,所以只留4个新连出的三角形
那么这样就不如建立一个临时的三角形边的表,用来帮助连线操作。(当然一条边就是两个点,所以形式是{{2,9},{-33,11}}这种的,这就表示一条边,当然现在需要的临时储存边的表要装很多条边,所以就是这种一条条边都装在那临时的边的表里)比如一个△ABC在插入了P以后判断出要重新连,那么就在边表里加入AB、BC、AC三条边,然后插入的点P再和这3条边连成3个三角形即可。所以这样,比如点P同时在三角形ABC的外接圆和三角形BCD的外接圆里的话,就先在边表里装入AB、BC、AC、BC、BD、DC这6条边,然后装好以后,再把重复的完全去掉,现在重复的是BC,那么边表里要完全去掉BC边,即得到AB、AC、BD、DC这4条边,这样再让这四条边和P连成4个三角形即可。(注意反复强调要完全去除重复,你看现在压根就没有BC边了)
③需要记录删除哪些三角形,因为之前说了,有删除三角形的操作,那么得知道要删除的三角形是装临时三角形的table中的第几个元素,要知道这个索引, 所以需要一个临时记录索引的表
好吧,再带着刚刚准备的工具,重新看一遍逐点插入法:
假设现在有4个离散点,如图1.21

假设离散点装在名叫r4m8的表里
(1).建立表u6h2用来存放delaunay三角形,这是最后要返回的表
建立表v9e5用来存放临时的三角形,这是中间要操作的,要删除和重连新三角形的表
现在一开始先在表v9e5里加入超级三角形,如图1.22

然后
因为是逐点插入,所以for循环遍历点集,即for i=1,#r4m8 do
每次循环时都重新建立一个临时储存边的表,和临时储存要删除的临时三角形的编号(索引)的表,即local edges={} local idx={} idx就是index的缩写,表示用来储存相应索引的(即下标)
然后
因为对于每一个插入的点,都需要检查它和每一个临时三角形的外接圆的关系,所以要遍历临时三角形的表,一个个判断其外接圆和插入点的关系,所以再for i1=1,#v9e5 do
现在判断临时三角形中的第一个三角形,△PQR的外接圆是否包含插入的点A,包含,则记录要删除的三角形△PQR的索引1、由于要重新连线,所以边表edges里加入三边
然后
内层for循环结束,即for i1=1,#v9e5 do结束,此时做完了判断工作和准备工作,如图1.23

然后
利用收集到的idx,讲临时三角形中需要删除的三角形删除掉。因为要将idx里记录的每一个三角形删除掉,所以for i2=#idx,1,-1 do,用上remove函数来删除,即table.remove(v9e5,idx[i2])
这样就删除掉了临时三角形表v9e5里所有需要删除的不合法三角形
然后
对边表edges去除重复(刚刚讲了,要完全去除重复),去重以后,用此时插入的点和边表中的边连成一个个新的三角形并把重新连的三角形添加入临时三角形表v9e5中,如图1.24

如图1.24,边表和idx表只是每次辅助你连个线、删除不合法三角形的,现在已经把临时三角形表里原来的△PQR删除了,并利用edges的边连接了新的三角形,然后加进了临时三角形表里了,所以edges和idx在本次循环(即外层的for i=1,#r4m8 do)里,已经出色的完成了任务
然后
插入第二点,即for i=1,#r4m8 do已经进行到i=2了,每次外层循环重新建立边表、idx表,如图1.25

对于插入的第二个点(即外层循环for i=1,#r4m8 do的i=2时),遍历临时三角形表v9e5中的每一个三角形,检查其外接圆和插入点的关系,如图1.26

显然B在△APQ的外接圆外面,不做任何改变,B在△APR的外接圆外部,不做任何改变,而B在△AQR的外接圆内部,所以判定△AQR需要被删除、记录下它的下标、放进idx,并在edges里添加相应的边,如图1.27

然后
此时已经遍历完临时三角形表v9e5里的所有三角形。所以还是利用获得的idx来删除需要删除的不合法的临时三角形。然后将edges里的边完全去重,再利用edges边表进行连线、连成新的三个三角形,并把新连的三角形放进临时三角形表中,如图1.28

然后
插入第三个点,即外层for循环到第3个点了,重新建立边表、idx表,如图1.29

对于插入的点,遍历临时三角形表中的每个三角形,检查其外接圆和插入点的关系
所以,C在△APQ的外接圆右侧,故△APQ是合法的德劳内三角形,直接将其添加入德劳内三角形的表u6h2里面,即u6h2[#u6h2+1]=v9e5[i1],并在idx里记录下该三角形△APQ的下标,因为这个三角形已经判定是德劳内△了,所以等会要从临时三角形表里删除
C在△APR的外接圆内部,在idx里记录下不合法三角形的编号(索引),并在边表edges里添加相应的边,如图1.30

C在△ABQ的外接圆右侧,故△ABQ为德劳内三角形,将其直接添加进德劳内三角形表u6h2里,并在idx里记录下该三角形△ABQ的下标,因为这个三角形已经判定是德劳内△了,所以等会要从临时三角形表里删除
C在△ABR的外接圆内部,在idx中记录下这个需要删除的三角形的下标,在边表edges里添加进△ABR的三条边,如图1.31

C在△BQR的外接圆外部,所以当然不用做任何操作与改变
然后
遍历完了临时三角形表的所有三角形,利用得到的idx对临时三角形表里的三角形进行相应删除,此时删掉第1、2、3、4个三角形。对边表edges进行完全去重,如图1.32

对edges去重以后,再用它和插入点连线。现在就是C和AP连成△ACP、C和PR连成△CPR、C和AB连成△ABC、C和BR连成△CBR,连好以后,当然还是一样添加进临时三角形表里面,如图1.33

然后
插入第4个点,即外层for循环到第4个点了,重新建立边表、idx表,如图1.34

又来又来,对于插入的点,遍历临时三角形表中的每个三角形,检查其外接圆和插入点的关系
插入点D在△BQR的外接圆外部,所以不进行任何改动
D在△ACP的外接圆外部,所以同样不做任何改动
继续,D在△CPR的外接圆外部,又是不做任何改动
D在△ABC的外接圆内部!所以△ABC一定不是德劳内△,在idx里记录下要删除的该三角形的编号(索引),并在edges里添加相应的边(即AB、BC、AC)
D在△BCR的外接圆内部,所以在idx里记录下要删除的该三角形的编号(索引),并在edges里添加相应的边(即BC、BR、CR)
如图1.35

然后利用idx收集到的索引,对临时三角形表中的元素进行删除,然后将边表edges进行完全去重,如图1.36

然后再用去重以后的边表edges里的边来和D连线。D和AB连成△ABD、D和AC连成△ACD、D和BR连成△BCR、D和CR连成△CDR,并将连好的几个三角形放进临时三角形表里,如图1.37

然后
现在已经插入完所有的离散点了,也就是外层for循环都已经结束了(即for i=1,#r4m8 do出色的完成了它的任务!辛苦了! 啊不对,是我写这么一大堆废话才辛苦,啊,这是什么,噢不, 我竟然吐血了)
那么现在就像前文中讲过的那样,在把所有点插入完以后,最后还要将临时三角形表里剩下的三角形都加入进德劳内三角形表里
来看看这两个表现在本来有哪些三角形,首先临时三角形表里的三角形就如图1.37所示,而德劳内三角形表u6h2里有两个三角形,如图1.38,当然这在前面也说得很清楚了

所以,最后把临时三角形表中的三角形添加进德劳内三角形表中,就得到图1.39

那么当然,我们需要的是离散点构成的德劳内三角网,所以这个三角网里不能有超级三角形的顶点参与其中,所以还是像刚刚前面说的,将德劳内三角形表u6h2中的和特大三角形的顶点有关系的三角形全部删除掉,便能得到只用我们的离散点连成的delaunay三角网,所以最后德劳内三角形表中就只剩下如图1.40这两个三角形了,最后返回你得到的这个表u6h2即可

在这么详细的讲了大概的算法以后,简单地提一下一些没讲的细节:
1、随机的点集怎么按x坐标升序排序:用table.sort
2、如何将edges里的边和插入点连成△的:遍历edges的每条边,每条边就是两个点一起的,比如{{-23,66},{55,7}},那么假设它要和点{0,0}连成三角形,就连成{{-23,66},{55,7},{0,0}}即可
3、如何将一个三角形拆成3条边(因为edges的边都是一个个的三角形的3条边加入进edges表里的):就比如三角形{{-23,66},{55,7},{0,0}}你只要拆成{{-23,66},{55,7}}、{{55,7},{0,0}}、{{0,0},{-23,66}}即可
4、如何删除德劳内三角形表中和特大三角形有关的三角形:只需要遍历德劳内三角形表的每个三角形,检查这三角形的顶点有没有和特大三角形的顶点一样的。比如三角形{{-23,66},{55,7},{0,0}},然后假设超级三角形的顶点是{55,7}那么就检查三角形{{-23,66},{55,7},{0,0}}里面有没有{55,7}即可,有的话,就要删除掉三角形{{-23,66},{55,7},{0,0}}
然后,我写的具体的代码就在视频里介绍了。