第09章 - 叠加分析详解

叠加分析（overlay analysis）是 GIS 里出现频率最高的一类几何运算：把两层矢量数据按空间位置叠在一起，求交集、并集、差集，从而回答“两块地重合的部分在哪”“这条管线穿过哪些行政区”“这个园区去掉水域之后剩多少面积”。在 JTS 中，这一族操作由 Geometry.intersection、Geometry.union、Geometry.difference、Geometry.symDifference 四个方法对外暴露，内部则牵出两套完全不同的算法实现：2016 年前后的旧引擎 OverlayOp，以及 1.18 引入、持续演进的 OverlayNG。

本章的目标有三个层次。第一层是会用：搞清四个方法的语义、返回类型、空几何与几何集合的处理规则，知道什么调用会抛异常、什么调用会静默给出错误结果。第二层是懂原理：拆开 OverlayNG 的流水线，看 noding、成边、建图、标注、构面这五步各自解决什么问题，理解 OverlayLabel 为什么要把边分成边界边、塌缩边、线边、非本几何边四种状态。第三层是能做工程决策：在浮点精度与固定精度之间取舍，知道 OverlayNGRobust 的重试阶梯在什么情况下救得了你，知道 CoverageUnion 什么时候能把并集从秒级压到毫秒级、又会在什么输入上悄悄算错。

阅读本章需要的前置知识：第 04 章的几何对象模型（Geometry、Polygon、MultiPolygon、GeometryFactory）、第 05 章的几何构造与 PrecisionModel、第 08 章的空间关系谓词与 DE-9IM。如果你只想先把功能跑起来，可以重点看 9.2、9.9、9.10 三节；如果你正被 TopologyException 折磨，直接跳到 9.7 与 9.11。

版本说明：本章以 JTS 1.21.0-SNAPSHOT 源码为准，示例代码在 1.18 以上均可编译。文中标注“1.18 新增”“1.20 新增”的特性，在更早版本中不存在。

9.1 叠加分析的数学语义

9.1.1 从集合论看四种运算

JTS 把几何对象看作平面上的点集（point set）。一个 Polygon 是它闭包内的所有点，一条 LineString 是其上所有点构成的连续点集，一个 Point 是单点集。在这个视角下，四个叠加运算就是集合的四种基本运算：

• 交（intersection）：属于 A 且属于 B 的点，记作 A ∩ B。
• 并（union）：属于 A 或属于 B 的点，记作 A ∪ B。
• 差（difference）：属于 A 但不属于 B 的点，记作 A − B。注意它不对称，A − B 与 B − A 是两个不同结果。
• 对称差（symDifference）：属于 A 或属于 B，但不同时属于两者的点，记作 A △ B，等价于 (A − B) ∪ (B − A)。

下面这张图用两个边长为 10 的正方形说明四种运算的几何形态。A 的坐标范围是 [0,10]×[0,10]，B 是 [5,15]×[5,15]，两者重叠一个 5×5 的方块：

        ┌───────────────┐
        │       A       │
        │      ┌────────┼────────┐
        │      │  A ∩ B │   B    │
        │      │   25   │        │
        └──────┼────────┘        │
               │                 │
               └─────────────────┘

  A            = 100        B    = 100
  A ∩ B        =  25        A ∪ B = 175
  A − B        =  75        B − A =  75
  A △ B        = 150（两个互不相连的 75 方块）

9.1.2 四种运算的语义对照

下表把数学定义、JTS 方法、结果维度语义与常见用途放在一起，方便对照查阅：

注意“结果维度”一列用的是“小于等于”。原因在于结果可能发生降维：两个多边形只在一条边上相接，其交集是一条 LineString；两个多边形只在一个角点相接，其交集是一个 Point。这就是所谓拓扑塌缩（topology collapse），是叠加分析里大量“意外结果”的根源。

9.1.3 代数不变量与它们为何会失效

在理想数学世界里，四个运算满足一组漂亮的恒等式：

• 交换律：A ∩ B = B ∩ A，A ∪ B = B ∪ A，A △ B = B △ A。
• 结合律：(A ∪ B) ∪ C = A ∪ (B ∪ C)，(A ∩ B) ∩ C = A ∩ (B ∩ C)。
• 分配律：A ∩ (B ∪ C) = (A ∩ B) ∪ (A ∩ C)。
• 对偶：A − (B ∪ C) = (A − B) ∩ (A − C)；A − (B ∩ C) = (A − B) ∪ (A − C)。

但在浮点坐标上，这些恒等式只是“近似成立”。原因是每做一次叠加，输出顶点都会被重新计算——要么保留原始坐标（浮点精度），要么四舍五入到网格（固定精度）。一旦顶点坐标变化，下一次运算看到的就是另一个几何，于是出现一些工程上很恼人的现象：

1. (A ∪ B) − B 不一定等于 A − B。并集可能把 B 的某个顶点吸附到 A 的边上，使得后续求差时多切掉一小块。
2. 交换律在数值层面成立，在坐标输出顺序上不成立。A.union(B) 与 B.union(A) 的顶点集合通常相同，但环的起点、方向、以及 MultiPolygon 内各元素的顺序可能不同。做数据比对时不能直接比较 WKT 字符串，应先 norm() 再比较，或比较面积、周长、顶点集合等不变量。
3. 面积不严格可加。对浮点结果，(A ∪ B).getArea() 与 A.getArea() + B.getArea() - A.intersection(B).getArea() 只在小数点后若干位上相等。

因此，工程上评价叠加结果时，应优先检查下面这组更稳健的不变量，而不是逐字符比较几何：

9.1.4 边界归属约定

JTS 遵循 OGC Simple Features 的边界约定：一个几何由内部（interior）、边界（boundary）、外部（exterior）三部分组成，三者互不相交。对叠加运算是封闭点集运算，边界点的归属由运算本身决定：

• 交集包含 A 与 B 的公共边界（例如两个相邻多边形的公共边属于交集，哪怕交集面积为 0，结果会是一条线）。
• 并集只包含边界一次，公共边在结果中不会重复出现。
• 差集保留 A 的完整闭包，再减去 B 的闭包；A 与 B 的公共边界会从结果中移除。

这解释了为什么 A.difference(B) 在 A、B 相邻时得到的是“去掉公共边之后”的多边形——公共边既不属于 A − B 的内部，也不属于它的边界。若需要保留这条边（例如做地块合并时保留历史界址线），应改用 LineMerger 或在线状图层上另行处理，而不能指望差集。

9.2 顶层 API：Geometry 的四个方法

9.2.1 方法签名与返回约定

四个叠加方法都定义在 Geometry 抽象基类上，签名一致——接受一个 Geometry 参数，返回一个新的 Geometry，不修改任何一方：

import org.locationtech.jts.geom.Geometry;

// 交集
public Geometry intersection(Geometry other);
// 并集（二元）
public Geometry union(Geometry other);
// 差集：this − other
public Geometry difference(Geometry other);
// 对称差
public Geometry symDifference(Geometry other);

另有三个一元/便捷方法容易混淆，先在此厘清：

• union() 无参版本：对本几何的所有组成元素求并，允许输入是 GeometryCollection，返回结果一定是最高的输入维度（多边形输入一定得到 Polygonal）。详见 9.9。
• buffer(double distance)：缓冲区，见第 10 章。它常被当成“修复工具”使用（geom.buffer(0)），但不是叠加运算。
• convexHull()：凸包，见第 12 章。

返回值类型遵循“最小充分类型”原则：结果是单个多边形就返回 Polygon，多个不相连的多边形返回 MultiPolygon，混合维度（例如一个面加一条线）返回 GeometryCollection。写代码时不要假设返回类型，用 getGeometryType() 或 instanceof 判断。

9.2.2 第一个完整示例

下面这段代码演示四种运算，并打印结果的类型、面积、周长与元素个数。输出值来自实际运行，可直接对照：

import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.io.WKTReader;

public class OverlayBasic {
    public static void main(String[] args) throws Exception {
        WKTReader reader = new WKTReader();
        Geometry a = reader.read("POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0))");
        Geometry b = reader.read("POLYGON ((5 5, 15 5, 15 15, 5 15, 5 5))");

        report("A ∩ B", a.intersection(b));
        report("A ∪ B", a.union(b));
        report("A − B", a.difference(b));
        report("B − A", b.difference(a));
        report("A △ B", a.symDifference(b));
    }

    static void report(String label, Geometry g) {
        System.out.println(label + " => " + g);
        System.out.printf("    type=%s dim=%d area=%.1f length=%.1f num=%d%n",
                g.getGeometryType(), g.getDimension(), g.getArea(),
                g.getLength(), g.getNumGeometries());
    }
}

运行结果：

A ∩ B => POLYGON ((5 10, 10 10, 10 5, 5 5, 5 10))
    type=Polygon dim=2 area=25.0 length=20.0 num=1
A ∪ B => POLYGON ((10 5, 10 0, 0 0, 0 10, 5 10, 5 15, 15 15, 15 5, 10 5))
    type=Polygon dim=2 area=175.0 length=60.0 num=1
A − B => POLYGON ((10 5, 10 0, 0 0, 0 10, 5 10, 5 5, 10 5))
    type=Polygon dim=2 area=75.0 length=40.0 num=1
B − A => POLYGON ((5 10, 5 15, 15 15, 15 5, 10 5, 10 10, 5 10))
    type=Polygon dim=2 area=75.0 length=40.0 num=1
