12块高程瓦片造一层等高线：DEM 转换实测

技术专栏 · 数据处理 · 2026-10-14

本站一键地形图选完范围几秒钟就吐出一张带等高线的地形图，这条流水线的每个环节其实都不神秘。这篇把它拆开在本地跑了一遍：拉 12 块公开高程瓦片，解出 78.6 万个高程点，再用 marching squares 按 50 米间隔抽线——641 毫秒出 4,421 条等高线。中间还踩出一个教科书级 bug：插值取错了端点，让等值线碎成 5 万段。前后两版的数字都在这里。

一张 PNG 里藏着一座山

实验区选在天津蓟州北部山地（上篇天津村界里村级单元最多的区），z12 缩放级下 4×3 共 12 块瓦片，文件合计约 0.45MB。这类高程瓦片用 terrarium 编码：红绿蓝三个通道各存高度的一部分，逐像素套一行公式就还原出米制高程：

height = R * 256 + G + B / 256 - 32768;   // terrarium 编码，单位米

PNG 解码在 Node 侧没有顺手的零依赖方案，本文借 playwright 的 canvas 读像素，12 块瓦片拼成 1,024×768 的网格——786,432 个高程点，实测高程 3 米到 839 米，无一处无效值。按 100 米一档数一遍，平原比想象中还要平：

一半的地面在 100 米以下，700 米以上只占千分之一——这张表直接决定了后面等高线的密度分布：山地等高线密、平原等高线稀，不是画出来的风格，是地形本身。

marching squares：16 种组合切出等值线

等值线的数学定义简单：高程恰好处在某个值的点连成的线。网格数据上找它用 marching squares——对每个 2×2 单元，四个角与等值线水平比大小，大于记 1、小于记 0，四种角拼出 16 种组合；每种组合要么不穿线，要么穿一条，要么（对角相反时）穿两条，交点位置在边上按线性插值求出：

const a = h(x, y) - level, b = h(x+1, y) - level,     // 四角减去等值线水平
      c = h(x+1, y+1) - level, d = h(x, y+1) - level;  // 得到正负号组合
const idx = (a>0)*8 | (b>0)*4 | (c>0)*2 | (d>0)*1;
// 交点必须从线段的同一端起算（t = -首端值 / (末端值 - 首端值)）：
const top    = () => [x - a / (b - a), y];        // 上边：从 (x,y) 起
const right  = () => [x + 1, y - b / (c - b)];    // 右边：从 (x+1,y) 起
const bottom = () => [x - d / (c - d), y + 1];    // 下边：从 (x,y+1) 起
const left   = () => [x, y - a / (d - a)];        // 左边：从 (x,y) 起
switch (idx) {
  case 1: case 14: push(left(), bottom()); break;
  case 2: case 13: push(bottom(), right()); break;
  case 3: case 12: push(left(), right()); break;
  // ……其余 9 种组合同理
}

786,432 个单元 × 16 个级别，纯 JavaScript 跑完 641 毫秒。把所有单元的线段按端点拼成折线（端点按 0.05 像素量化配对），统计表如下：

折线数从低海拔向中山递减、顶点数却在 200–250 米档冲上峰值——山麓一圈圈小环多而短，中山带的大环每条都长。这张表本身就是地形剖面。

一个真 bug：5 万条碎线是怎么来的

这组干净的数字来之不易。第一版跑出来的是 49,657 条折线、87.5% 不足 5 个点、最长一条只有 99 个顶点——等值线全裂成了渣。原因是插值取错了端点：上边的交点写成了 x + b/(b-a)，而正确值是 x - a/(b-a)。两者在数值上关于单元中点镜像，只有当上下两个角值对称时才相等。相邻单元共享同一条边，一个从左端算、一个等效地从右端算，接缝处错开最多一个像元——每条等值线走到单元边界就断，全图碎成渣。

修正成「从线段同一端起算」之后，折线数从 49,657 降到 4,421，最长一条 6,499 个顶点——那是绕着主山体转了一整圈的长线；百点以上的长线 113 条。剩下 3,735 个短环也不是噪声：缓坡上高程在等值线水平附近抖动，一个单元自成一环是正常产物，出图时按长度过滤即可。这个 bug 教训和边界数据拓扑检查里的教训是同一条：相邻要素必须对共享边界上的同一个点给出同一个答案，顶点如此，插值也如此。

想自己跑一遍，路线有三条

脚本路线就是本文的全套：拉瓦片 → 套 terrarium 公式 → marching squares → 拼折线，一张山地图的等高线层 641 毫秒，全部代码不到两百行。桌面路线在 QGIS 里对高程栅格直接用「提取等值线」工具，间隔与属性字段都是现成选项。产品路线最省事——一键地形图把瓦片解码、等高线、山体渲染都搬进了浏览器，框选范围秒级出图，还能导 PNG。生成的等高线若是嫌顶点太多，瘦身方法在《GeoJSON 太大？先抽点，再谈小数位》里实测过；给高程分层配色的思路与《填色图配色一次做对：分档、选色与对比度》完全同源。

常见问题

高程瓦片里的高度是怎么编码的？

本文用的 terrarium 编码把高度拆进 PNG 的三个颜色通道：高度 = R×256 + G + B/256 − 32768。一块 z12 瓦片 256×256 像素就是 65,536 个高程点，12 块拼出 1,024×768 的网格。解码只要逐像素套公式，本文 78 万个点无一处无效值。

marching squares 生成等高线的原理是什么？

对每个高程网格单元，比较四个角与等值线水平的大小关系，16 种组合里最多穿过两条等值线段，交点位置按线性插值求出；把所有单元的线段按端点拼接，就是完整的等高线。本文 78.6 万个点的网格跑 16 个级别，641 毫秒完成，纯 JavaScript 无依赖。

等高线间隔选多少合适？

跟着用途和数据分辨率走。z12 瓦片一格约 30 米，50 米间隔能在山地出清晰纹理、在平原不至于过密，本文实测 16 个级别共 4,421 条折线。做展示可加大到 100 米只留计曲线；做地形分析可以加密，但要接受短碎环变多。

生成的等高线全是短碎段怎么办？

先查插值。本文第一版把交点比例从线段的错误一端取值，相邻单元在共享边上各算各的，接缝错位让等值线碎成 49,657 段、87.5% 不足 5 个点；修正为同一端起算后降到 4,421 条，最长一条 6,499 个顶点。修完仍存在的短环多为缓坡高程抖动的正常产物，可按长度过滤。

不想自己写代码，怎么快速得到等高线图？

用本站的一键地形图：框选范围即可生成带等高线的高程图，支持 15 套配色、乡镇标注与 PNG 导出，解码渲染均在浏览器本地完成。想控制每个细节再回到本文的脚本路线，两条路用的都是同一套高程数据源。