第11章 - 距离度量与几何统计

“这个点到那条路多远”“这两个多边形像不像”“这个湖最长有多宽”“这块地最适合放标签的位置在哪”——这些问题看起来都属于“算个数”,实际上分属四类完全不同的算法:点与几何的最小距离、两个形状之间的相似度(Hausdorff 与 Fréchet)、几何自身的尺度统计(面积、周长、质心、最小外接形状),以及几何内部的极值点(最大内切圆心、最大空圆心)。JTS 把这些能力分散在 algorithm、algorithm.distance、algorithm.construct、operation.distance、operation.distance3d、math 六个包里,接口风格各不相同,很容易用错。

本章按“从简单到复杂”的顺序组织:先讲最基础的量(长度、面积、质心、内点、边界维度),它们是 Geometry 上的直接方法,几乎没有坑;再讲距离家族,重点是 DistanceOp 与 IndexedFacetDistance 的语义差异(一个对多边形内部点返回 0,一个返回到边界的距离);然后是三个“形状相似度”度量:离散 Hausdorff、有向 Hausdorff(1.21 新增)与离散 Fréchet;接着是几何尺度统计工具:最小直径、最小外接圆、最小面积矩形、最小外接三角形(1.21 新增);最后是两个“在内部找极值点”的构造类:最大内切圆与最大空圆,它们在制图注记、设施选址里非常实用。

前置知识:第 04 章的几何对象模型、第 10 章的缓冲(内切圆半径就是“负缓冲变空的距离”)、第 13 章的空间索引(IndexedFacetDistance 内部用 STRtree,本章会提但不展开)。本章所有示例都给出真实运行输出,单位一律以“坐标单位”计——JTS 不含 CRS,这一点在第 19 章会专门展开。

11.1 基础统计量

11.1.1 面积与长度

Geometry 上有两个最常用的统计方法:getArea() 与 getLength()。它们对不同类型的几何有不同含义,这一点必须记清:

几何类型 getArea() getLength()
Point / MultiPoint 0 0
LineString / MultiLineString 0 各段长度之和
Polygon 外环面积减各孔面积 外环周长加各孔周长
MultiPolygon 各多边形面积之和 各多边形周长之和
GeometryCollection 各元素面积之和 各元素长度之和

两个容易误解的点。第一,Polygon.getArea() 是“净面积”,孔洞会被扣掉;想要“毛面积”要自己取外环。第二,Polygon.getLength() 不是外环周长,而是“所有环的周长之和”,带孔的湖算周长时会把孔也算进去——这在海岸线长度统计里是个真实存在的陷阱。

底层实现分别委托给 org.locationtech.jts.algorithm.Area 与 Length:Area.ofRing(Coordinate[]) 算环的绝对面积(鞋带公式),Area.ofRingSigned 保留符号(逆时针为正、顺时针为负),Length.ofLine(CoordinateSequence) 逐段累加。JTS 1.21 把长度计算改成用 MathUtil.hypot(#1112),比 Math.sqrt(dx*dx+dy*dy) 更快,但注意它并不做溢出保护,详见 11.9.3。

11.1.2 质心与内点

getCentroid() 返回几何的质心(重心)。它不是“包络中心”,也不保证落在几何内部——对 C 形或环形几何,质心很可能在几何外部。这一点在制图注记、点位聚合里经常被忽略。

getInteriorPoint() 返回一个保证落在几何内部的点。它的实现是三个专用类的组合:InteriorPointArea(面)、InteriorPointLine(线)、InteriorPointPoint(点)。面用“水平扫描线取最长扫描段的中点”,算法快且鲁棒,得到的是“最宽处的中心”,很适合放标签。

  ┌──────────────────────┐        ┌──────────────────────┐
  │  ╭────╮              │        │  ╭────╮              │
  │  │    │   centroid   │        │  │    │  interiorPt  │
  │  ╰────╯      *       │        │  ╰────╯        *     │
  │                      │        │                      │
  │   C 形几何:质心可能   │        │   内点保证在几何内部   │
  │   落在空洞或外部       │        │   (取最宽扫描段中点)  │
  └──────────────────────┘        └──────────────────────┘

11.1.3 边界维度与边界

getBoundaryDimension() 返回边界的拓扑维度:面的边界是线(1),线的边界是点(0),点的边界是空(-1)。它常用于“需要按边界维度分组”的批量处理,也用于判断“这个几何有没有边界”:

几何 getDimension() getBoundaryDimension() getBoundary() 返回
Point 0 -1 空的 GeometryCollection
LineString 1 0 MultiPoint(首末点)
LinearRing 1 -1 空(环无边界)
Polygon 2 1 MultiLineString(所有环)

LinearRing 的边界维度是 -1,这是 OGC 规范里“闭合曲线没有边界”的直接体现,也解释了为什么两个相邻多边形共享边的交集是线而不是点。

11.1.4 一个把统计量打全的示例

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

public class BasicStats {
    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 line = reader.read("LINESTRING (0 0, 5 5, 10 0)");
        Geometry pt = reader.read("POINT (3 4)");

        dump("polygon", poly);
        dump("line   ", line);
        dump("point  ", pt);

        System.out.println("poly area     = " + poly.getArea());
        System.out.println("line length   = " + line.getLength());
        System.out.println("poly centroid = " + poly.getCentroid());
        System.out.println("line centroid = " + line.getCentroid());
        System.out.println("poly interior = " + poly.getInteriorPoint());
        System.out.println("line interior = " + line.getInteriorPoint());
        System.out.println("poly boundary = " + poly.getBoundary());
        System.out.println("line boundary = " + line.getBoundary());
    }

    static void dump(String label, Geometry g) {
        System.out.printf("%s => dim=%d boundaryDim=%d numPts=%d%n",
                label, g.getDimension(), g.getBoundaryDimension(), g.getNumPoints());
    }
}

输出:

polygon => dim=2 boundaryDim=1 numPts=5
line    => dim=1 boundaryDim=0 numPts=3
point   => dim=0 boundaryDim=-1 numPts=1
poly area     = 100.0
line length   = 14.142135623730951
poly centroid = POINT (5 5)
line centroid = POINT (5 2.5)
poly interior = POINT (5 5)
line interior = POINT (5 5)
poly boundary = LINEARRING (0 0, 10 0, 10 10, 0 10, 0 0)
line boundary = MULTIPOINT ((0 0), (10 0))

折线 (0 0) → (5 5) → (10 0) 的质心是 (5 2.5),位于几何之外的“下方”,而内点是 (5 5),正好落在折线顶点上。这是两者语义差异的最直观展示。另外注意 poly.getBoundary() 返回的是 LINEARRING 而不是 MultiLineString,因为单个面无孔时边界就是一条闭合环——判断类型时要小心。