A △ B => MULTIPOLYGON (((10 5, 10 0, 0 0, 0 10, 5 10, 5 5, 10 5)), ((10 5, 10 10, 5 10, 5 15, 15 15, 15 5, 10 5)))
    type=MultiPolygon dim=2 area=150.0 length=80.0 num=2

可以验证两个细节：并集的周长 60 而不是 40+40=80，因为重叠部分的边被溶解了；对称差的周长 80 正好是 A、B 两个“L 形”的边界之和，因为两者在 (10,5)–(5,10) 处仅以点相接，不共享边。

9.2.3 空几何的短路规则

GeometryOverlay 在进入算法前先处理空几何，避免无谓的图构建。规则可以总结为：

“副本”意味着返回的是 copy() 而非原对象引用，但坐标与结构完全相同。空结果的维度规则由 OverlayOp.createEmptyResult / OverlayUtil.resultDimension 决定，例如 POLYGON EMPTY 与 POINT EMPTY 求交，得到的是 LINESTRING EMPTY 之外的哪一种，取决于输入维度而非运算。实际写代码时建议统一用 result.isEmpty() 判断，不要依赖空结果的类型。

9.2.4 GeometryCollection 的限制

四个方法对 GeometryCollection 的态度并不一致，这是实践中很容易踩的坑。用第 04 章的思路理解：JTS 的叠加引擎要求输入“同质”（homogeneous），即一个 MultiPolygon 里的元素必须都是多边形，一个 MultiLineString 里必须都是线。

• 二元并、差、对称差：遇到非空 GeometryCollection 直接抛 IllegalArgumentException，消息为 Operation does not support GeometryCollection arguments。
• 交集：GeometryOverlay.intersection 对 GeometryCollection 做了特殊处理——把集合拆开，逐元素与另一个几何求交，再用 GeometryCollectionMapper 合并。所以它“能跑”，但结果是逐元素相交的拼接，可能包含大量低维碎片。
• 一元并集 union()：完整支持异构 GeometryCollection，是清理多类型数据的推荐入口（见 9.9）。

下面这段代码演示三种行为。第二个 catch 分支的异常消息是实际输出，可以据此在日志里快速识别问题类型：

import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.io.WKTReader;

public class OverlayWithCollection {
    public static void main(String[] args) throws Exception {
        WKTReader reader = new WKTReader();
        Geometry poly = reader.read("POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0))");
        Geometry gc = reader.read(
                "GEOMETRYCOLLECTION (POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0)), LINESTRING (20 20, 30 30))");

        // 交集：逐元素映射，能算
        System.out.println("intersection => " + poly.intersection(gc));

        // 一元并集：支持异构集合
        System.out.println("unary union  => " + gc.union());

        // 二元并集：拒绝 GeometryCollection
        try {
            poly.union(gc);
        } catch (RuntimeException ex) {
            System.out.println("union failed => " + ex.getClass().getSimpleName()
                    + ": " + ex.getMessage());
        }
    }
}

输出：

intersection => POLYGON ((10 0, 0 0, 0 10, 10 10, 10 0))
unary union  => GEOMETRYCOLLECTION (LINESTRING (20 20, 30 30), POLYGON ((0 0, 0 10, 10 10, 10 0, 0 0)))
union failed => IllegalArgumentException: Operation does not support GeometryCollection arguments

一个实用的规避手法是：先按维度把集合拆成面层、线层、点层（PolygonExtracter、LineStringExtracter、PointExtracter 或 GeometryExtracter），分别在各自图层上做叠加，最后按需合并。这既符合业务语义，也避开了异构输入的限制。

9.3 两套引擎与全局开关

9.3.1 为什么会有两套引擎

JTS 的叠加分析最初只有一套实现，就是今天的“旧引擎”OverlayOp。它建立在 geomgraph 包的平面图结构上：把所有输入的边插入一个带标签的平面图，做 noding（求交打断），再按标签规则遍历子图、收集结果环。这套代码从 JTS 1.x 一直用到 1.17，支撑了 PostGIS、GeoTools、Hibernate Spatial 等一大票下游项目，稳定性经过了长期验证。

但它有两个绕不过去的短板。第一，精度问题：旧引擎的 noding 在双精度浮点上做，两条线段“数学上应该相交在一点”，浮点计算却可能给出两个略有差异的交点，于是 noding 不完整（non-noded intersection），最后在构面时抛 TopologyException。第二，性能与可维护性：平面图结构推理困难，遇到塌缩边、重合边时要靠大量启发式补丁。

1.18 引入的 OverlayNG（NG = Next Generation）从设计上就是把“精度”摆在第一位的重写。它的核心决策是：noding 必须在一个明确的精度模型下进行。若精度模型是固定的（PrecisionModel.FIXED），就用 Snap-Rounding Noder，把交点吸附到网格上，从根上保证 noding 完整；若精度模型是浮点的，就退回到普通的 MCIndexNoder，用校验器兜底，仍可能抛异常。

9.3.2 GeometryOverlay：隐藏的分发器

Geometry 的四个方法并不直接调用算法，而是委托给包级私有类 org.locationtech.jts.geom.GeometryOverlay。这个类做三件事：空几何短路、GeometryCollection 特殊处理、按全局开关选择引擎。它的核心分支如下：

                   Geometry.union(a, b)
                            │
                            ▼
                ┌───────────────────────┐
                │   GeometryOverlay     │
                │  （geom 包，包级私有） │
                └───────────┬───────────┘
                            │  isOverlayNG ?
              ┌─────────────┴──────────────┐
        false │                            │ true
              ▼                            ▼
  ┌──────────────────────┐      ┌───────────────────────┐
  │ SnapIfNeededOverlayOp│      │   OverlayNGRobust     │
  │  ├─ OverlayOp（旧）  │      │   └─ OverlayNG + 重试 │
  │  └─ SnapOverlayOp    │      └───────────────────────┘
  └──────────────────────┘

开关由系统属性 jts.overlay 控制，取值只有两个：old（默认）与 ng。源码中对应 GeometryOverlay.OVERLAY_NG_DEFAULT = false，也就是说在 1.21.0-SNAPSHOT 中，高层 union/intersection/difference/symDifference 默认仍然走旧引擎。这一点与很多教程中“JTS 已全面切换到 OverlayNG”的说法不一致，写代码时必须自己确认。

可以用一行代码在启动时切换：

# 让所有 Geometry 高层叠加方法使用 OverlayNG（带鲁棒重试）
java -Djts.overlay=ng -cp jts-core-1.20.0.jar:app.jar com.example.Main

# 显式回退到旧引擎（也是默认值）
java -Djts.overlay=old -cp jts-core-1.20.0.jar:app.jar com.example.Main

如果你不想依赖全局状态（在多租户服务、并发测试里尤其危险），就不要碰这个属性，而是在需要的地方显式调用 OverlayNG 或 OverlayNGRobust 的静态方法。这样每个调用点的行为都是确定的。

9.3.3 引擎对照表

9.3.4 该选哪个引擎

实践中的建议可以归纳成三条。第一，新代码优先显式使用 OverlayNGRobust：它把重试逻辑内置好了，输入再脏也有很大概率算出一个合法结果，代价是失败时比一次 OverlayNG 慢（因为要重试）。第二，批量数据处理优先考虑固定精度：把坐标按业务精度（例如 0.001 米）建立 PrecisionModel，一次性舍入后再算，既快又可复现。第三，不要混用：一个数据集里如果一半走旧引擎一半走新引擎，两者对塌缩边的处理规则不同，合并结果时会出现难以解释的缝隙与重叠。

至于结果差异，后面 9.12 节会用同一组输入给出对照。这里先记住一个结论：对“正常”的输入（本身合法、无近重合边），两套引擎的结果在拓扑上等价，差异主要体现在顶点坐标的舍入与环的起点方向；对“病态”输入，差异可能直接表现为一个能算、一个抛异常。

9.4 旧引擎 OverlayOp 的流程

9.4.1 流水线总览

虽然新项目不建议主动使用旧引擎，但读懂它的流程有两个好处：一是看到 TopologyException 时能判断异常从哪一步冒出来；二是理解 SnapIfNeededOverlayOp 为什么那样重试。旧引擎的主体在 org.locationtech.jts.operation.overlay.OverlayOp，流程大致如下：

 输入 geom0, geom1, opCode
        │
        ▼
 ① OverlayNodeFactory + PlanarGraph
        │  把两侧的边与节点加入图
        ▼
 ② EdgeNodingBuilder / EdgeSetNoder
        │  求交打断（noding），生成带 Label 的边
        ▼
 ③ OverlayOp.computeOverlay()
        │  按 opCode 计算每个节点的标签
        ▼
 ④ 遍历 MaximalEdgeRing → MinimalEdgeRing
        │  收集属于结果的面
        ▼
 ⑤ PointBuilder / LineBuilder / PolygonBuilder
        │  分别构造结果中的点、线、面
        ▼
   buildGeometry(...) → 结果

其中第 ② 步是异常高发区。旧引擎用 RobustLineIntersector 计算线段交点，再把交点插入到边链中；如果同一对线段被算出两个“几乎相同但不相等”的交点，插入后图里就会出现两条本该相交却错开的边，第 ④ 步构环时便无法闭合，于是抛出：

org.locationtech.jts.geom.TopologyException: found non-noded intersection between
  LINESTRING ( 0 0, 10 10 ) and LINESTRING ( 10 0, 0 10 ) [ (5.0, 5.0, NaN) ]

