点击蓝字关注我们
编者按:本文作者奇舞团前端开发工程师魏川凯。
狄洛尼三角剖分(Delaunay triangulation)是指离散分布的点集P的其中一种三角剖分DT(P),可以使得点集P中没有任意一个点严格处于任意一个三角形的外接圆的内部。狄洛尼三角剖分构建的三角网可以尽量避免狭长三角形的出现,有效提高逼近精度,使得网格整体质量保持最优,在三维显示,图像处理,人工智能等领域被广泛运用。
上图为利用三角部分构成三角网拟合我国某山脉地形。
定义
空圆性:DT(P)是唯一的(任意四点不能共圆),在DT(P)中,任意三角形的外接圆范围内不会有其它点存在。
最大化最小角:在点集所有可能的三角剖分中,狄洛尼三角剖分所形成的三角形的最小角最大,也可以说在两个相邻的三角形构成凸四边形的对角线,在相互交换后,两个内角的最小角不再增大。
性质
最外部的三角形边的集合是点集的凸包。
以最接近的三个点形成三角形,且各线段(三角形的边)皆不相交。
不论从区域何处开始构建,最终都将得到一致的结果(点集中任意四点不能共圆)。
任意两个相邻三角形构成的凸四边形的对角线如果可以互换的话,那么两个三角形六个内角中最小角度不会变化。
在点集所有可能的三角剖分中,狄洛尼三角剖分所形成的三角形的最小角最大,至少狄洛尼三角剖分中的最小角与其他三角剖分中的最小角相等。
新增、删除、移动某一个顶点只会影响邻近的三角形。
与泰森多边形(Voronoi图)对偶。
算法
首先回顾一下前端开发模式的演进,我觉得主要有四个阶段。实现狄洛尼三角剖分的算法有很多,我们这里介绍其中一种:分治算法。这种算法主要过程是递归的分割点集直到子集大小不超过三,然后在合并的过程中逐渐选取最优点连接左右子集,最终完成三角剖分。
一般分为一下几个步骤完成。
将所有点的数组按照x坐标升序排序,如下图是排好序的点集。
将有序的点集递归的分成两个部分,直到子集大小不超过3。然后将这些子点集剖分成为一个三角形或者一条线段。
把已经剖分好的左右子点集可以依次合并。合并后的剖分包含左子点集的边(红色),右子点集的边(绿色),连接左右剖分产生的新的边(蓝色)。对于合并后的三角剖分,为了维持狄洛尼三角剖分原则,我们可能需要删除部分左子集边和右子集边。
合并左右子集首先是要确定一条基线(base edge),基线是最底部的且不与任何左右子集边相交的一条边。
接着就是确定紧接当前基线的下一条基线。以右子集为例,下一条基线为为右子集中的可能点与当前基线左端点构成。对于这些可能的端点,我门有两个检验标准:当前基线与可能点的夹角小于180°,当前基线两个端点与这个可能点的外接圆不能包好其他的可能点。
如上图所示,点I包含是其他可能点不能作为下一个基线的端点,点G不包含其他可能点可以作为下一个基线的端点。同时,左子集也做同样的处理。
当左右子集都不存在符合标准的可能点时,合并完成。在添加下一条基线的时,同时删除与其相交的左右子集的边。
当左右子集均存在可能点时,判断左子集的可能点的外接圆是否包含右子集的点,如果包含则不符合,右子集同理。如下图,点D的外接圆包含右子集的点,点G的外接圆不包含左子集的点,所以点G作为下一个基线的端点。
重复以上步骤,直到合并完成。
以下为主要步骤的代码实现:
/** * 下面有几个工具函数没有列出,这里说明一下 * dist(a, b) - 计算两点之间的距离 * cross(0, a, b) - 计算向量OA和OB的叉积 * inCircle(a, b, c, p) - 点P是否在点a,b,c三点的外接圆内 * interp(a, b, c, d) - 判断线段ab,cd是否相交 */ function delaunayTriangulation(l: number, r: number, points: Point[], edges: Edge[]) { if (r - l <= 2) { // 将所有的点两两相连 for (let i = l; i <= r; i++) { for (let j = i + 1; j <= r; j++) { addEdge(i, j, edges); } } return; } let mid = ~~((l + r) / 2); delaunayTriangulation(l, mid, points, edges); delaunayTriangulation(mid + 1, r, points, edges); let founded = false; let currentL = l, currentR = r; // 查找基线(base edge) while (!founded) { founded = true; const pL = points[currentL], pR = points[currentR]; // 查找左子集最底部的点 for (let i = 0; i < edges[currentL].next.length; i++) { const tIdx = edges[currentL].next[i]; const t = points[tIdx]; const d = cross(pR, pL, t); if (d > 0 || (d === 0 && dist(pR, t) < dist(pR, pL))) { currentL = tIdx; founded = false; break; } } if (!founded) continue; // 查找右子集最底部的点 for (let i = 0; i < edges[currentR].next.length; i++) { const tIdx = edges[currentR].next[i]; const t = points[tIdx]; const d = dir(pL, pR, t); if (d < 0 || (d === 0 && dist(pL, t) < dist(pR, pL))) { currentR = tIdx; founded = false; break; } } } // 添加base edge addEdge(currentL, currentR, edges); while (true) { const pL = points[currentL], pR = points[currentR]; let potential = -1, side = 0; // 查找可能点, 可能点与base edge小于180°, 可能点和base edge两个端点构成的圆内不能包含其他的点 for (let i = 0; i < edges[currentL].next.length; i++) { const tIdx = edges[currentL].next[i]; const t = points[tIdx]; if ( cross(pL, pR, t) > 0 && (potential == -1 || inCircle(pL, pR, points[potential], points[tIdx]) < 0) ) { (potential = tIdx), (side = -1); } } for (let i = 0; i < edges[currentR].next.length; i++) { const tIdx = edges[currentR].next[i]; const t = points[tIdx]; if ( cross(pR, t, pL) > 0 && (potential == -1 || inCircle(pL, pR, points[potential], points[tIdx]) < 0) ) { (potential = tIdx), (side = 1); } } // 不存在任何可能点,合并完成 if (potential == -1) break; if (potential !== -1) { // 删除与下一条base edge相交的边 if (side == -1) { for (let i = 0; i < edges[currentL].next.length; i++) { const tIdx = edges[currentL].next[i]; const t = points[tIdx]; const d = edges[tIdx]; if (interp(pL, t, pR, points[potential])) { d.next.splice(d.next.indexOf(currentL), 1); edges[currentL].next.splice(i, 1); i--; } } currentL = potential; } else { for (let i = 0; i < edges[currentR].next.length; i++) { const tIdx = edges[currentR].next[i]; const t = points[tIdx]; const d = edges[tIdx]; if (interp(pR, t, pL, points[potential])) { d.next.splice(d.next.indexOf(currentR), 1); edges[currentR].next.splice(i, 1); i--; } } currentR = potential; } // 添加下一个base edge addEdge(currentL, currentR, edges); } } }
相关链接
https://en.wikipedia.org/wiki/Delaunay_triangulation
http://www.geom.uiuc.edu/~samuelp/del_project.html
往期精彩回顾
苹果强制要求更换启动方式的解决方案
实习招聘|360云平台火热招聘中
360Stack裸金属服务器部署实践
360技术公众号
技术干货|一手资讯|精彩活动
扫码关注我们