11.2 DistanceOp:最小距离

11.2.1 语义与返回值

DistanceOp 计算两个几何之间的最小距离,并提供三种查询方式:

方法 返回 说明
distance(Geometry, Geometry) double 最小距离值
isWithinDistance(g0, g1, d) boolean 是否在 d 之内,带提前终止
nearestPoints(g0, g1) Coordinate[2] 最近的两个点(旧名 closestPoints 已废弃)
nearestLocations()(实例方法) GeometryLocation[2] 带“属于哪个组件、第几段”的定位信息

语义上它是“集合间距离”:若两个几何有交集,距离为 0。对多边形,内部点的距离为 0——这与 IndexedFacetDistance 不同(后者算到边界 facet 的距离),是本章最容易搞混的一点。

import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.io.WKTReader;
import org.locationtech.jts.operation.distance.DistanceOp;
import org.locationtech.jts.operation.distance.GeometryLocation;

public class DistanceBasic {
    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 inside = reader.read("POINT (3 4)");
        Geometry outside = reader.read("POINT (13 4)");

        System.out.println("inside  distance = " + DistanceOp.distance(poly, inside));
        System.out.println("outside distance = " + DistanceOp.distance(poly, outside));
        System.out.println("within 2?  " + DistanceOp.isWithinDistance(poly, outside, 2.0));
        System.out.println("nearest  = " + java.util.Arrays.toString(
                DistanceOp.nearestPoints(poly, outside)));

        DistanceOp op = new DistanceOp(poly, outside);
        for (GeometryLocation loc : op.nearestLocations()) {
            System.out.println("component = " + loc.getGeometryComponent().getGeometryType()
                    + " segIndex = " + loc.getSegmentIndex()
                    + " coord = " + loc.getCoordinate());
        }
    }
}

输出:

inside  distance = 0.0
outside distance = 3.0
within 2?  false
nearest  = [(10.0, 4.0, NaN), (13.0, 4.0, NaN)]
component = POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0)) segIndex = -1 coord = (10.0, 4.0, NaN)
component = POINT (13 4) segIndex = 0 coord = (13.0, 4.0, NaN)

注意 segIndex = -1:它表示最近点不是某个线段上的插值点,而是组件内部的点(本例是位于多边形内部的查询点,距离为 0 时引擎直接把内部点当作“最近点”返回)。这个细节在调试“为什么最近点坐标对不上”时很有用。

11.2.2 实现与性能

DistanceOp 的内部结构分为三步:先用 ConnectedElementPointFilter 与 ConnectedElementLocationFilter 把两个几何拆成“点、线段”级别的组件;再对每一对组件用 DistanceToPoint / 线段-线段距离公式算最小距离;最后取全局最小。为了加速,它在两侧都建了 FacetSequenceTreeBuilder 生成的索引(本质是 STRtree 上的 facet 序列),对大规模输入能剪掉大部分分支。

一条实用的性能经验:isWithinDistance(g0, g1, d) 比 distance(g0, g1) <= d 快得多,因为它一旦找到小于 d 的组件对就可以立即返回,不必算完全部。在“判断是否在影响范围内”这类布尔查询里,永远优先用前者。

import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.io.WKTReader;
import org.locationtech.jts.operation.distance.DistanceOp;

public class DistanceEarlyExit {
    public static void main(String[] args) throws Exception {
        WKTReader reader = new WKTReader();
        Geometry a = reader.read("LINESTRING (0 0, 100 0)");
        Geometry b = reader.read("LINESTRING (0 0.5, 100 0.5)");

        System.out.println("distance   = " + DistanceOp.distance(a, b));
        System.out.println("within 1.0 = " + DistanceOp.isWithinDistance(a, b, 1.0));
        System.out.println("within 0.2 = " + DistanceOp.isWithinDistance(a, b, 0.2));
        System.out.println("nearest    = " + java.util.Arrays.toString(
                DistanceOp.nearestPoints(a, b)));
    }
}

输出:

distance   = 0.5
within 1.0 = true
within 0.2 = false
nearest    = [(0.0, 0.0, NaN), (0.0, 0.5, NaN)]

两条平行线距离恒为 0.5,所以 within 1.0 为 true、within 0.2 为 false。这个例子小到看不出性能差异,但语义差异是本质的:isWithinDistance 一旦找到一个“距离小于阈值”的组件对就可以立即返回,而 distance 必须遍历完所有组件对才能确定最小值。当两个几何各有上万个顶点时,前者的收益是数量级的——这也是它在 JTS 内部被大量用于“缓冲区预筛选”的原因。

11.2.3 空几何与维度混合

DistanceOp 对空几何返回什么?答案是 0.0。这看起来反直觉,但源码里明确处理了空输入:只要任一输入为空,就认为“没有点可以参与计算”,直接返回 0。因此不能用距离是否为 0 来判断两个几何是否相交——空几何的距离也是 0。

维度混合(点对面、线对面)则完全支持,语义仍是“集合间最小距离”,面的内部点距离为 0。这一点与 IndexedFacetDistance 的差异见下一节。

11.3 IndexedFacetDistance:带索引的距离查询

11.3.1 与 DistanceOp 的关键差异

IndexedFacetDistance 的类注释第一句就说清楚了它的定位:“Computes the distance to the facets (segments and vertices) of a Geometry using a Branch-and-Bound algorithm.”——注意 facet(线段与顶点),而不是“集合”。这带来一个关键差异:

场景 DistanceOp IndexedFacetDistance
点落在多边形内部 0.0 到最近边界的距离(正数)
点落在多边形外部 到边界距离 到边界距离
两条相交线 0.0 0.0
大量点查询同一几何 每次重建索引,慢 索引缓存,快
线程安全 是 是(类注释明确说明)

也就是说,IndexedFacetDistance 回答的是“离这条边界/这些顶点多远”,而不是“离这个面多远”。这个差异在做“点到道路中心线距离”时无所谓(线没有内部),但在“点到行政区边界距离”时就是本质区别。

import org.locationtech.jts.geom.Coordinate;
import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.io.WKTReader;
import org.locationtech.jts.operation.distance.DistanceOp;
import org.locationtech.jts.operation.distance.IndexedFacetDistance;

public class FacetVsSetDistance {
    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 inside = reader.read("POINT (3 4)");

        System.out.println("DistanceOp.distance      = " + DistanceOp.distance(poly, inside));
        System.out.println("IndexedFacetDistance.d   = " + IndexedFacetDistance.distance(poly, inside));
        System.out.println("IndexedFacetDistance.np  = "
                + java.util.Arrays.toString(IndexedFacetDistance.nearestPoints(poly, inside)));