这条消息是经典款，几乎每个用过 JTS 的人都见过。它告诉你：在 (5,5) 附近，两条线段本应相交但 noding 没成功。注意方括号里给出的是引擎认为“应该存在交点”的位置估计，实际计算出的交点可能差了几个 ulp。

9.4.2 Self-intersecting 输入触发的异常

最容易被复现的异常来源是自相交多边形（俗称蝴蝶结，bow-tie）。下面这段代码直接对蝴蝶结多边形做交运算，异常消息就是上一节展示的那条：

import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.io.WKTReader;

public class BowtieOverlay {
    public static void main(String[] args) throws Exception {
        WKTReader reader = new WKTReader();
        Geometry bowtie = reader.read("POLYGON ((0 0, 10 10, 10 0, 0 10, 0 0))");
        Geometry square = reader.read("POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0))");

        System.out.println("bowtie isValid = " + bowtie.isValid());
        try {
            System.out.println(bowtie.intersection(square));
        } catch (RuntimeException ex) {
            System.out.println("intersection failed: " + ex.getMessage());
        }

        // buffer(0) 是常用的“合法化”手段
        Geometry fixed = bowtie.buffer(0);
        System.out.println("buffer(0) => " + fixed);
        System.out.println("fixed isValid = " + fixed.isValid());
    }
}

输出：

bowtie isValid = false
intersection failed: found non-noded intersection between LINESTRING ( 0 0, 10 10 ) and LINESTRING ( 10 0, 0 10 ) [ (5.0, 5.0, NaN) ]
buffer(0) => POLYGON ((5 5, 10 10, 10 0, 5 5))
fixed isValid = true

这里有两个值得记住的点。其一，buffer(0) 只保留“最大面积的那一瓣”，所以蝴蝶结被削成了一个三角形（面积 25），另一瓣被丢弃——这不是通用的修复方案，第 17 章会给出更合适的 GeometryFixer。其二，SnapIfNeededOverlayOp 的重试在这里没能救场：它会用 GeometrySnapper 的容差重新吸附后再算一遍，但这个容差（由 GeometrySnapper.computeSizeBasedSnapTolerance 给出，等于输入包络短边 × 1e-9，本例为 1e-8）远小于蝴蝶结两瓣的交叠尺度（这里是 5 个单位），吸附后依然是同一条自交边，所以第二次尝试照样失败，最终把原始异常抛给调用者。

9.4.3 SnapIfNeededOverlayOp 的重试机制

SnapIfNeededOverlayOp 的逻辑非常直白，可以概括为“一次正式 + 一次补救”：

// 伪代码，对应 SnapIfNeededOverlayOp.getResultGeometry
try {
    result = OverlayOp.overlayOp(g0, g1, opCode);   // 直接用原始坐标算
    isSuccess = true;
} catch (RuntimeException ex) {
    savedException = ex;                           // 记住原始异常
}
if (!isSuccess) {
    try {
        result = SnapOverlayOp.overlayOp(g0, g1, opCode);  // 吸附后再算
    } catch (RuntimeException ex) {
        throw savedException;                      // 还是失败就抛最初那个
    }
}

SnapOverlayOp 内部做了三件事：先用 CommonBitsRemover 去掉两个几何的公共高位（把大坐标平移成小坐标，减少精度损失），再用 GeometrySnapper.snap 把两边的顶点互相吸附，然后交给 OverlayOp 计算，最后把公共位加回去。吸附容差由 GeometrySnapper.computeOverlaySnapTolerance(g0, g1) 给出：取两个几何各自的“基于尺寸容差”中较小者，而单个几何的基于尺寸容差是 min(width, height) × 1e-9；若几何使用固定精度模型，则容差至少为 (1/scale) × 2 / 1.415。

这解释了为什么旧引擎对“差一点点就重合”的输入很有效（例如两个测绘单位分别数字化同一条界址线，顶点偏差 1e-7 米），而对“尺度级别”的病态输入无能为力——容差是相对尺度的，不是相对错误的。

这条重试链的能力上限，正好解释了 OverlayNGRobust 存在的意义：它把“吸附一次”扩展为“吸附多次、容差递增、最后退回 Snap-Rounding”，覆盖了旧引擎够不着的场景。

9.5 OverlayNG 架构

9.5.1 五阶段流水线

OverlayNG 把叠加拆成五个阶段，每个阶段职责单一、数据形态明确。理解这条流水线，是读懂 TopologyException、读懂“为什么固定精度更慢但更稳”的关键。

        geom0, geom1, pm, opCode
                 │
                 ▼
  ┌──────────────────────────────────────────────┐
  │ 1. nodeEdges()  —— EdgeNodingBuilder         │
  │    · 抽取边（RingClipper/LineLimiter 可裁剪） │
  │    · 按 pm 选 Noder：SnapRounding / MCIndex   │
  │    · 丢弃完全塌缩的边，成 Edge 并按 EdgeKey 合并│
  └───────────────────┬──────────────────────────┘
                      ▼  List<Edge>（唯一、已 noding）
  ┌──────────────────────────────────────────────┐
  │ 2. buildGraph() —— OverlayGraph               │
  │    每条 Edge → 一对反向 OverlayEdge（共享 Label）│
  └───────────────────┬──────────────────────────┘
                      ▼
  ┌──────────────────────────────────────────────┐
  │ 3. labelGraph() —— OverlayLabeller             │
  │    computeLabelling() 求各边两侧 Location      │
  │    markResultAreaEdges(opCode) 标记结果边     │
  └───────────────────┬──────────────────────────┘
                      ▼
  ┌──────────────────────────────────────────────┐
  │ 4. extractResult() —— 按维度分层构结果         │
  │    PolygonBuilder → LineBuilder → PointBuilder│
  └───────────────────┬──────────────────────────┘
                      ▼
  ┌──────────────────────────────────────────────┐
  │ 5. ElevationModel.populateZ()（仅当有 Z）      │
  └───────────────────┬──────────────────────────┘
                      ▼
                   result

阶段 1 的产物是一组唯一的边：两条完全重合的边（比如两个相邻多边形共享的公共边）会被 EdgeMerger 合并成一条，并在标签里记录“它同时来自几何 0 与几何 1”。这一步决定了结果的正确性下限——如果 noding 没把该断开的边断开，后面所有步骤都白搭。

阶段 3 的产物是每条边两侧的 Location（INTERIOR / BOUNDARY / EXTERIOR）。OverlayLabeller 从一个已知位置的边出发，沿着节点周围的边按角度顺序推进，逐步推导出整张图所有边的位置。这是叠加分析里最“几何”的一步，也是性能瓶颈之一。

阶段 4 按运算类型决定要构哪些维度的结果。源码中 extractResult 的规则值得记住：先构面，若没有面再构线，若连线都没有且运算是交集，才构点。也就是说，交集的面优先于线、线优先于点，这是 isAllowMixedIntResult = !isStrictMode 与“最高维度优先”语义共同作用的结果。

9.5.2 Edge 与 OverlayLabel：四种边的状态

Edge 是“可见的一段线”，OverlayLabel 是“这段线相对于两个输入几何的拓扑身份”。标签用 dim 字段记录边的状态，只有四种取值：

这四种状态之所以够用，靠两条前提：输入必须同质（一条边不可能在同一个几何里既是塌缩边又是线边）；两个输入几何各自的维度可以不同，所以同一条边可以在几何 0 里是边界边、在几何 1 里是线边。

 一条边在两个输入几何中的标签组合（示意）
 ┌──────────────┬──────────────┬─────────────────────────┐
 │ 相对 geom0   │ 相对 geom1   │ 在结果中的角色           │
 ├──────────────┼──────────────┼─────────────────────────┤
 │ BOUNDARY/壳  │ NOT_PART     │ A 独有的外边 → 并集保留  │
 │ BOUNDARY/壳  │ BOUNDARY/壳  │ 重合边 → 按 opCode 决定  │
 │ BOUNDARY/壳  │ COLLAPSE     │ A 的边+B 的塌缩 → 需判定 │
 │ COLLAPSE     │ COLLAPSE     │ 双塌缩 → 通常丢弃        │
 │ LINE         │ NOT_PART     │ B 输入是面时 → 按位置判  │
 └──────────────┴──────────────┴─────────────────────────┘

OverlayEdge 则是一对方向相反的半边的容器：一条 Edge 生成两条 OverlayEdge（正向与反向），它们共享同一个 OverlayLabel。这样“边的左侧”与“边的右侧”就分别对应两条半边各自的视角，避免了左右信息在合并时丢失。

9.5.3 裁剪与短路优化

OverlayNG 默认开启优化（isOptimized = true）。优化分两类。

第一类是包络裁剪。交集的结果不可能超出两个输入包络的交集，差集的结果不可能超出 A 的包络，所以在这两种运算下，OverlayUtil.clippingEnvelope 会算出一个安全的裁剪包络，交给 RingClipper 把面裁剪进去、交给 LineLimiter 把线截断。安全外扩量是浮点精度下包络短边的 10%、固定精度下 3 倍网格边长，确保舍入后的顶点仍落在包络内。并集与对称差的结果范围就是两输入的全范围，无法裁剪，此时返回 null，不做优化。

第二类是结果面积启发式检查。浮点精度下，OverlayNG 在返回结果前会用 OverlayUtil.isResultAreaConsistent 校验一次：交集的面积不应大于任一输入，并集的面积不应小于任一输入，差集的面积应满足上下界关系，容差系数是 0.1。若明显违反，说明 noding 让图“翻转”了，直接抛 TopologyException("Result area inconsistent with overlay operation")。这条检查是 1.18 之后新增的，它把一部分原本静默的错误结果变成了显式异常，是好事。

