简介面向TIN不规则三角网的高效构建需求这套基于VC6.0的完整工程实现了被称为“最快”的Delaunay三角网生成算法提供可直接编译运行的C源码及演示程序适合GIS、测绘、地形可视化与计算机图形学方向的开发者学习或二次复用。资源共88个文件压缩包仅769KB核心代码由13个头文件与13个C源文件组成涵盖CTriangle、CEdge、CVertex等关键数据结构同时包含obj、sbr、pdb等编译过程文件以及exe可执行演示便于在Visual C 6.0中打开工程并对照调试。从输入对话框、绘制信息模块可以看出项目不仅实现了算法内核还附有点集输入、三角网渲染与结果显示等完整示例有助于理解Delaunay剖分中顶点插入、边翻转、三角网更新等核心步骤并可在Debug/Release两种模式下比较运行效率。目前已有1318人浏览学习适合欲深究快速TIN算法原理、或需要将算法移植到新项目的开发者参考借鉴。1. 最快的 TIN 三角网生成算法先别急着写代码无人机航测或地面激光扫描得到的点云动辄几十万上百万个点让这些点变成 TIN 三角网是等高线绘制、填挖方计算、地质建模共同的第一道工序。所谓最快的 TIN 三角网生成算法不是一个固定的程序名而是选型、数据结构和浮点处理的组合拳。同一个数据集有人十几秒出结果有人要等两分钟差的就是这套组合拳里的细节。这篇笔记先讲怎么选算法再给出一个最小实现骨架最后把我在实践里反复踩过的坑一起说透。适合手里已经有点云数据、正被生成耗时折磨的开发者如果你的数据只有几千个点读到这里可以先跳过十万点以上后面内容会直接变成你的参考路径。2. 分治 vs 增量先立住“最快”的算法选型2.1 空外接圆判据一条公式决定正确与速度TIN 三角网的标准骨架是 Delaunay 三角网任意一个三角形的外接圆内不包含点集中的其他点。这个“空外接圆判据”是几乎所有生成算法的总开关。判断流程通常是两步先确定三个顶点是逆时针排列再做圆内/圆外检测。别看它只有一行行列式它同时决定结果正确性和速度——点接近共圆时符号出错的代价是整个后续逻辑全乱。我在实际项目里见过最快也最惨的错误判据里用了不精确的比较导致该翻转的边没翻不该翻的边翻了。于是一条错误的边一直留在三角网里后续合并和查找全部错位速度随之崩坏。很多号称 O(n log n) 的算法跑成 O(n²)不是算法不行是判据没有给算法提供稳定输入。所以结论先放在这判据的正确性先于速度。先让空外接圆判据在“点接近共圆”时仍能给出正确符号再去谈算法快不快。这一步省不了否则后面所有优化都是给错误的网格加速。2.2 分治法为什么它是十万点以上最稳的“快”分治法的思路很直白点集先按 x 或 (x, y) 排序从中间切一刀左右两半分别三角化再把两块沿交界缝合成一个完整的 Delaunay 三角网。缝合要做两件事找到左右两个子集之间的公共切线然后沿切线逐一加入新三角形每加一步就查一次空外接圆性质不满足就立即翻转。因为每一层每个点只参与常数次翻转整体复杂度能做到严格的 O(n log n)。工程上分治法另一个好处是内存访问模式稳定。点按排序后左右两半天然落在连续区间合并时数据局部性好缓存命中率高在百万点随机分布数据上稳定实现能压到秒级。这也是我常用它作基准的原因不管最后选哪个算法先用一个标准分治实现测出基线后面的优化才有参照。分治的代价是代码量。递归、合并、切线、翻转每一部分都有边界条件。十万点以下它不一定赢过实现良好的增量法常系数摆在那但数据一旦过了百万复杂度的下限优势就会体现出来增量的最坏情况是 O(n²)随机分布也很难保证每次定位都在常数步内。2.3 Bowyer-Watson 增量法什么时候它能反超增量法的流程更直观先构造一个足够大的超级三角形把全部点罩住每插入一个新点先定位它落在哪个三角形里把那个三角形拆成三个新三角形再从这三个三角形出发检查相邻三角形是否违反空外接圆判据违反就翻转直到局部满足 Delaunay 性质。它快不快几乎全看“定位”这一步。常见做法是把三角形按网格或四叉树索引从点出发找三角形当输入数据本身有空间顺序比如扫描线顺序输出的点云定位几乎就是 O(1)整体速度能非常接近分治法。我实际见到的反超场景往往是数据点在空间上有强相关性增量算法利用上一轮留下的上下文移动窗口式插入缓存友好常数小得惊人。但增量法的翻车点也在定位和回溯。如果点坐标分布很随机又只用一个不划分的全局三角形表去线性扫每插一个点都扫大片三角形百万点数据会慢到让人怀疑人生。即使有索引当翻转链跨度很大时递归回溯的代价也远高于分治的局部优化。总结成一条经验十万点以下、数据有序增量法更省事百万点以上、随机分布我宁愿选分治。2.4 选型决策表什么数据选什么算法把上面几条放进一张表边界一眼就能扫出来。数据特征推荐算法理由1 万点以下增量法或直接枚举常数小实现短10 万~100 万随机分布分治法严格 O(n log n)行为稳定10 万~100 万扫描线有序增量法配网格索引定位 O(1)缓存命中好大量共线、重复点先排序去重再走主算法判据稳定性是前提需要增量更新结果增量法局部修改不用全量重建这套选型不是拍脑袋。理论复杂度相同时常数、内存布局和退化点处理决定胜负。算法选型这一关过了后面实现才有意义选错了后面再怎么调参都是给一辆跑偏的车换轮胎。3. 把分治法写成能跑的最小 TIN 生成器数据结构、递归与归并这里给一个最小实现骨架。我常用 Python 讲逻辑因为可读性好、能直接跑生产环境建议换 C 或 Rust数据结构思路完全一致后面我也会点出工程版要改的地方。3.1 数据结构点数组与三角形数组别用对象链表三角网生成里最隐性的坑是把每个三角形建成独立对象、再用指针互相指。每次翻转要改三四个指针合并时还要从链表头找邻居日积月累就是慢。数组顺序存储才是正解。class Vert: __slots__ (x, y, idx) def __init__(self, x, y, idx): self.x x self.y y self.idx idx class Tri: __slots__ (v0, v1, v2, nb0, nb1, nb2) def __init__(self, a, b, c): self.v0, self.v1, self.v2 a, b, c self.nb0 self.nb1 self.nb2 None每个三角形存三个顶点和三个邻居。邻居不带左右区分按顺序对应nb0是v0-v1这条边对面的三角形nb1对应v1-v2nb2对应v2-v0。合并操作里最高频的动作是“拿到一个三角形沿某条边找邻居再走过去”这个结构一次就能查到。参数说明顶点idx保留原始点编号最后写回结果时要靠它。三角形类的六个引用不要省很多实现因为“只存顶点、需要时全图搜索”而直接拖慢 LOP 阶段。工程版我建议把六个引用改成三个整数索引数组顶点坐标再用两个单独的浮点数组存放缓存能再上一个台阶。3.2 预处理排序、去重、再谈递归分治的第一步不是递归而是预处理。我在同一段逻辑里做排序和去重一并返回去重后的新数组。def prepare(pts): pts sorted(pts, keylambda p: (p.x, p.y)) # x 主序y 次序 uniq [] last None for p in pts: if last and last.x p.x and last.y p.y: continue # 精确比较去重 uniq.append(p) last p return uniq按 (x, y) 字典序排序能让后续递归切分天然产生空间连续性。分治所有性能都建立在“切开的两个子集各占连续内存段”之上排序这一步宁可多花几十毫秒也要保证后面的局部性。去重是关键中的关键。重复点进入三角网会生成零面积三角形空外接圆判据在这种面上退化翻转逻辑会变成来回循环的鬼打墙。这里用精确比较不吃 epsilon坐标来自采集端时通常已经整理过直接按明文去重即可。3.3 切分与基线条件小于等于多少点就暴力递归终止条件不是 n1而是一个阈值。点足够少时直接枚举所有三点组合用空外接圆判据筛出合法三角网比继续递归更快。MIN_T 12 def triangulate(pts): n len(pts) if n MIN_T: return brute_force_tri(pts) # 枚举三点组合保留 Delaunay 边 k n // 2 left triangulate(pts[:k]) right triangulate(pts[k:]) return merge(left, right, pts, k)阈值取 12 是我常用的起点。太小递归层数多函数调用和数组切片开销浪费太大暴力枚举的排列组合又会吞掉收益。我试过 8 到 20 之间差别不大从 12 开始调没问题。这段代码还要注意切片开销pts[:k]在 Python 里是复制数组并未共享内存。工程版我只传(start, end)区间避免每层复制。分治递归深度大约是 log2(百万)≈20 层栈没有风险真正的风险在归并层下面细说。3.4 合并与 LOP快慢分水岭在这里归并函数是分治实现真正的难点。两半递归各自生成了三角网合并要做三件事找下切线、沿切线推进、对每条新边做局部优化。合并骨架的示意代码如下。def merge(left, right, tris): L, R lower_tangent(left, right) # 左右子集的下侧公切线 while True: cl next_cand_left(left, L, R) # 左侧候选点 cr next_cand_right(right, L, R) # 右侧候选点 if cl is None and cr is None: break # 到达上切线合并结束 # 右候选若落入左候选的外接圆右候选更优否则左候选更优 if cr is not None and (cl is None or inside_circle(R, L, cl, cr)): add_triangle(L, R, cr) R cr else: add_triangle(L, R, cl) L cl local_optimize(L, R, tris) # 对邻接边做 LOP 翻转inside_circle(R, L, cl, cr)的含义是以 R、L、cl 三点做外接圆判断 cr 是否落在圆内。若 cr 在圆内说明 cl 对应的三角形不满足空外接圆判据要用 cr 这一侧前进。这个方向在教科书里经常一句话带过实际编码时最容易写反写反的后果就是合并不收敛程序卡死在同一组边上。local_optimize就是 LOP加入新三角形后检查它对侧的三角形如果对侧顶点落入新三角形外接圆就翻转公共边翻转后可能波及下一条邻接边于是循环检查直到局部全部合法。LOP 是 O(1) 的局部操作实现时别把它写成全图搜索否则分治的复杂度优势就白搭了。注意上面是教学骨架去掉了几何边界判断。生产实现里next_cand_left和next_cand_right要顺着凸包链走每走一步都要用一次方向判据确认“线两侧不再有点”这是合并步最容易漏的环节。4. 三角网生成避坑笔记五个必踩的坑和排查方法4.1 外接圆漏判与三角形撕裂浮点误差在作怪现象一百万个点跑完局部检查发现少了几条边三角形互相穿过面积出现负值。原因空外接圆判据用 double 直接算行列式共圆或近共圆时符号被舍入误差反转。该翻的边没有翻后续邻接关系就按错误网格传播越传越乱。解决把空外接圆判据统一成一个独立函数不内联散落。先用 double 算当行列式绝对值接近零比如小于 1e-12或相对机器精度的数十倍落回扩展精度重新算一遍。这个分支在普通数据上基本进不去几乎不影响性能却能把稳定性拉回正轨。4.2 递归栈溢出百万点退化成一条直线现象普通数据跑得好好的换成沿道路扫描的窄长条点云后程序直接栈溢出。原因分治按 x 切分但窄长条点云排序后“左半”和“右半”贴在同一侧子集仍然严重共线三角形退化合并时切线找来找去找不到归并又套递归栈就越压越深。解决排序时用 (x, y) 字典序切分阈值往上调到 16 或 32让退化点集在更早的暴力阶段结束平台栈浅就改成显式工作栈。这个坑一踩就是大半天建议开发期先拿窄长数据测一轮再谈百万点。4.3 重复点与共线点现象随机、原因明确现象同样一版数据有时少几个三角形有时多出几条零面积边结果不稳定。原因重复点没被预处理干掉共线点让多个三角形同时符合“空外接圆”判据同时进入候选集翻转循环在等价解之间来回试错。解决预处理阶段严格去重遇到共线点按“只选第一个进入决策”的策略处理LOP 阶段用“严格大于 0 才翻转”来终止循环。这个“严格大于 0”是许多实现翻车的静默点用的话边界会反复触发翻转结果就是看着正常、实际慢一倍。4.4 合并桥架错起始点选错导致死循环现象合并函数卡住日志停在同一条切线边上反复打印手动中断才发现程序没退出。原因找下切线时起点选了左子集的任意边界点或右子集的任意边界点两点连起来的那条线并没有真正把两个点集分开候选点始终落在线的错误侧循环来回换边永远不收敛。解决找切线时先各自取凸包边界链从最左、最右端点出发往下走每一步用一次方向判据确认“所有点都在连线同侧”这个方向判断单独写成一个函数。死循环九成出在切线方向先打印 L、R 坐标再打印方向判据结果五个排查能命中三个。4.5 比增量法还慢把分治写成了伪分治现象理论上是分治实际测出来比 Bowyer-Watson 慢一倍百万点跑了四十秒。原因切片用复制、每层都拷贝数组合并时又全图扫描寻找候选点LOP 也是全局搜索。本质上复杂度已经回到 O(n²)只是披着递归的外壳。解决用偏移量加长度传子区间候选点只从前一步邻接三角形里找别扫全子集LOP 每翻一条边都要快速找“对面”翻得慢就回头查邻接字段是否维护错位。分治法的面子是递归里子是局部性局部性丢了不如改用增量法。5. 验证阶段用三张表检验“快”和“对”再谈落地5.1 对照验证先证明没生成错跑通第一版后别急着测时间。我习惯与现成三角化实现做对照三角形数量、凸包顶点顺序、随机采样点的邻接关系。两份输出可能不完全一致但三角形总数和凸包必须一致。若数量差一两个多半是退化点处理差异若数量级都对不上说明结构有根本问题。对照这一步看起来笨却最省时间因为后面所有速度测试都建立在“正确”的前提上。5.2 规模曲线从 1 万到 100 万点的耗时趋势测速只看一组数据没有意义。我按 1 万、5 万、20 万、100 万四个规模测同一份随机数据看时间曲线是否近似 n log n若 100 万点出现悬崖式拉长大概率是索引或数据结构在某个规模上退化。数据点数分治实现耗时参考健康度判断1 万约几十毫秒正常5 万约几百毫秒正常20 万约 1~2 秒留意是否陡增100 万约 5~10 秒应接近 n log n 趋势这个表的绝对数值在每台机器上都不一样别把它当标准重点看趋势从 5 万到 100 万耗时不应出现几十倍的跳变。5.3 不变量三角形数量守恒与空圆抽查我每次调参后都会跑两条不变量。第一核对欧拉关系去重后点数 n、凸包顶点数 h三角形数 t 必须满足 t 2n - 2 - h差一个数都说明邻接没缝上。第二随机抽十万条边对每条边取对侧三角形顶点重跑空外接圆判据违规比例应接近零。这两条能兜住绝大多数“看着对、实际错”的结果。5.4 落到生产前的三个小改造和习惯验证之后真要做生产版本我还会做三件事把坐标减掉包围盒左下角转成相对值让数值量级集中在 1e-5 到 1e6 之间把递归里的数组切片改成(start, end)区间把三角形邻居的六七个引用换成三个整数索引。这三件事不改算法逻辑但能让 C 或 Rust 版本里的缓存命中和向量化顺利很多。到这一步才能真正放心去优化“最快”先有正确基线再有规模曲线最后做硬件级优化。我自己的习惯是把对照脚本、规模曲线和不变量检查写成同一个测试程序跑一次出一张报告一切调参都用同一套报告比对。这个习惯让我少走了很多弯路也希望帮到你。本文还有配套的精品资源点击获取