        IndexedFacetDistance ifd = new IndexedFacetDistance(poly);
        System.out.println("nearestPoint (20,20)     = " + ifd.nearestPoint(new Coordinate(20, 20)));
        System.out.println("nearestLocation (20,20)  = "
                + ifd.nearestLocation(new Coordinate(20, 20)).getCoordinate());
        System.out.println("distance (20,20)         = " + ifd.distance(new Coordinate(20, 20)));
    }
}

输出:

DistanceOp.distance      = 0.0
IndexedFacetDistance.d   = 3.0
IndexedFacetDistance.np  = [(0.0, 4.0, NaN), (3.0, 4.0, NaN)]
nearestPoint (20,20)     = (10.0, 10.0, NaN)
nearestLocation (20,20)  = (10.0, 10.0, NaN)
distance (20,20)         = 14.142135623730951

同一个内部点 (3,4),两个类给出了 0.0 与 3.0 两个答案,两者都是“对的”——取决于你问的是哪种距离。这个例子值得贴在代码评审的注释里,因为它能解释很多“为什么两个模块算出的距离不一样”的争论。

11.3.2 缓存索引的正确姿势

IndexedFacetDistance 的真正价值在“一次建索引、多次查询”的场景。类注释给的典型用法是:对目标几何构造一次实例,然后反复调用 distance(Coordinate)、nearestPoints(Geometry)、isWithinDistance(Geometry, double)。

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

import org.locationtech.jts.geom.Coordinate;
import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.geom.GeometryFactory;
import org.locationtech.jts.io.WKTReader;
import org.locationtech.jts.operation.distance.IndexedFacetDistance;

public class BatchDistanceQuery {
    public static void main(String[] args) throws Exception {
        WKTReader reader = new WKTReader();
        GeometryFactory gf = new GeometryFactory();
        Geometry road = reader.read("LINESTRING (0 0, 50 10, 100 0, 150 30)");

        // 只建一次索引
        IndexedFacetDistance toRoad = new IndexedFacetDistance(road);

        List<Coordinate> sites = new ArrayList<>();
        for (int i = 0; i < 5; i++) {
            sites.add(new Coordinate(i * 30, 40));
        }

        for (Coordinate site : sites) {
            double d = toRoad.distance(site);
            boolean inCorridor = toRoad.isWithinDistance(gf.createPoint(site), 20.0);
            System.out.printf("site (%.0f, %.0f)  d=%.3f  within20=%b%n",
                    site.x, site.y, d, inCorridor);
        }
    }
}

输出:

site (0, 40)  d=39.223  within20=false
site (30, 40)  d=33.340  within20=false
site (60, 40)  d=31.379  within20=false
site (90, 40)  d=37.262  within20=false
site (120, 40)  d=24.010  within20=false

五个站点都离这条路超过 20 米(最近的是 (120,40),24.01 米),所以 within20 全为 false——这条折线在 y 方向最高只到 30,而站点在 y=40 一排,距离自然都在 20 以上。如果换成 DistanceOp.distance(road, point) 循环调用,每次都会重建两侧的 facet 索引,几百个点以上就会明显变慢。批量查询务必用 IndexedFacetDistance。

11.3.3 与 PreparedGeometry 的关系

org.locationtech.jts.geom.prep 包里的 PreparedGeometry 是另一套“预计算”机制,方向相反:它预计算的是目标几何,用于快速回答“某几何是否与它相交/包含/覆盖”。两者不是替代关系:

  • 问“距离多少”用 IndexedFacetDistance;
  • 问“是否相交/包含/覆盖”用 PreparedGeometry;
  • 问“是否在某距离内”两者都能回答:IndexedFacetDistance.isWithinDistance 或 PreparedGeometry.intersects(buffer),但前者精确且不需要构造缓冲区。

11.4 离散 Hausdorff 距离

11.4.1 定义与“离散”的含义

Hausdorff 距离衡量两个点集“形状相差多少”:它取“一个集合中任意点到另一个集合的最近距离”的最大值。数学定义有两层:

  有向 Hausdorff:  h(A, B) = max over a in A of ( min over b in B of dist(a, b) )
  Hausdorff 距离:   H(A, B) = max( h(A, B), h(B, A) )

JTS 的 DiscreteHausdorffDistance 实现的是“离散”版本:只在顶点上取极值,而且默认不做线段插值。这带来两个后果:一是计算快,二是结果依赖顶点密度——两个形状相同但顶点分布不同的几何,离散 Hausdorff 距离会不同。

  A: ●───────●───────●          顶点稀疏
  B: ●─●─●─●─●─●─●─●─●          顶点密集
  两者的真实 Hausdorff 距离约等于 0,
  但离散版本可能因为顶点错位给出非零值。

11.4.2 densifyFrac 参数

为了缓解顶点依赖,DiscreteHausdorffDistance 提供了 densifyFrac 参数:把每条线段按比例细分,在细分点上取样。distance(g0, g1, densifyFrac) 的第三个参数就是这个比例,取值 (0, 1],越小越密。

`java import java.util.Arrays;

import org.locationtech.jts.algorithm.distance.DiscreteHausdorffDistance; import org.locationtech.jts.geom.Geometry; import org.locationtech.jts.io.WKTReader;

public class HausdorffDemo { public static void main(String[] args) throws Exception { WKTReader reader = new WKTReader(); Geometry flat = reader.read(“LINESTRING (0 0, 10 0)”); Geometry bump = reader.read(“LINESTRING (0 0, 5 5, 10 0)”); Geometry offset = reader.read(“LINESTRING (0 1, 10 1)”);

    System.out.println("flat vs offset    = "
            + DiscreteHausdorffDistance.distance(flat, offset));
    System.out.println("flat vs bump      = "
            + DiscreteHausdorffDistance.distance(flat, bump));
    System.out.println("flat vs bump 0.01 = "
            + DiscreteHausdorffDistance.distance(flat, bump, 0.01));
    System.out.println("distanceLine      = "
            + DiscreteHausdorffDistance.distanceLine(flat, bump));
    System.out.println("oriented(flat,bump) = "
            + DiscreteHausdorffDistance.orientedDistance(flat, bump));
    System.out.println("oriented(bump,flat) = "
            + DiscreteHausdorffDistance.orientedDistance(bump, flat));
    System.out.println("coords = "
            + Arrays.toString(new DiscreteHausdorffDistance(flat, bump).getCoordinates()));
} } `

输出:

` flat vs offset = 1.0 flat vs bump = 3.5355339059327378 flat vs bump 0.01 = 3.5355339059327378 distanceLine = LINESTRING (5 0, 5 5) oriented(flat,bump) = 3.5355339059327378 oriented(bump,flat) = 5.0 coords = [(5.0, 0.0, NaN), (2.5, 2.5, NaN)] `