FastOverlayFilter 是另一层短路：若两个几何的包络根本不相交，或其中一个完全包含另一个，则不必进入完整流水线。它只在“结果必然为空或必然等于某一输入”时才生效，属于纯加速手段。

9.5.4 用代码观察中间产物

OverlayNG 提供了三个调试用输出开关：setOutputNodedEdges(true) 输出 noding 之后的全部边，setOutputEdges(true) 输出进入图的边，setOutputResultEdges(true) 输出被标记为结果的边。它们返回的是 MultiLineString 而不是最终几何，非常适合排查“为什么这里多出一条线”。

import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.io.WKTReader;
import org.locationtech.jts.operation.overlayng.OverlayNG;

public class OverlayDebugEdges {
    public static void main(String[] args) throws Exception {
        WKTReader reader = new WKTReader();
        Geometry a = reader.read("POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0))");
        Geometry b = reader.read("POLYGON ((5 5, 15 5, 15 15, 5 15, 5 5))");

        OverlayNG union = new OverlayNG(a, b, null, OverlayNG.UNION);
        union.setOutputNodedEdges(true);
        System.out.println("noded edges  = " + union.getResult());

        OverlayNG inter = new OverlayNG(a, b, null, OverlayNG.INTERSECTION);
        inter.setOutputResultEdges(true);
        System.out.println("result edges = " + inter.getResult());

        // 关掉裁剪优化，对比结果是否变化
        OverlayNG noOpt = new OverlayNG(a, b, null, OverlayNG.INTERSECTION);
        noOpt.setOptimized(false);
        System.out.println("no-opt result = " + noOpt.getResult());
    }
}

输出（节选，noded edges 共 12 段）：

noded edges  = MULTILINESTRING ((0 0, 10 0, 10 5), (10 5, 10 0, 0 0), (10 5, 10 10, 5 10), ... , (5 5, 5 10))
result edges = MULTILINESTRING ((5 10, 10 10, 10 5), (10 5, 5 5), (5 5, 5 10))
no-opt result = POLYGON ((5 10, 10 10, 10 5, 5 5, 5 10))

注意 noded edges 里每条公共边都出现了两次、方向相反——这正是“一对反向 OverlayEdge”在调试输出里的直观体现。而 result edges 只剩三条边，正好围成那个 5×5 的交集方块。关掉优化后结果与开启优化完全一致，说明裁剪只影响性能、不影响语义（前提是安全包络计算正确）。

9.6 精度模型与 PrecisionReducer

9.6.1 浮点与固定精度

PrecisionModel 有两种类型：FLOATING（默认，直接用 double）与 FIXED（按固定比例尺舍入到网格）。两者的差异在叠加分析里会被放大，因为它决定了 noding 策略：

scale 的含义是“每单位长度分成多少格”，例如 new PrecisionModel(1000) 表示保留到 0.001。源码注释里有一条重要限制：固定精度的可靠计算上限约 13 位十进制有效数字（PrecisionUtil.MAX_ROBUST_DP_DIGITS = 14，留了一位余量）。因为坐标量级越大，可用于小数的有效位越少——量级 1e6 的坐标再乘 1e9 的 scale，就已经逼近 double 的 15 到 16 位有效数字极限。

9.6.2 三个精度工具函数

PrecisionUtil 提供了一组静态方法，作用是把“经验”变成“可计算的数”。最容易混淆的是 safeScale 与 robustScale：

对本章那对 0 到 15 量级的正方形，safeScale 给出 1e12，意味着“在这个量级下，保留 12 位小数仍然安全”；而 inherentScale 只给出 1.0，因为所有坐标都是整数——数据本身没有小数位，没必要用更高精度。robustScale 取两者中更保守的一个，所以得到 1.0。

这三个函数的取舍逻辑值得展开。safeScale 只关心“会不会溢出精度”，所以对整数坐标也会给出 1e12；inherentScale 只关心“数据长什么样”，所以浮点噪声也会被当成有效精度；robustScale 则是“在安全的前提下尽量保留输入精度”。实际工程里，如果你明确知道业务精度是 0.001 米，直接用 new PrecisionModel(1000) 比任何启发式都可靠。

9.6.3 固定精度叠加的实际效果

下面这段代码对同一对“歪斜”多边形分别做浮点叠加与两种固定精度叠加，并打印面积。输入刻意带了 9 位小数，方便观察舍入效果：

import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.geom.PrecisionModel;
import org.locationtech.jts.io.WKTReader;
import org.locationtech.jts.operation.overlayng.OverlayNG;

public class OverlayPrecision {
    public static void main(String[] args) throws Exception {
        WKTReader reader = new WKTReader();
        Geometry c = reader.read(
                "POLYGON ((0.123456789 0.1, 10.987654321 0.2, 10.5 10.4, 0.3 9.9, 0.123456789 0.1))");
        Geometry d = reader.read(
                "POLYGON ((5.111111111 5.2, 15.2 5.3, 15.4 15.6, 5.5 15.1, 5.111111111 5.2))");

        show("float    ", OverlayNG.overlay(c, d, OverlayNG.INTERSECTION));
        show("fixed 1  ", OverlayNG.overlay(c, d, OverlayNG.INTERSECTION, new PrecisionModel(1.0)));
        show("fixed 10 ", OverlayNG.overlay(c, d, OverlayNG.INTERSECTION, new PrecisionModel(10.0)));
    }

    static void show(String label, Geometry g) {
        System.out.println(label + " => " + g);
        System.out.printf("            area=%.15f numPts=%d%n", g.getArea(), g.getNumPoints());
    }
}

输出：

float     => POLYGON ((5.305373344890099 10.145361438475005, 10.5 10.4, 10.745937852375098 5.2558518069072875, 5.111111111 5.2, 5.305373344890099 10.145361438475005))
            area=27.319979782567877 numPts=5
fixed 1   => POLYGON ((5 10, 11 10, 11 5, 5 5, 5 10))
            area=30.000000000000000 numPts=5
fixed 10  => POLYGON ((5.3 10.1, 10.5 10.4, 10.7 5.3, 5.1 5.2, 5.3 10.1))
            area=27.000000000000000 numPts=5

三条结论可以直接写进工程规范。其一，固定精度的输出坐标全部落在网格上（5.3、10.7），这正是“结果可复现”的基础：同样的输入、同样的 scale，任何机器上都得到同一份 WKT。其二，精度越低，面积误差越大：scale=1（整米网格）把面积从 27.32 拉到了 30.00，误差接近 10%；scale=10（0.1 米）得到 27.00，误差约 1.2%。其三，浮点结果虽然“最接近真值”，但顶点坐标带着 15 位小数的噪声，做后续比较、入库、切片都会引入不确定性。

9.6.4 PrecisionReducer：只降精度、不改变拓扑意图

PrecisionReducer.reducePrecision(geom, pm) 是“把几何舍入到指定精度”的专用工具。它内部并不是简单地遍历坐标四舍五入，而是用 OverlayNG 的单输入模式（new OverlayNG(geom, pm)，opCode 固定为 UNION）在舍入后的精度下重新计算一遍——这样即使舍入导致自相交或相邻边粘连，输出也是合法几何。对面积维度的输入，它会设置 setAreaResultOnly(true)，避免塌缩出线状碎片。

import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.geom.PrecisionModel;
import org.locationtech.jts.io.WKTReader;
import org.locationtech.jts.operation.overlayng.PrecisionReducer;

public class ReducePrecisionDemo {
    public static void main(String[] args) throws Exception {
        WKTReader reader = new WKTReader();
        Geometry g = reader.read(
                "POLYGON ((0.2 0.3, 9.7 0.6, 9.9 9.8, 0.1 9.4, 0.2 0.3))");

        Geometry p1 = PrecisionReducer.reducePrecision(g, new PrecisionModel(1.0));
        Geometry p10 = PrecisionReducer.reducePrecision(g, new PrecisionModel(10.0));

        System.out.println("pm=1   => " + p1 + "  area=" + p1.getArea());
        System.out.println("pm=10  => " + p10 + "  area=" + p10.getArea());
        System.out.println("origin area=" + g.getArea());
    }
}

输出：

pm=1   => POLYGON ((0 0, 0 9, 10 10, 10 1, 0 0))  area=90.0
pm=10  => POLYGON ((0.2 0.3, 0.1 9.4, 9.9 9.8, 9.7 0.6, 0.2 0.3))  area=88.28
origin area=88.28

pm=10 时面积没变（88.28 恰好等于原值，因为原坐标本来就只有一位小数），只是环的方向与起点变了——又一次说明“面积、顶点集合”比“WKT 字符串”更适合做断言。pm=1 时面积从 88.28 变成 90，四个顶点被分别吸到最近的整数格点，这是可接受的（也是可预期的）精度损失。

注意：org.locationtech.jts.precision.GeometryPrecisionReducer 是另一个类，功能更全（支持点、线、面，支持按元素处理），1.18 起内部同样改用 OverlayNG 与 Snap-Rounding。两者不要混用：PrecisionReducer 只处理线面输入，遇到非法输入抛 IllegalArgumentException。

9.7 OverlayNGRobust 的重试阶梯

9.7.1 为什么需要“重试”

OverlayNG 的浮点模式并不鲁棒：它用的 MCIndexNoder 在 noding 出错时只会让 ValidatingNoder 抛异常，不会自动修复。而现实数据里，noding 出错是常态——不同来源的图层边界共线但顶点不重合、坐标经过多次投影转换后带着细微抖动、从 Shapefile 读进来的双精度值本身就不可靠。于是 JTS 提供了 OverlayNGRobust：一个把“多次尝试、逐级放宽”封装好的门面类。

它的算法思想非常朴素：先按最保真的方式算，算不出来就牺牲一点精度再算，直到算出来为止。源码的类注释把策略列为五条，对应到实现里就是三级阶梯：

 输入 geom0, geom1, opCode
        │
        ▼
 ┌─────────────────────────────────────────┐
 │ 第 1 级：OverlayNG + 浮点 noding          │  快，坐标改动最小
 │   失败 → 捕获原始异常保存起来             │
 └───────────────────┬─────────────────────┘
                     ▼
 ┌─────────────────────────────────────────┐
 │ 第 2 级：SnappingNoder，最多 5 轮         │
 │   容差 tol = max(|ordinate|) / 1e12       │
 │   每轮：① 两几何互相吸附后叠加            │
 │         ② 各自先自吸附、再叠加            │
 │   失败 → tol *= 10，重试                  │
 └───────────────────┬─────────────────────┘
                     ▼
 ┌─────────────────────────────────────────┐
 │ 第 3 级：SnapRoundingNoder + 安全比例尺   │
 │   scale = PrecisionUtil.safeScale(a,b)    │
 │   失败 → 抛出第 1 级保存的原始异常        │
 └───────────────────┬─────────────────────┘
                     ▼
                  结果 / 异常

9.7.2 三级策略的细节

第 1 级直接调用 OverlayNG.overlay(geom0, geom1, opCode)，用几何自身的精度模型（通常是浮点）。成功就返回，失败则把异常暂存起来——注意是“暂存”而不是“记录并继续”，因为如果后面所有策略都失败，抛给调用者的应该是这个最原始的异常，它的消息里带着最接近真实的交点位置，最有诊断价值。

第 2 级用 SnappingNoder，容差由 snapTolerance 计算：取两个几何包络中绝对值最大的 ordinate，除以 SNAP_TOL_FACTOR = 1e12。以本章的正方形为例，最大 ordinate 是 15，容差就是 1.5e-11。每一轮做两次尝试：先让两个几何互相吸附后叠加（overlaySnapping），再让每个几何先自吸附、然后叠加（overlaySnapBoth）。都失败就把容差乘 10 进入下一轮，最多 NUM_SNAP_TRIES = 5 轮。容差从 1e-11 一路放大到 1e-7，覆盖了“浮点噪声”到“数据录入误差”的大部分场景。

第 3 级退回 Snap-Rounding：用 PrecisionUtil.safeScale(a, b) 算出一个保证不超过 14 位有效数字的比例尺，建 PrecisionModel，再调一次固定精度的 OverlayNG。这一级会把坐标吸附到网格上，是最“暴力”也最可靠的一级。

值得强调的是：第 2 级的自吸附（snapSelf）内部用 OverlayNG(geom, null) 加 setStrictMode(true)，也就是“把几何与空几何做并集”。这相当于一次“自我 noding 加去重”，能把几何内部那些“几乎重合的相邻边”先粘在一起，从而消除大部分自相交类问题。这个技巧在工程上也可以单独使用：对一个脏的多边形执行 new OverlayNG(geom, null).getResult()，往往就能得到合法几何。

9.7.3 什么情况下它救不了你

OverlayNGRobust 并非万能。它救不了的三类情况：

1. 输入真正无效且错误尺度大于容差。两个多边形大面积重叠、自相交的蝴蝶结、环方向错误导致的“面积在外部”，这些不是精度问题而是数据问题，容差放大一万倍也没用。
2. 输入维度混合。OverlayNG 明确拒绝异构 GeometryCollection，异常消息是 Overlay input is mixed-dimension；这类参数错误不会因为重试而消失。
3. 拓扑语义本身不成立。例如对一对“一边完全在另一个内部”的多边形求差，结果为空是正确行为，不应期待异常；反之，如果你期待得到非空结果而得到空，那是业务假设错了，不是引擎错了。

   9.7.4 用代码观察三个入口的差异

OverlayNGRobust 没有暴露“第几级成功”的回调，但可以从结果结构上间接观察。下面这段代码对一对“顶点有微小偏差”的相邻多边形做并集，对比三个入口的返回结构：

import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.io.WKTReader;
import org.locationtech.jts.operation.overlayng.OverlayNG;
import org.locationtech.jts.operation.overlayng.OverlayNGRobust;

public class RobustRetryDemo {
    public static void main(String[] args) throws Exception {
        WKTReader reader = new WKTReader();
        // 两个本应相邻的多边形，接缝处顶点偏差 1e-9
        Geometry a = reader.read("POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0))");
        Geometry b = reader.read(
                "POLYGON ((10.000000001 0, 20 0, 20 10, 10.000000001 9.999999999, 10.000000001 0))");

        report("OverlayNG.union      ", OverlayNG.overlay(a, b, OverlayNG.UNION));
        report("OverlayNGRobust.union", OverlayNGRobust.overlay(a, b, OverlayNG.UNION));
        report("Geometry.union       ", a.union(b));
    }

    static void report(String label, Geometry g) {
        System.out.printf("%s => num=%d area=%.9f%n",
                label, g.getNumGeometries(), g.getArea());
    }
}

输出：

OverlayNG.union       => num=1 area=200.000000000
OverlayNGRobust.union => num=1 area=200.000000000
Geometry.union        => num=2 area=200.000000001

三个入口都能给出结果，但 Geometry.union 走的是旧引擎（默认 jts.overlay=old），结果是两个互不相连的多边形（num=2，因为 1e-9 的缝隙没有被识别为重合边）；而 OverlayNG 系列把缝隙正确地消除成了 num=1。面积差异只在小数点后第 9 位，肉眼几乎不可见，但几何结构完全不同——这正是“精度决策会改变拓扑结论”的典型例证。

如果想手工复现第 3 级的效果，可以直接用 PrecisionUtil.safeScale 构造安全比例尺：

import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.geom.PrecisionModel;
import org.locationtech.jts.io.WKTReader;
import org.locationtech.jts.operation.overlayng.OverlayNG;
import org.locationtech.jts.operation.overlayng.PrecisionUtil;

public class SafeScaleOverlay {
    public static void main(String[] args) throws Exception {
        WKTReader reader = new WKTReader();
        Geometry a = reader.read("POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0))");
        Geometry b = reader.read("POLYGON ((5 5, 15 5, 15 15, 5 15, 5 5))");

        double scale = PrecisionUtil.safeScale(a, b);
        System.out.println("safeScale   = " + scale);
        System.out.println("robustScale = " + PrecisionUtil.robustScale(a, b));

        Geometry g = OverlayNG.overlay(a, b, OverlayNG.INTERSECTION, new PrecisionModel(scale));
        System.out.println("result = " + g);
    }
}

输出：

safeScale   = 1.0E12
robustScale = 1.0
result = POLYGON ((5 10, 10 10, 10 5, 5 5, 5 10))

9.8 严格模式与结果维度

9.8.1 非严格模式下的混合维度结果

默认情况下（isStrictMode = false），OverlayNG 允许结果混合维度。最典型的例子是交集：两个 MultiPolygon 相交，既有一块面积重合的区域，又有一对只在角点相接的多边形，那么结果就是一个 GeometryCollection，里面既有 Polygon 又有 Point：

  A = MULTIPOLYGON( [0,10]×[0,10] , [20,30]×[0,10] )
  B = MULTIPOLYGON( [5,15]×[5,15] , [30,40]×[10,20] )

  A ∩ B ：第一对相交于 5×5 的面（面积 25）
          第二对仅在点 (30,10) 相接
          → 非严格模式：GEOMETRYCOLLECTION (POLYGON(25), POINT(30 10))
          → 严格模式  ：POLYGON(25)

这个行为继承了旧引擎的语义，对某些分析场景有用（例如你想同时知道“面接触”和“点接触”的位置），但对数据处理管线是灾难：一个 GEOMETRYCOLLECTION 传进下一个算子，多数算子会抛异常或返回空。

9.8.2 严格模式的三条规则

调用 setStrictMode(true) 后，结果语义收紧为三条：

1. 因拓扑塌缩产生的线、点不进入结果。例如两个多边形只共享一条边，交集在非严格模式下是 LINESTRING，严格模式下是空。
2. INTERSECTION 与 DIFFERENCE 的结果一定是同质的（homogeneous）。
3. 若两个输入维度相同，UNION 与 SYMDIFFERENCE 的结果也同质。

严格模式带来的好处是“可链式调用”：结果永远是干净的面或干净的线，不需要每次都判断“这次会不会又有低维碎片”。代价是丢失信息，因此严格模式适合做数据加工流水线，不适合做空间关系诊断。

9.8.3 代码对比

import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.io.WKTReader;
import org.locationtech.jts.operation.overlayng.OverlayNG;