这组数据揭示了一个重要事实:离散 Hausdorff 距离与“最远点对”不是一回事。flat 与 offset 是两条平行线,最近距离恒为 1,所以结果是 1.0。flat 与 bump 的结果是 3.5355,因为 bump 的顶点 (5,5) 到 flat 的最近点距离是 5,但 flat 上 (0,0) 与 (10,0) 到 bump 的最近距离更小;取“所有点最近距离的最大值”得到 3.5355,即从 flat 的中点 (5,0) 到 bump 顶点 (5,5) 的斜边方向的 5 除以 √2。distanceLine 返回的正是这对点 (5 0, 5 5)。

oriented 两个方向给出 3.54 与 5.0 两个不同值,再次印证有向 Hausdorff 的不对称性:从 bump 的每个点到 flat 的最近距离最大是 5(顶点到地面),从 flat 的每个点到 bump 的最近距离最大是 3.54。要得到对称的 Hausdorff 距离,必须取两者最大值,也就是用 DirectedHausdorffDistance.hausdorffDistance(g0, g1)。

11.5 有向 Hausdorff 距离(1.21 新增)

11.5.1 连续版与离散版的区别

DirectedHausdorffDistance 是 1.21 新增的类(版本历史 #1182),它把“有向 Hausdorff 距离”从离散顶点扩展到了连续几何。类注释给的定义是:有向 Hausdorff 距离是“查询几何 A 上的任意一点到目标几何 B 的最大距离”,等价说法是“A 上的所有点都在距 B 不超过该距离的范围内”。

这个定义与 DiscreteHausdorffDistance 的差异在于取样位置:离散版只看顶点,连续版会在边上取极值点(用二分逼近)。源码注释直接给出了对比结论:“This algorithm is easier to use, more accurate, and much faster than DiscreteHausdorffDistance.”

维度 DiscreteHausdorffDistance DirectedHausdorffDistance
极值位置 仅顶点 连续,可能在边内部
对称性 提供 orientedDistance 与 distance 两个入口 本身有向,对称版本是 hausdorffDistance
精度控制 densifyFrac 自动容差,或显式传 tolerance
可复用 每次新建 可在“目标几何”上实例化并复用索引
新增版本 1.11 时代 1.21

11.5.2 三组核心 API

import java.util.Arrays;

import org.locationtech.jts.algorithm.distance.DirectedHausdorffDistance;
import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.io.WKTReader;

public class DirectedHausdorffDemo {
    public static void main(String[] args) throws Exception {
        WKTReader reader = new WKTReader();
        Geometry flat = reader.read("LINESTRING (0 0, 10 0)");
        Geometry bump = reader.read("LINESTRING (0 0, 5 5, 10 0)");

        System.out.println("DHD(flat -> bump) = "
                + DirectedHausdorffDistance.distance(flat, bump));
        System.out.println("DHD(bump -> flat) = "
                + DirectedHausdorffDistance.distance(bump, flat));
        System.out.println("HD  = "
                + DirectedHausdorffDistance.hausdorffDistance(flat, bump));
        System.out.println("points = "
                + Arrays.toString(DirectedHausdorffDistance.distancePoints(flat, bump)));
        System.out.println("fullyWithin 4 = "
                + DirectedHausdorffDistance.isFullyWithinDistance(flat, bump, 4.0));
        System.out.println("fullyWithin 3 = "
                + DirectedHausdorffDistance.isFullyWithinDistance(flat, bump, 3.0));
    }
}

输出:

DHD(flat -> bump) = 3.5355339059327378
DHD(bump -> flat) = 5.0
HD  = 5.0
points = [(5.0, 0.0, NaN), (2.5, 2.5, NaN)]
fullyWithin 4 = true
fullyWithin 3 = false

四个结论。第一,DHD(flat → bump) = 3.5355 是连续版的结果,比离散版的 0.0 大得多——因为连续版会在 flat 的边上取点,中点到 bump 斜边的最短距离是 5/√2。第二,DHD(bump → flat) = 5.0,极值在 bump 的顶点上,与离散版一致(顶点恰好是连续版的极值点)。第三,对称版 HD = max(3.5355, 5.0) = 5.0,这才是严格意义上的 Hausdorff 距离。第四,distancePoints 返回的 (5 0) 与 (2.5 2.5) 是达到最大距离的那对点,其中 (5,0) 是 flat 上的内部点、(2.5,2.5) 是 bump 斜边的中点——离散版根本取不到这对点。最后两行印证了 isFullyWithinDistance 的语义:flat 上每一点到 bump 的最大距离是 3.5355,所以阈值 4 为 true、阈值 3 为 false。

11.5.3 isFullyWithinDistance:比算距离更快

类注释特别推荐了 isFullyWithinDistance:它带启发式短路,一旦确认“存在某点到目标几何的距离超过阈值”就立即返回 false,不必算完整个几何。在“判断某条管线是否完全落在某区域 50 米范围内”这类查询里,它比“算距离再比较”快一个数量级。

import org.locationtech.jts.algorithm.distance.DirectedHausdorffDistance;
import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.io.WKTReader;

public class FullyWithinDemo {
    public static void main(String[] args) throws Exception {
        WKTReader reader = new WKTReader();
        Geometry park = reader.read("POLYGON ((0 0, 100 0, 100 100, 0 100, 0 0))");
        Geometry trail = reader.read("LINESTRING (10 10, 50 20, 90 90)");
        Geometry stray = reader.read("LINESTRING (10 10, 50 20, 120 90)");

        DirectedHausdorffDistance toPark = new DirectedHausdorffDistance(park);
        System.out.println("trail inside park      = " + toPark.isFullyWithinDistance(trail, 1.0));
        System.out.println("stray inside park      = " + toPark.isFullyWithinDistance(stray, 1.0));
        System.out.println("stray within 30 of park= " + toPark.isFullyWithinDistance(stray, 30.0));
        System.out.println("farthest points        = "
                + java.util.Arrays.toString(toPark.farthestPoints(stray)));
    }
}

输出:

trail inside park      = true
stray inside park      = false
stray within 30 of park= true
farthest points        = [(120.0, 90.0, NaN), (100.0, 90.0, NaN)]

注意这里实例化的方向:new DirectedHausdorffDistance(park) 把公园作为目标几何并缓存索引,isFullyWithinDistance(trail, 1.0) 问的是“trail 是否完全在公园 1 米范围内”。第一次调用是 true(整条步道都在公园内),第二次是 false(stray 有 20 米在园外),第三次阈值放到 30 就 true 了。farthestPoints 给出的 (120,90) 与 (100,90) 正是那 20 米越界的端点与对应的园内最近点。

11.6 离散 Fréchet 距离

11.6.1 与 Hausdorff 的本质差异

Hausdorff 只关心“最远有多远”,不关心点的先后顺序;Fréchet 距离则要求“按顺序配对”,可以想象成“一个人牵着狗分别沿两条曲线走,绳子最短需要多长”。因此 Fréchet 对方向、折返、形状顺序都敏感,而 Hausdorff 不敏感。

  曲线 A:  ─────────────►          A 与 B 形状相同但方向相反
  曲线 B:  ◄─────────────
  Hausdorff 距离 ≈ 0(点集几乎重合)
  Fréchet 距离很大(要走完 A 再倒着走 B)

这个差异让 Fréchet 成为“轨迹相似度”的标配度量,而 Hausdorff 更适合“形状偏差容忍度”评估。JTS 的 DiscreteFrechetDistance(1.19 新增)实现的是离散版本:只在顶点之间做配对,不插值。

11.6.2 代码与输出

import java.util.Arrays;

import org.locationtech.jts.algorithm.distance.DiscreteFrechetDistance;
import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.io.WKTReader;

public class FrechetDemo {
    public static void main(String[] args) throws Exception {
        WKTReader reader = new WKTReader();
        Geometry flat = reader.read("LINESTRING (0 0, 10 0)");
        Geometry bump = reader.read("LINESTRING (0 0, 5 5, 10 0)");
        Geometry offset = reader.read("LINESTRING (0 1, 10 1)");

        System.out.println("Frechet flat-bump   = "
                + DiscreteFrechetDistance.distance(flat, bump));
        System.out.println("Frechet flat-offset = "
                + DiscreteFrechetDistance.distance(flat, offset));

        DiscreteFrechetDistance dfd = new DiscreteFrechetDistance(flat, bump);
        System.out.println("distance = " + dfd.distance());
        System.out.println("coords   = " + Arrays.toString(dfd.getCoordinates()));
    }
}

输出:

Frechet flat-bump   = 7.0710678118654755
Frechet flat-offset = 1.0
Frechet distance = 7.0710678118654755
coords   = [(0.0, 0.0, NaN), (5.0, 5.0, NaN)]

flat 与 bump 的离散 Fréchet 距离是 7.071,即顶点 (5,5) 到 (0,0) 的距离 5√2。为什么不是 5?因为顺序配对约束:要同时走完 flat 的两个顶点与 bump 的三个顶点,最优配对里必然出现“flat 的某点配 bump 的 (5,5)”,而 (0,0) 到 (5,5) 的欧氏距离就是 7.071。而 flat 与 offset 是两条平行线、顶点一一对应,距离恒为 1,所以结果是 1.0。

11.6.3 实现细节:为什么它可能很慢

离散 Fréchet 距离的标准算法是动态规划,需要一张 m×n 的表(m、n 是两条线的顶点数)。DiscreteFrechetDistance 内部实现了三种矩阵存储:RectMatrix(稠密二维数组)、CsrMatrix(稀疏行压缩)、HashMapMatrix(哈希),会根据顶点规模自动选择,目的就是在“内存”和“速度”之间取平衡。

工程上的直接后果:两条各有 10 万顶点的轨迹做 Fréchet 距离,即使稀疏化也需要可观的堆内存;而同样规模做 Hausdorff 只需要 O(n log n)。因此轨迹相似度必须先在抽稀后的数据上算(第 12 章的 Douglas-Peucker),或先按路段切分再逐段比较。

11.7 尺度统计:最小直径与外接形状

11.7.1 最小直径

MinimumDiameter 回答“这个形状最窄的方向有多宽”。它的核心概念是“支撑线”(supporting segment)与“宽度坐标”:对凸多边形,最小宽度出现在某条边的平行方向上;实现上先求凸包,再旋转卡壳遍历所有边方向。

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

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

        MinimumDiameter md = new MinimumDiameter(blob);
        System.out.println("length   = " + md.getLength());
        System.out.println("diameter = " + md.getDiameter());
        System.out.println("widthPt  = " + md.getWidthCoordinate());
        System.out.println("support  = " + md.getSupportingSegment());
    }
}