public class StrictModeDemo {
    public static void main(String[] args) throws Exception {
        WKTReader reader = new WKTReader();
        Geometry mp = reader.read(
                "MULTIPOLYGON (((0 0, 10 0, 10 10, 0 10, 0 0)), ((20 0, 30 0, 30 10, 20 10, 20 0)))");
        Geometry mp2 = reader.read(
                "MULTIPOLYGON (((5 5, 15 5, 15 15, 5 15, 5 5)), ((30 10, 40 10, 40 20, 30 20, 30 10)))");

        Geometry loose = OverlayNG.overlay(mp, mp2, OverlayNG.INTERSECTION);
        System.out.println("non-strict => " + loose);
        System.out.println("            type=" + loose.getGeometryType()
                + " num=" + loose.getNumGeometries());

        OverlayNG ng = new OverlayNG(mp, mp2, null, OverlayNG.INTERSECTION);
        ng.setStrictMode(true);
        Geometry strict = ng.getResult();
        System.out.println("strict     => " + strict);
        System.out.println("            type=" + strict.getGeometryType()
                + " num=" + strict.getNumGeometries());
    }
}

输出：

non-strict => GEOMETRYCOLLECTION (POLYGON ((5 10, 10 10, 10 5, 5 5, 5 10)), POINT (30 10))
            type=GeometryCollection num=2
strict     => POLYGON ((5 10, 10 10, 10 5, 5 5, 5 10))
            type=Polygon num=1

注意严格模式并不是“扔掉所有低维元素”那么简单：如果交集中只有线、没有面（例如两条线交叉），严格模式下的交集结果会是空的 LineString，而不是保留那个交点。做“线网相交检测”时不要用叠加结果去判断，应该用第 08 章的 RelateNG 或 IntersectionMatrix。

9.8.4 与旧引擎语义的对照

如果不方便构造 OverlayNG 实例，也可以全局切换：-Djts.overlay=ng 之后，Geometry 的高层方法会走 OverlayNGRobust，但它仍然是非严格模式。严格模式只能通过显式构造 OverlayNG 并调用 setStrictMode(true) 获得，没有全局开关。

9.9 一元并集家族

9.9.1 为什么需要专门的一元并集

把 1000 个相邻宗地合并成一个街道，最直接的做法是循环调用 a = a.union(b)。这个朴素实现有两个致命问题：第一，每次二元并集都会重建整张拓扑图，总代价是 O(n) 次图构建，远高于必要的量级；第二，中间结果越来越复杂，后续每次并集的输入规模都在膨胀，而实际上大部分内部边在最终结果里会被溶解掉。

JTS 为此提供了“一元并集”这一族 API，核心思路是先收集所有边界边，一次性做 noding，再统一成面，避免反复构建中间几何。几个类的分工如下：

9.9.2 UnaryUnionOp：默认入口

UnaryUnionOp.union(geom) 与 Geometry.union() 等价，是“把所有东西并起来”的默认选择。它的内部流程分三步：先用 InputExtracter 把输入按维度拆成面、线、点三组；面组交给 CascadedPolygonUnion，线组走 noding 加 LineMerger 合并，点组做去重；最后按维度从高到低拼装结果。

它的默认并集函数是 CascadedPolygonUnion.CLASSIC_UNION，也就是旧引擎。这一点在源码里写得很清楚：字段声明为 private UnionStrategy unionFunction = CascadedPolygonUnion.CLASSIC_UNION;，而 CLASSIC_UNION 的实现就是调 OverlayOp.overlayOp(g0, g1, OverlayOp.UNION)。所以，Geometry.union() 默认也走旧引擎。若要让一元并集走新引擎，有两个办法：显式调 OverlayNGRobust.union(geom)，或者设置 -Djts.overlay=ng（此时 GeometryOverlay.union(a) 会转发到 OverlayNGRobust.union）。

9.9.3 级联并集的空间索引思想

CascadedPolygonUnion 的“级联”不是简单的两两迭代，而是把输入按空间邻近性组织：每次从集合里取一个多边形，在 STRtree 里查与它相交的候选，把候选取出来先合并成一个小结果，再放回集合继续。这样每次合并只处理“真正相邻”的几何，避免了大范围无效计算。

 输入：N 个多边形
        │
        ▼
 ┌────────────────────────────┐
 │ STRtree（node capacity = 4）│
 └─────────────┬──────────────┘
               │ 取一个多边形 p
               ▼
 ┌────────────────────────────┐
 │ 查询与 p 包络相交的候选集合  │
 └─────────────┬──────────────┘
               │ 候选逐个与 p 合并
               ▼
 ┌────────────────────────────┐
 │ 合并结果放回集合，重复直到 1 │
 └─────────────┬──────────────┘
               ▼
            最终并集

这段逻辑在源码里的常量是 STRTREE_NODE_CAPACITY = 4，与第 13 章要讲的 STRtree 默认值一致。级联策略对小规模数据（几十个面）优势不明显，但在上万面的批量溶解中通常能带来数倍加速。

9.9.4 四者结果一致、性能不同

下面这段代码用同一组 400 个互相重叠的小方块，对比四个入口的结果面积与耗时。四个入口都给出面积 78400（因为方块按 1.4 倍间距铺排，重叠部分被正确溶解），差别只在耗时：

import java.util.ArrayList;
import java.util.List;

import org.locationtech.jts.geom.Envelope;
import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.geom.GeometryFactory;
import org.locationtech.jts.geom.PrecisionModel;
import org.locationtech.jts.operation.overlayng.OverlayNGRobust;
import org.locationtech.jts.operation.overlayng.UnaryUnionNG;
import org.locationtech.jts.operation.union.CascadedPolygonUnion;
import org.locationtech.jts.operation.union.UnaryUnionOp;

public class UnaryUnionCompare {
    static final GeometryFactory GF = new GeometryFactory();

    static List<Geometry> squares(int n, double step, double offset) {
        List<Geometry> out = new ArrayList<>();
        int side = (int) Math.ceil(Math.sqrt(n));
        for (int i = 0; i < n; i++) {
            double x = (i % side) * step + offset;
            double y = (i / side) * step + offset;
            out.add(GF.toGeometry(new Envelope(x, x + step * 1.4, y, y + step * 1.4)));
        }
        return out;
    }

    public static void main(String[] args) {
        List<Geometry> geoms = squares(400, 10, 0);

        time("UnaryUnionOp        ", () -> UnaryUnionOp.union(geoms));
        time("CascadedPolygonUnion", () -> CascadedPolygonUnion.union(geoms));
        time("OverlayNGRobust     ", () -> OverlayNGRobust.union(geoms));
        time("UnaryUnionNG pm=1   ", () -> UnaryUnionNG.union(geoms, new PrecisionModel(1.0)));
    }

    static void time(String label, java.util.function.Supplier<Geometry> op) {
        op.get(); // 预热，避免把 JIT 时间算进去
        long t0 = System.nanoTime();
        Geometry g = op.get();
        long ms = (System.nanoTime() - t0) / 1_000_000;
        System.out.printf("%s => area=%.0f num=%d  %d ms%n",
                label, g.getArea(), g.getNumGeometries(), ms);
    }
}

一台普通笔记本上的实测（耗时为第二次运行，仅供量级参考）：

UnaryUnionOp         => area=78400 num=1  14 ms
CascadedPolygonUnion => area=78400 num=1  16 ms
OverlayNGRobust      => area=78400 num=1  19 ms
UnaryUnionNG pm=1    => area=78400 num=1  14 ms

结论与直觉略有出入：在这个规模下四个入口相差不大，新引擎甚至略慢（因为它多了校验与重试框架的开销）。级联策略的收益要到大几千个面才明显。因此选型时不要“唯新是从”，而应按数据规模与脏净程度决定：小数据用默认入口，大数据优先试 CascadedPolygonUnion 与 UnaryUnionNG，脏数据用 OverlayNGRobust。