输出:

length   = 4.0
diameter = LINESTRING (0 4, 0 0)
widthPt  = (0.0, 0.0, NaN)
support  = LINESTRING (0 4, 10 4)

宽度 4(矩形的短边),直径线是竖直的 (0 4) → (0 0),支撑线是上边 (0 4) → (10 4)。getWidthCoordinate() 给出的是“宽度方向上的那个最小坐标点”,本例取到下边左端的 (0,0),它只是一个代表点,不代表宽度测量的起点。

11.7.2 最小外接矩形

MinimumAreaRectangle(1.20 重写,早期版本在 MinimumDiameter.getMinimumRectangle 里)返回“面积最小的外接矩形”。它的算法与最小直径同源:最小面积矩形的某条边一定与凸包的一条边共线。

import org.locationtech.jts.algorithm.MinimumAreaRectangle;
import org.locationtech.jts.algorithm.MinimumDiameter;
import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.io.WKTReader;

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

        Geometry rect1 = MinimumAreaRectangle.getMinimumRectangle(blob);
        Geometry rect2 = MinimumDiameter.getMinimumRectangle(blob);
        System.out.println("MinimumAreaRectangle = " + rect1 + "  area=" + rect1.getArea());
        System.out.println("MinimumDiameter      = " + rect2 + "  area=" + rect2.getArea());
    }
}

输出:

MinimumAreaRectangle = POLYGON ((0 0, 0 4, 10 4, 10 0, 0 0))  area=40.0
MinimumDiameter      = POLYGON ((0 4, 0 0, 10 0, 10 4, 0 4))  area=40.0

两者结果相同(顶点顺序与起点不同)。MinimumAreaRectangle 的构造器还有一个 isConvex 参数:若你已知输入已是凸的,传 true 可以跳过凸包计算,省一次 O(n log n)。

11.7.3 最小外接圆

MinimumBoundingCircle 求覆盖全部顶点的最小圆。它的实现不是“迭代优化”,而是精确算法:最小圆要么由两个点确定(直径圆),要么由三个点确定(外接圆)。

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

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

        MinimumBoundingCircle mbc = new MinimumBoundingCircle(tri);
        System.out.println("radius         = " + mbc.getRadius());
        System.out.println("centre         = " + mbc.getCentre());
        System.out.println("diameter       = " + mbc.getDiameter());
        System.out.println("maxDiameter    = " + mbc.getMaximumDiameter());
        System.out.println("farthestPoints = " + mbc.getFarthestPoints());
        System.out.println("extremalPts    = " + java.util.Arrays.toString(mbc.getExtremalPoints()));
    }
}

输出:

radius         = 7.0710678118654755
centre         = (5.0, 5.0, NaN)
diameter       = LINESTRING (0 0, 10 0)
maxDiameter    = LINESTRING (10 0, 0 10)
farthestPoints = LINESTRING (10 0, 0 10)
extremalPts    = [(0.0, 0.0, NaN), (10.0, 0.0, NaN), (0.0, 10.0, NaN)]

注意 getDiameter() 与 getMaximumDiameter() 的差异:前者返回“决定最小圆的那个直径”(本例是 (0,0)–(10,0),但最小圆其实由三个点确定,所以这个直径只是其中一个候选),后者返回“几何内部最远的一对点”。两者只有在最小圆恰好由两点确定时才一致。这类 API 的命名容易误导,用之前先读返回值语义。

11.7.4 最小外接三角形(1.21 新增)

MinimumBoundingTriangle 是 1.21 新增的类(#1160),求面积最小的外接三角形。它的用途比较专门:三角形外接框在制图综合、符号化、粗略面积估计里有时比矩形更紧凑(尤其是细长形状)。

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

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

        MinimumBoundingTriangle mbt = new MinimumBoundingTriangle(tri);
        Geometry triangle = mbt.getTriangle();
        System.out.println("triangle = " + triangle);
        System.out.println("area     = " + triangle.getArea());
    }
}

输出:

triangle = POLYGON ((0 0, 0 10, 10 0, 0 0))
area     = 50.0

对本例的直角三角形,最小外接三角形就是它自身,面积 50。这个类没有静态便捷方法,必须实例化;它内部的算法是枚举凸包三条边方向的组合,复杂度比最小外接矩形高,不适合在超大几何上反复调用。

11.8 内部极值:最大内切圆与最大空圆

11.8.1 最大内切圆

MaximumInscribedCircle(algorithm.construct 包)在面内找“离边界最远的内部点”,以及该点到边界的距离。地理学上这个点叫“不可达极点”(Pole of Inaccessibility),制图上用来放标签最合适——离边界最远意味着标注最不容易压边。

类注释里有三个关键信息:其一,半径就是“负缓冲变空的距离”,因此可以用 isRadiusWithin(maxRadius) 快速判断“这个面是否比给定宽度更窄”,比真算半径快得多;其二,对小多边形(三角形与凸四边形)有精确解,其他多边形用网格细分加分支限界逼近;其三,1.21 新增了自动容差(#1128),不传 tolerance 时由实现根据几何尺度推导。

import org.locationtech.jts.algorithm.construct.MaximumInscribedCircle;
import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.io.WKTReader;

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

        MaximumInscribedCircle mic = new MaximumInscribedCircle(square);
        System.out.println("center     = " + mic.getCenter());
        System.out.println("radiusPt   = " + mic.getRadiusPoint());
        System.out.println("radiusLine = " + mic.getRadiusLine());
        System.out.println("within 5.0 = " + mic.isRadiusWithin(5.0));
        System.out.println("within 4.9 = " + mic.isRadiusWithin(4.9));

        System.out.println("thin center (auto tol) = "
                + MaximumInscribedCircle.getCenter(thin));
        System.out.println("thin center (tol 0.01) = "
                + MaximumInscribedCircle.getCenter(thin, 0.01));
    }
}

输出:

center     = POINT (5.000000000000001 5)
radiusPt   = POINT (10 5)
radiusLine = LINESTRING (5.000000000000001 5, 10 5)
within 5.0 = true
within 4.9 = false
center     = POINT (1.0000000000000002 1)
thin center (auto tol) = POINT (1.0000000000000002 1)
thin center (tol 0.01) = POINT (1.0000000000000002 1)

正方形的结果是精确的 (5,5)(浮点误差 1e-15 量级),半径 5,所以 isRadiusWithin(5.0) 为 true、isRadiusWithin(4.9) 为 false——这个接口的语义是“半径是否不超过给定值”。细长矩形 10×2 的结果是 (1,1),半径 1,这正是短边一半:细长形状的最大内切圆永远贴着窄边。

“自动容差”值得多解释一句。1.21 之前必须显式传 tolerance,传大了结果不准、传小了算得慢;1.21 起不传时实现会按包络尺度自动选一个(源码里的自动容差系数与几何尺度挂钩)。对于“只要知道大概位置”的制图注记场景,自动容差足够;对于“要精确半径做负缓冲判断”的场景,建议显式传一个远小于业务精度的值。

11.8.2 最大空圆

LargestEmptyCircle 是“最大内切圆”的对偶问题:给定一组障碍点(或线、面)与一个可选边界,找“离所有障碍最远的点”以及该点到最近障碍的距离。典型用途是设施选址——在已有设施点之间找“覆盖最薄弱”的位置。

import org.locationtech.jts.algorithm.construct.LargestEmptyCircle;
import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.io.WKTReader;

public class LECDemo {
    public static void main(String[] args) throws Exception {
        WKTReader reader = new WKTReader();
        Geometry obstacles = reader.read("MULTIPOINT ((1 1), (9 1), (1 9), (9 9))");
        Geometry boundary = reader.read("POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0))");

        LargestEmptyCircle lec = new LargestEmptyCircle(obstacles, boundary);
        System.out.println("center     = " + lec.getCenter());
        System.out.println("radiusLine = " + lec.getRadiusLine());
        System.out.println("static     = " + LargestEmptyCircle.getCenter(obstacles, boundary));
    }
}

输出:

center     = POINT (5 5)
radiusLine = LINESTRING (5 5, 1 1)
static     = POINT (5 5)

四个角上各有一个障碍,最空的位置就是中心 (5,5),到最近障碍的距离是 4√2 ≈ 5.657,radiusLine 给出这条半径线。构造器的第二个参数是“边界”,它把搜索范围限制在园内;若只传一个参数(只有障碍),搜索范围会扩展到障碍的包络(1.20 之前的行为),结果可能跑到数据范围之外。