补充：1.18 的版本历史里有一条 “Improve performance of UnaryUnionOp by removing OverlayUnion optimization (#644)”，说明早期 UnaryUnionOp 里曾有一条更激进的优化路径被移除，因为带来的边界情况比收益多。读老版本源码时不要被那部分代码误导。

9.10 CoverageUnion：覆盖面专用溶解

9.10.1 前提条件与收益

CoverageUnion（1.20 新增）针对的是“多边形覆盖”（polygonal coverage）这种特殊数据：一组多边形，彼此不重叠、不留缝隙，拼起来正好铺满一片区域。行政区划、地籍图、土地利用图都属于这一类。

对这类数据，溶解内部边其实不需要完整的叠加计算——只要把所有多边形的边界边收集起来，去掉那些“被两个多边形共享”的边，剩下的边围成的面就是并集。这正是 CoverageUnion 的做法：它用一个 BoundaryChainNoder 把边界链切开，按边的共享关系去重，再用 OverlayNG 的严格模式拼面。省掉了 noding 求交与标签推导两大开销，因此比通用并集快得多。

它的入口非常简洁，只有一个条带参数的静态方法：

// 签名：public static Geometry union(Geometry coverage)
Geometry dissolved = CoverageUnion.union(coverageMultiPolygon);

9.10.2 用错前提的后果

致命的一点是：CoverageUnion 不校验输入是否为合法覆盖。如果传进去的多边形互相重叠，它不会报错，而是按“边的共享次数”这个纯拓扑规则硬算，结果可能把重叠区域重复计入。下面这段代码就是反例：

import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.io.WKTReader;
import org.locationtech.jts.operation.overlayng.CoverageUnion;

public class CoverageUnionPitfall {
    public static void main(String[] args) throws Exception {
        WKTReader reader = new WKTReader();
        // 两个重叠的多边形：不是合法覆盖
        Geometry badCoverage = reader.read(
                "MULTIPOLYGON (((0 0, 10 0, 10 10, 0 10, 0 0)), ((5 5, 15 5, 15 15, 5 15, 5 5)))");

        Geometry result = CoverageUnion.union(badCoverage);
        System.out.println("result = " + result);
        System.out.println("num    = " + result.getNumGeometries());
        System.out.println("area   = " + result.getArea());
        System.out.println("expected union area = 175");
    }
}

输出（注意这是错误结果，但没有任何异常）：

result = MULTIPOLYGON (((0 0, 0 10, 10 10, 10 0, 0 0)), ((5 5, 5 15, 15 15, 15 5, 5 5)))
num    = 2
area   = 200.0
expected union area = 175

两个多边形只是被原样保留了下来，面积 200 而不是正确的 175（重叠的 25 被重复计算）。所以工程上必须遵守一条硬规则：调用 CoverageUnion 之前，先用 CoverageValidator.isValid(coverage) 校验。校验通过再溶解；校验失败则改用 UnaryUnionOp（慢但正确），或先走第 17 章的清洗流程。

9.10.3 与通用并集的对比

9.10.4 组合使用的最佳实践

覆盖面的完整处理流程通常是：先用 CoverageValidator.validate 找出问题，用 CoverageCleaner（1.21 新增）或 GeometryFixer 修掉，再用 CoverageUnion 溶解，最后对结果做一次 PrecisionReducer 统一精度。第 17 章会展开质检细节，这里给一段骨架代码：

import java.util.Arrays;

import org.locationtech.jts.coverage.CoverageCleaner;
import org.locationtech.jts.coverage.CoverageValidator;
import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.geom.MultiPolygon;
import org.locationtech.jts.geom.Polygon;
import org.locationtech.jts.operation.overlayng.CoverageUnion;

public class CoveragePipeline {
    public static void main(String[] args) {
        Geometry[] polygons = loadCoverage();

        if (!CoverageValidator.isValid(polygons)) {
            System.out.println("覆盖不合法，先清洗");
            // 第二个参数为负表示由实现自动推导吸附容差，第三个参数为可填补的最大缝隙宽度
            polygons = CoverageCleaner.clean(polygons, -1, 1.5);
        }

        MultiPolygon coverage = polygons[0].getFactory()
                .createMultiPolygon(Arrays.stream(polygons)
                        .map(g -> (Polygon) g)
                        .toArray(Polygon[]::new));

        Geometry dissolved = CoverageUnion.union(coverage);
        System.out.println("dissolved area = " + dissolved.getArea());
        System.out.println("dissolved num  = " + dissolved.getNumGeometries());
    }

    static Geometry[] loadCoverage() {
        throw new UnsupportedOperationException("替换为真实读取逻辑");
    }
}

9.11 TopologyException 的成因与容差策略

9.11.1 三类异常消息

叠加分析抛出的异常几乎全部是 TopologyException（继承自 RuntimeException），消息可以归成三类，每一类指向不同的根因：

第一类最常见，也最容易误判成“JTS 的 bug”。实际上它几乎总是“输入本身不可靠”或“浮点精度不足”导致的。第二类更严重：它通常意味着输入几何无效——例如两个多边形互相重叠（自相交的多边形集合本身就违反“面内部不相交”的定义），这时无论怎么重试都不会得到正确结果。第三类是保护性检查，它宁可抛异常也不返回错误几何，遇到时直接降精度重算即可。

9.11.2 一次完整的“踩坑到解决”过程

下面这段代码构造了一对“大量共享边且互相重叠”的方块集合，它们在旧引擎与浮点新引擎下都会抛出 side location conflict。这不是构造的假想，而是本章基准程序里真实发生的一次失败：

import java.util.ArrayList;
import java.util.List;

import org.locationtech.jts.geom.Envelope;
import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.geom.GeometryFactory;
import org.locationtech.jts.geom.PrecisionModel;
import org.locationtech.jts.operation.overlay.OverlayOp;
import org.locationtech.jts.operation.overlay.snap.SnapIfNeededOverlayOp;
import org.locationtech.jts.operation.overlayng.OverlayNG;

public class OverlayFailureDemo {
    static final GeometryFactory GF = new GeometryFactory();

    static Geometry grid(int n, double step, double size, double offset) {
        List<Geometry> parts = new ArrayList<>();
        int side = (int) Math.ceil(Math.sqrt(n));
        for (int i = 0; i < n; i++) {
            double x = (i % side) * step + offset;
            double y = (i / side) * step + offset;
            parts.add(GF.toGeometry(new Envelope(x, x + size, y, y + size)));
        }
        return GF.buildGeometry(parts);
    }

    public static void main(String[] args) {
        Geometry a = grid(400, 10, 14, 0);
        Geometry b = grid(400, 10, 14, 5);

        attempt("old OverlayOp   ", () -> OverlayOp.overlayOp(a, b, OverlayOp.INTERSECTION));
        attempt("SnapIfNeeded    ", () -> SnapIfNeededOverlayOp.intersection(a, b));
        attempt("OverlayNG float ", () -> OverlayNG.overlay(a, b, OverlayNG.INTERSECTION));
        attempt("OverlayNG fixed1", () ->
                OverlayNG.overlay(a, b, OverlayNG.INTERSECTION, new PrecisionModel(1.0)));
    }

    static void attempt(String label, java.util.function.Supplier<Geometry> op) {
        try {
            Geometry g = op.get();
            System.out.println(label + " => OK area=" + g.getArea()
                    + " num=" + g.getNumGeometries());
        } catch (RuntimeException ex) {
            System.out.println(label + " => FAILED " + ex.getMessage());
        }
    }
}

输出（本机实测）：

old OverlayOp    => FAILED side location conflict [ (0.0, 10.0, NaN) ]
SnapIfNeeded     => FAILED side location conflict [ (0.0, 10.0, NaN) ]
OverlayNG float  => FAILED side location conflict: arg 0 [ (50.0, 114.0, NaN) ]
OverlayNG fixed1 => OK area=39601.0 num=1

只有固定精度活了下来。原因是这些方块大量共享边、且互相重叠，浮点 noding 会在共享边的端点处产生方向与位置矛盾，而 Snap-Rounding 把所有顶点先吸到整数网格上，让共享边真正重合，矛盾随之消失。这是前文“固定精度更稳”的最直接证据。

9.11.3 容差策略的四条经验

结合源码里的容差计算与上面的实验，可以总结出四条可操作的经验：

第一条，先判断是数据问题还是算法问题。用 geom.isValid() 与 IsValidOp 检查输入；如果输入本身重叠或自相交，先去修数据（第 17 章），不要指望换引擎或加容差。

第二条，对“本应重合但不完全重合”的边，主动吸附。可以用 GeometrySnapper.snapToSelf(geom, tol, true) 对每个输入自吸附，容差取业务精度的 1 到 10 倍；也可以直接用固定精度叠加，让引擎替你吸。

第三条，容差要与坐标量级匹配。GeometrySnapper 的默认容差是包络短边乘 1e-9，OverlayNGRobust 的初始吸附容差是最大 ordinate 乘 1e-12；这两个值在 1e2 量级的数据上分别是 1e-7 与 1e-10，都很小。如果数据是经纬度（量级 1e2），1e-10 度约等于 1e-5 米，勉强够用；如果是投影坐标（量级 1e7 米），1e-12 乘 1e7 等于 1e-5 米，也够用。但如果业务精度是厘米级而输入噪声是米级，这些容差完全不够，必须显式指定固定精度。

第四条，捕获异常后要有降级路径。生产代码建议按“OverlayNGRobust → 固定精度重算 → 记录并跳过”三级降级，而不是让一个要素的拓扑异常拖垮整个批处理：

import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.geom.PrecisionModel;
import org.locationtech.jts.operation.overlayng.OverlayNG;
import org.locationtech.jts.operation.overlayng.OverlayNGRobust;
import org.locationtech.jts.operation.overlayng.PrecisionUtil;

public class OverlayWithFallback {
    public static Geometry safeIntersection(Geometry a, Geometry b) {
        try {
            return OverlayNGRobust.overlay(a, b, OverlayNG.INTERSECTION);
        } catch (RuntimeException first) {
            double scale = PrecisionUtil.safeScale(a, b);
            try {
                return OverlayNG.overlay(a, b, OverlayNG.INTERSECTION,
                        new PrecisionModel(scale));
            } catch (RuntimeException second) {
                throw new IllegalStateException(
                        "叠加失败: " + first.getMessage() + " / " + second.getMessage(), first);
            }
        }
    }
}

9.12 新旧引擎的结果与性能对比

9.12.1 结果层面的差异

对合法、无近重合边的输入，两套引擎的结果在拓扑上等价，差异只体现在三处：顶点坐标（新引擎在浮点模式下通常保留更多有效位）、环的起点与方向、MultiPolygon 内部元素顺序。用 norm() 规范化后再比较，通常可以判定为“完全相同”。

对“接缝有微小偏差”的输入，两者会给出结构不同的结果。9.7.4 的例子已经展示过：旧引擎把 1e-9 的缝隙当成真实缝隙，返回两个多边形；新引擎把它当作重合边，返回一个多边形。这不是谁对谁错，而是精度假设不同——旧引擎的吸附只在抛异常之后才发生（SnapIfNeededOverlayOp），而新引擎在固定精度下的 Snap-Rounding 是默认行为。

对“病态”输入，9.11.2 给出了更极端的情况：旧引擎与新引擎的浮点模式都失败，只有固定精度成功。

9.12.2 性能层面的差异

性能是选型的另一个维度。9.9.4 的数据显示：在小规模、干净的输入上，旧引擎往往更快，因为它的数据结构更紧凑、没有校验与重试框架；在需要大量 noding 的大输入上，新引擎的裁剪优化与更轻的边结构开始占优。

用锯齿多边形做基准（两个各含 8005 个顶点的“梳齿”多边形，互相错位 0.5 重叠）可以得到一组有代表性的对比。构造方法如下：

static Geometry comb(int teeth, double toothW, double toothH, double offset) {
    List<Coordinate> cs = new ArrayList<>();
    cs.add(new Coordinate(0, 0));
    for (int i = 0; i < teeth; i++) {
        double x0 = offset + i * toothW * 2;
        cs.add(new Coordinate(x0, 0));
        cs.add(new Coordinate(x0, toothH));
        cs.add(new Coordinate(x0 + toothW, toothH));
        cs.add(new Coordinate(x0 + toothW, 0));
    }
    cs.add(new Coordinate(offset + teeth * toothW * 2, 0));
    cs.add(new Coordinate(offset + teeth * toothW * 2, -5));
    cs.add(new Coordinate(0, -5));
    cs.add(new Coordinate(0, 0));
    return new GeometryFactory().createPolygon(cs.toArray(new Coordinate[0]));
}

本机实测（每个操作先预热一次，取第二次耗时；绝对值随机器变化，只比较同一台机器上的相对关系）：

三点解读。第一，浮点新引擎在中大型输入上比旧引擎快 30% 到 40%，收益来自更轻的边结构与裁剪优化。第二，固定精度慢得多（交集 488 ms），因为 Snap-Rounding 要把每个交点吸附到网格并重建拓扑；但它换来的是“绝不失败”，这笔账是否划算取决于业务。第三，OverlayNGRobust 的交集只要 36 ms，比纯浮点路径更快——因为它内部第一级就用浮点算成功了，而它的调用路径上少了一层旧引擎的平面图结构开销。

重要提醒：固定精度下若 scale 选得过大（网格过粗），薄特征会塌缩。上例中 0.5 宽的梳齿在 scale=1.0（整米网格）下全部被吸成零面积，交集结果从 1 个面积 70000 的多边形变成 4000 个面积 20000 的碎片。选 scale 时必须保证“最窄业务特征远大于 1/scale”。

9.12.3 决策清单

| 你的情况 | 建议 | | — | — | | 数据来历不明、批量处理、不能中断 | OverlayNGRobust.overlay 加固定精度兜底 | | 需要结果可复现、可入库比对 | 固定精度 OverlayNG，scale 取业务精度的 10 倍 | | 单次小规模计算、追求低延迟 | 默认 Geometry.union（旧引擎）即可 | | 合法覆盖面批量溶解 | CoverageUnion（务必先跑 CoverageValidator） | | 大量多边形溶解 | CascadedPolygonUnion 或 UnaryUnionNG | | 需要混合维度结果做诊断 | OverlayNG 非严格模式 | | 下游算子要求同质输入 | OverlayNG 加 setStrictMode(true) |

9.13 常见问题与避坑

问题 1：两个多边形明明相交，intersection 却返回 POLYGON EMPTY。 先确认两个几何的包络是否相交（getEnvelopeInternal().intersects(other)）。若包络相交但结果为空，多半是坐标系不一致（一个经纬度、一个投影坐标），导致两者在数值上相距极远。JTS 不含 CRS 概念，这类错误不会被引擎发现，必须由上层校验。

问题 2：结果面积比两个输入之和还大。 典型的重复计算，常见于三种情况：对互相重叠的多边形用 CoverageUnion（见 9.10.2）；对同一个几何做了两次并集且中间结果未被正确去重；用了 GeometryCollection 输入而交集被逐元素映射后简单拼接。解决方法是改用 UnaryUnionOp.union，并对结果做 isValid() 检查。

问题 3：union 报 Operation does not support GeometryCollection arguments。 二元并、差、对称差不接受非空 GeometryCollection。按维度拆分输入（PolygonExtracter、LineStringExtracter、PointExtracter），或改用一元并集 geom.union()。

问题 4：同样的代码，本地结果和服务器结果不一致。 两种可能：一是全局属性不同（服务器上设了 -Djts.overlay=ng），二是浮点运算在不同 JIT 与 CPU 上的差异。要彻底消除，必须用固定精度；只靠浮点无法保证跨平台逐位一致。

问题 5：TopologyException 捕获后重试仍然失败。 先看消息类型。side location conflict 与 non-noded intersection 若在固定精度下依旧出现，说明输入几何本身无效（自相交、环方向错误、坐标含 NaN），应先跑 GeometryFixer.fix。IsValidOp 能定位到具体的问题点，比盲目重试高效得多。

问题 6：叠加之后顶点暴增，内存吃紧。 交集与差集在裁剪优化下只处理相关区域，但并集必须处理全部边。若结果顶点过多，可对结果做 PrecisionReducer.reducePrecision 降精度，或用第 12 章的 DouglasPeuckerSimplifier 化简（注意化简必须在叠加之后，且要考虑拓扑保持）。

问题 7：Z 值在结果里丢了。 WKTWriter 默认输出二维，看起来像丢了；实际上 OverlayNG 的 ElevationModel 会沿边线性插值 Z。确认方式是用 new WKTWriter(3) 输出，或在代码里读 getCoordinate().getZ()。旧引擎不插值，只保留来源坐标，所以来回切换引擎会看到 Z 值差异。

问题 8：difference 之后公共边消失了。 这是 OGC 语义决定的，不是 bug（见 9.1.4）。需要保留界址线时，应在线图层单独求交：boundaryA.intersection(boundaryB)。

问题 9：symDifference 的结果总是 MultiPolygon，无法直接入库。 对称差天然可能产生多个不相连的块，返回 MultiPolygon 是正确行为。若业务表只接受 Polygon，需要按面积排序取主块，或改用两次 difference 分别处理。

问题 10：大量小多边形求并集时 OOM。 改用 CascadedPolygonUnion 分批处理；或先把输入按网格分块（STRtree 查询），块内溶解后再块间溶解。不要一次性把上百万个多边形塞进 UnaryUnionOp。

9.14 本章小结

本章从数学语义一路讲到源码实现与工程决策，要点可以浓缩成十条：

1. 四个叠加运算在数学上是集合的交、并、差、对称差，JTS 用 Geometry 的四个方法对外暴露，返回值遵循“最小充分类型”，可能是 Polygon、MultiPolygon 或 GeometryCollection。
2. 结果维度可能低于输入维度（拓扑塌缩），这是语义的一部分而不是异常；需要同质结果时用 OverlayNG.setStrictMode(true)。
3. 代数恒等式在浮点坐标下只是近似成立；断言结果时应比较面积、周长、顶点集合等不变量，而不是 WKT 字符串。
4. 空几何与 GeometryCollection 有专门的短路与限制规则：交集支持 GC（逐元素映射），二元并、差、对称差不支持。
5. JTS 有两套引擎：旧 OverlayOp 与新 OverlayNG；高层方法默认仍走旧引擎，开关是系统属性 jts.overlay，取值为 old 或 ng。
6. OverlayNG 的流水线是“noding → 建图 → 标注 → 按维度构结果 → 插值 Z”，其中 noding 的精度模型决定鲁棒性，标注阶段决定结果语义。
7. OverlayLabel 用边界边、塌缩边、线边、非本几何边四种状态描述一条边的身份，这个设计使“重合边合并”与“塌缩边丢弃”成为统一规则的特例。
8. 精度是可控的：PrecisionUtil 给出安全比例尺与固有比例尺，PrecisionReducer 负责降精度且保证输出合法；固定精度更稳但更慢，且可能让薄特征塌缩。
9. OverlayNGRobust 用“浮点 → 五轮递增容差吸附 → Snap-Rounding”三级阶梯兜底，仍失败则抛出最初的异常；它救不了非法输入与混合维度输入。
10. 一元并集是一族工具而非一个方法：UnaryUnionOp 通用、CascadedPolygonUnion 适合大批量、CoverageUnion 只对合法覆盖正确且最快；用错前提会静默给出错误结果。

下一章进入缓冲区与偏移曲线：叠加回答“交集在哪”，缓冲区回答“距离它多远的范围内”。两者在实现上共享 noding 与构面技术，但在参数化与误差控制上有完全不同的关注点。

相关资源

• JTS 官方 Javadoc：https://locationtech.github.io/jts/javadoc/ ，重点读 OverlayNG、OverlayNGRobust、PrecisionUtil、CoverageUnion 的类注释，本章很多结论直接从这些注释与源码得来。
• 源码位置：modules/core/src/main/java/org/locationtech/jts/operation/overlayng/（新引擎）、同目录下的 operation/overlay/（旧引擎）、operation/union/（一元并集）。
• 版本历史：doc/JTS_Version_History.md。OverlayNG 的引入见 1.18.0 条目（#599、#615、#645）；CoverageUnion 见 1.20.0 条目（#900）；OverlayNG 面积启发式见 1.18.1 条目（#812）。
• 本站第 08 章《空间关系谓词与 DE-9IM》：https://znlgis.github.io/gis/tutorial/jts/第08章-空间关系谓词与DE-9IM/
• 本站第 10 章《缓冲区与偏移曲线》：https://znlgis.github.io/gis/tutorial/jts/第10章-缓冲区与偏移曲线/
• 本站第 17 章《有效性与 Coverage 质检》：https://znlgis.github.io/gis/tutorial/jts/第17章-有效性与Coverage质检/
• 本站第 18 章《精度模型与数值鲁棒性》：https://znlgis.github.io/gis/tutorial/jts/第18章-精度模型与数值鲁棒性/