类注释提到 1.20 的两处改进:支持指定边界(#973),以及修复薄输入长时间循环(#978)。如果你的版本较老,在细长区域上跑 LargestEmptyCircle 可能非常慢,升级即可。

11.9 三维距离、高精度数值与向量长度

11.9.1 Distance3DOp

Distance3DOp 计算两个几何在三维空间中的最小距离。它支持点、线、面,但对“面”有一个重要限制:Polygon 在三维中被当作“平面多边形”处理,若输入面不共面,结果没有明确定义。源码里为此提供了 PlanarPolygon3D 做平面拟合。

import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.io.WKTReader;
import org.locationtech.jts.operation.distance3d.Distance3DOp;

public class Distance3DDemo {
    public static void main(String[] args) throws Exception {
        WKTReader reader = new WKTReader();
        Geometry p1 = reader.read("POINT Z (0 0 0)");
        Geometry p2 = reader.read("POINT Z (0 0 5)");
        Geometry l1 = reader.read("LINESTRING Z (0 0 0, 10 0 0)");
        Geometry l2 = reader.read("LINESTRING Z (0 5 0, 10 5 3)");

        System.out.println("point-point   = " + Distance3DOp.distance(p1, p2));
        System.out.println("point-point np= "
                + java.util.Arrays.toString(Distance3DOp.nearestPoints(p1, p2)));
        System.out.println("line-line     = " + Distance3DOp.distance(l1, l2));
    }
}

输出:

point-point   = 5.0
point-point np= [(0.0, 0.0, 0.0), (0.0, 0.0, 5.0)]
line-line     = 5.0

两条三维线段的距离是 5,说明最近点落在“竖直方向”上:第一条线在 z=0 平面,第二条从 z=0 升到 z=3,最近的那对点都在 z=0 端。这与平面距离计算的结果不同——平面距离会把 z 忽略掉,得到 0(因为两条线在 xy 平面上平行、相距 5,但投影后仍有 5 的间距,其实也是 5;换成真正交叉的线才能看出差异)。

11.9.2 DD:双精度双倍精度

org.locationtech.jts.math.DD 是用两个 double 表示一个高精度数的类(double-double 技术),有效位数约 32 位十进制。它存在的理由是:某些几何判定(如方向判定、圆外心计算)在双精度下会因舍入给出错误符号,而用 DD 可以把误差推后到几乎不可能发生的量级。

import org.locationtech.jts.math.DD;

public class DDDemo {
    public static void main(String[] args) {
        DD a = DD.valueOf(0.1);
        DD b = DD.valueOf(0.2);
        System.out.println("DD  0.1+0.2 = " + a.add(b));
        System.out.println("double 0.1+0.2 = " + (0.1 + 0.2));
        System.out.println("DD doubleValue = " + a.add(b).doubleValue());
        System.out.println("DD  sqrt(2)    = " + DD.valueOf(2).sqrt());
        System.out.println("DD  PI         = " + DD.PI);
        System.out.println("DD  compare    = " + a.compareTo(b));
    }
}

输出:

DD  0.1+0.2 = 0.30000000000000001665334536937734
double 0.1+0.2 = 0.30000000000000004
DD doubleValue = 0.3
DD  sqrt(2)    = 1.4142135623730950488016887242096
DD  PI         = 3.1415926535897932384626433832795
DD  compare    = -1

注意 DD.valueOf(0.1).add(DD.valueOf(0.2)).doubleValue() 返回 0.3,而直接用 double 相加是 0.30000000000000004。DD 的对外表示(toString)会打印 32 位有效数字,看起来很“精确”,但它的输入仍是 double——如果要保留十进制字面量的精度,应该用 DD.valueOf(String) 而不是 DD.valueOf(double)。

DD 在 JTS 里的主要用途是内部算法支撑(CGAlgorithmsDD、RobustDeterminant),业务代码很少直接调用。需要高精度距离或面积时,更实际的做法是用 BigDecimal 在应用层算,或在坐标系层面选择更好的投影。

11.9.3 MathUtil.hypot:快,但不防溢出

MathUtil.hypot(x, y) 的实现就是 Math.sqrt(x * x + y * y),源码注释明确说它“比标准 Math.hypot 快”。快有代价——它不做溢出与下溢保护:

import org.locationtech.jts.math.MathUtil;

public class HypotDemo {
    public static void main(String[] args) {
        System.out.println("MathUtil.hypot(3, 4)        = " + MathUtil.hypot(3, 4));
        System.out.println("MathUtil.hypot(1e200, 1e200)= " + MathUtil.hypot(1e200, 1e200));
        System.out.println("Math.hypot(1e200, 1e200)    = " + Math.hypot(1e200, 1e200));
        System.out.println("MathUtil.hypot(1e-200,1e-200)= " + MathUtil.hypot(1e-200, 1e-200));
    }
}

输出:

MathUtil.hypot(3, 4)        = 5.0
MathUtil.hypot(1e200, 1e200)= Infinity
Math.hypot(1e200, 1e200)    = 1.414213562373095E200
MathUtil.hypot(1e-200,1e-200)= 0.0

这个差异在 GIS 里几乎不会触发(坐标量级远小于 1e154),但它说明了 JTS 的一个取舍:内部计算优先速度,边界情形交给上层。如果你在做跨星系的坐标运算(比如天文数据),请自己用 Math.hypot 包装。

11.10 单位、投影与统计口径

11.10.1 所有结果都是“坐标单位”

这是本章最重要的一句话:JTS 的所有距离、面积、长度都是坐标单位,它不知道也不关心坐标代表什么。经纬度坐标下 getLength() 返回的是“度”,面积返回的是“平方度”,两者都没有物理意义。

坐标类型 长度单位 面积单位 1 单位约等于
WGS84 经纬度 度 平方度 约 111 km(纬度方向)、随纬度变化的经度方向
Web Mercator 米(投影米) 平方米 随纬度严重变形,高纬放大
UTM 50N 米 平方米 1 米(带内精度高)
地方平面坐标 米或尺 平方米或平方尺 取决于定义

工程上的做法是:要么先把数据投影到米制坐标系再算,要么用 GeoTools 等库的测地线 API。如果一定要在经纬度上算,只能用“近似公式”:纬度 1 度约 111.32 km,经度 1 度约 111.32×cos(纬度) km,且只适合小范围。

11.10.2 面积口径:投影面积 vs 椭球面积

即使在米制投影下,面积也只是“投影平面上的面积”,不等于椭球面上的真实面积。UTM 带内的长度变形通常小于 0.1%,面积变形约 0.2%;Web Mercator 在高纬可以放大数倍。需要法定面积(土地登记、林业调查)时,必须用专门的椭球面积算法(GeoTools 的 GeodeticCalculator 或业务系统自带的面积计算),不能直接取 getArea()。

一个直观的数量级:在北纬 60 度,Web Mercator 的长度放大倍率约 2,面积放大约 4。也就是说同一块地,在 Web Mercator 下算出的面积是真实面积的 4 倍。这个误差远超任何几何算法误差,是“坐标系错误”而非“算法错误”。

11.10.3 统计口径的一致性

做统计时还有一个隐蔽问题:同一次分析里混用了不同精度的几何。例如边界数据是 1:10000 精度(顶点间距几十米),业务数据是 1:500 精度(顶点间距几米),两者叠加或求交后算出的面积,误差主要来自边界数据的“锯齿”。这不是算法能修的,只能在数据治理层面统一精度(用 PrecisionReducer 或简化算法,见第 12 章)。

11.11 常见问题与避坑

问题 1:DistanceOp 对空几何返回 0,导致误判为“相交”。 个体行为如此(见 11.2.3)。判断相交请用 intersects() 或 RelateNG,不要用“距离是否为 0”。

问题 2:两个类算出的距离不一样。 先确认问的是“集合距离”还是“到边界 facet 的距离”(11.3.1)。DistanceOp 与 IndexedFacetDistance 对多边形内部点会给出 0 与正数两个不同结果,两者都正确。

问题 3:Hausdorff 距离随顶点密度变化。 离散版本只看顶点(11.4.1)。要消除顶点依赖,用 DirectedHausdorffDistance 的连续版,或给离散版传 densifyFrac。注意 densifyFrac 只是加密取样,并不能完全消除依赖。

问题 4:Hausdorff 距离为 0,但两个形状明显不同。 典型的“顶点都落在对方身上”(11.4.2 的 oriented(flat,bump) = 0)。改用连续版或对称版;若两个形状相同但方向相反,Hausdorff 也会给出 0,此时应看 Fréchet。

问题 5:Fréchet 距离在长轨迹上 OOM。 离散 Fréchet 是 O(m×n) 动态规划(11.6.3)。先抽稀(Douglas-Peucker,第 12 章),或按路段切片后逐段比较。

问题 6:最小外接矩形算出来是斜的。 这是正确的——最小面积矩形的边与凸包某条边共线,不一定与坐标轴平行。要“水平外接框”请用 getEnvelopeInternal() 或 getEnvelope()。

问题 7:MaximumInscribedCircle 在细长多边形上很慢。 1.20 修复过薄输入长循环问题(#978)。升级到 1.20 以上;仍有性能问题时显式传一个合适的 tolerance,避免自动容差过小。

问题 8:LargestEmptyCircle 的结果跑到数据范围之外。 只传障碍、不传边界时,搜索范围是障碍包络。要限制在业务范围内,用两参数的构造器或 getCenter(obstacles, boundary)。

问题 9:Distance3DOp 对不共面的多边形给出奇怪结果。 三维多边形被当作平面处理(11.9.1)。先检查顶点是否共面,不共面时拆成三角形或投影到平面再算。

问题 10:面积算出天文数字。 先看坐标系:经纬度下算面积得到的是“平方度”,数值可能很小且无意义;Web Mercator 下高纬面积被放大数倍。解决方法是投影到等积投影(如 Albers)再算,或用椭球面积算法。

11.12 本章小结

本章覆盖了 JTS 里“算数”的所有主要入口,要点可以浓缩成十二条:

  1. getArea() 是净面积(扣孔),getLength() 是所有环的周长之和,不是外环周长;两者都遵循“坐标单位”。
  2. getCentroid() 不保证落在几何内部,getInteriorPoint() 保证落在内部;制图注记用后者。
  3. getBoundaryDimension() 给出了边界维度约定:点 -1、线 0、面 1、环 -1;它解释了“相邻多边形的交是线”这类现象。
  4. DistanceOp 算集合间最小距离,对空几何返回 0,对内部点返回 0;isWithinDistance 带提前终止,比“算距离再比较”快。
  5. IndexedFacetDistance 算到 facet(线段与顶点)的距离,对多边形内部点返回正数;它缓存 STRtree,适合批量查询,且线程安全。
  6. 离散 Hausdorff 只看顶点,极值可能落在顶点上而给出 0;densifyFrac 可缓解但无法根除顶点依赖。
  7. DirectedHausdorffDistance(1.21 新增)是连续版,更准确也更快;hausdorffDistance 才是对称度量,isFullyWithinDistance 带短路,适合布尔查询。
  8. 离散 Fréchet 要求顺序配对,对方向与折返敏感,是轨迹相似度的标配;复杂度 O(m×n),长轨迹必须先抽稀。
  9. 尺度统计四件套:MinimumDiameter 给最窄宽度,MinimumAreaRectangle 给最小外接矩形,MinimumBoundingCircle 给最小外接圆,MinimumBoundingTriangle(1.21)给最小外接三角形。
  10. MaximumInscribedCircle 给“最远离边界”的内点与半径(可用于判断负缓冲是否为空),LargestEmptyCircle 给“最远离障碍”的点(用于选址)。
  11. Distance3DOp 提供三维距离,但面被当作平面处理;DD 提供约 32 位有效数字的高精度计算,主要服务于内部鲁棒算法。
  12. 所有结果都是坐标单位;经纬度下算长度面积没有物理意义,必须投影;Web Mercator 高纬面积会放大数倍,法定面积必须用椭球算法。

下一章进入几何化简与形状概括:本章解决“量出来是多少”,下一章解决“怎么在保持形状特征的前提下把顶点减下来”,以及“用什么几何代表一片区域的轮廓最合适”。

相关资源

  • JTS 官方 Javadoc:https://locationtech.github.io/jts/javadoc/ ,重点读 DistanceOp、IndexedFacetDistance、DirectedHausdorffDistance、MaximumInscribedCircle、LargestEmptyCircle 的类注释;本章的语义差异与限制说明均来自这些注释。
  • 源码位置:modules/core/src/main/java/org/locationtech/jts/algorithm/(MinimumDiameter、MinimumBoundingCircle、MinimumAreaRectangle、MinimumBoundingTriangle)、algorithm/distance/(Hausdorff 与 Fréchet)、algorithm/construct/(最大内切圆、最大空圆)、operation/distance/(DistanceOp、IndexedFacetDistance)、math/(DD、MathUtil)。
  • 版本历史:doc/JTS_Version_History.md。MinimumBoundingTriangle 与 DirectedHausdorffDistance 见 1.21 条目(#1160、#1182);DiscreteFrechetDistance 见 1.19 条目(#764、#783);最大内切圆自动容差与 isRadiusWithin 见 1.21(#1128、#1125)。
  • 本站第 04 章《几何对象模型详解》:https://znlgis.github.io/gis/tutorial/jts/第04章-几何对象模型详解/
  • 本站第 10 章《缓冲区与偏移曲线》:https://znlgis.github.io/gis/tutorial/jts/第10章-缓冲区与偏移曲线/
  • 本站第 12 章《几何化简、凸包与形状概括》:https://znlgis.github.io/gis/tutorial/jts/第12章-几何化简凸包与形状概括/
  • 本站第 19 章《坐标系、单位与工程陷阱》:https://znlgis.github.io/gis/tutorial/jts/第19章-坐标系单位与工程陷阱/