第24章 - 综合实战案例

package com.example.cases;

import java.util.ArrayList;
import java.util.HashMap;
import java.util.List;
import java.util.Map;

import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.geom.GeometryFactory;
import org.locationtech.jts.geom.PrecisionModel;
import org.locationtech.jts.index.strtree.STRtree;
import org.locationtech.jts.io.WKTReader;
import org.locationtech.jts.operation.overlayng.OverlayNG;

/**
 * 案例一：片区 × 行政区叠加统计。
 * <p>
 * 输出每个片区在不同行政区里的面积与占比，并做面积守恒校验。
 */
public class Case01ZonalStats {

  /** 叠加用的精度模型：1e6 表示保留 6 位小数，足够米制数据使用。 */
  private static final PrecisionModel PM = new PrecisionModel(1e6);

  private final GeometryFactory gf = new GeometryFactory(PM);
  private final WKTReader reader = new WKTReader(gf);

  /** 地区：名称 + 几何 */
  static class Zone {
    final String name;
    final Geometry geom;

    Zone(String name, Geometry geom) {
      this.name = name;
      this.geom = geom;
    }
  }

  /** 一条统计结果 */
  static class Stat {
    final String unit;
    final String zone;
    final double area;

    Stat(String unit, String zone, double area) {
      this.unit = unit;
      this.zone = zone;
      this.area = area;
    }
  }

  public List<Stat> compute(List<Zone> units, List<Zone> zones) throws Exception {
    // 1. 为行政区建 STRtree（按包络索引）
    STRtree index = new STRtree();
    for (Zone z : zones) {
      index.insert(z.geom.getEnvelopeInternal(), z);
    }
    index.build();

    List<Stat> stats = new ArrayList<Stat>();

    // 2. 逐片区计算
    for (Zone unit : units) {
      double unitArea = unit.geom.getArea();
      double sumArea = 0.0;

      @SuppressWarnings("unchecked")
      List<Zone> candidates = index.query(unit.geom.getEnvelopeInternal());

      for (Zone zone : candidates) {
        if (!unit.geom.intersects(zone.geom)) continue;

        Geometry inter = OverlayNG.overlay(unit.geom, zone.geom, OverlayNG.INTERSECTION, PM);
        if (inter.isEmpty()) continue;

        double area = inter.getArea();
        sumArea += area;
        stats.add(new Stat(unit.name, zone.name, area));
      }

      // 3. 面积守恒校验
      double diff = Math.abs(unitArea - sumArea);
      if (diff > 1e-6) {
        System.err.printf("WARN: unit %s area mismatch: total=%.9f sum=%.9f diff=%.9f%n",
            unit.name, unitArea, sumArea, diff);
      }
    }
    return stats;
  }

  public static void main(String[] args) throws Exception {
    Case01ZonalStats app = new Case01ZonalStats();

    List<Zone> zones = new ArrayList<Zone>();
    zones.add(new Zone("A区", app.reader.read(
        "POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0))")));
    zones.add(new Zone("B区", app.reader.read(
        "POLYGON ((10 0, 20 0, 20 10, 10 10, 10 0))")));

    List<Zone> units = new ArrayList<Zone>();
    units.add(new Zone("片区1", app.reader.read(
        "POLYGON ((0 0, 20 0, 20 10, 0 10, 0 0))")));

    List<Stat> stats = app.compute(units, zones);

    // 4. 按片区汇总，输出交叉表
    Map<String, Double> unitTotals = new HashMap<String, Double>();
    for (Stat s : stats) {
      Double v = unitTotals.get(s.unit);
      unitTotals.put(s.unit, (v == null ? 0.0 : v) + s.area);
    }

    System.out.println("unit\tzone\tarea\tshare");
    for (Stat s : stats) {
      double share = s.area / unitTotals.get(s.unit);
      System.out.printf("%s\t%s\t%.3f\t%.1f%%%n", s.unit, s.zone, s.area, share * 100);
    }
  }
}

package com.example.cases;

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

import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.geom.GeometryFactory;
import org.locationtech.jts.geom.LineString;
import org.locationtech.jts.io.WKTReader;
import org.locationtech.jts.operation.linemerge.LineMerger;
import org.locationtech.jts.operation.linemerge.LineSequencer;
import org.locationtech.jts.operation.union.UnaryUnionOp;

/**
 * 案例二：道路网连通性处理。
 * <p>
 * 演示节点化 → 线合并 → 欧拉路径判定的标准流水线。
 */
public class Case02RoadNetwork {

  private final GeometryFactory gf = new GeometryFactory();
  private final WKTReader reader = new WKTReader(gf);

  /** 节点化：在交叉点插入节点，并去掉重复线段。 */
  public Geometry node(Geometry linework) {
    return UnaryUnionOp.union(linework);
  }

  /** 线合并：在度 1 或度 >= 3 的节点处断开，度 2 节点合并。 */
  public List<LineString> merge(Geometry linework) {
    LineMerger merger = new LineMerger();
    merger.add(linework);

    List<LineString> result = new ArrayList<LineString>();
    Collection<?> merged = merger.getMergedLineStrings();
    for (Object o : merged) {
      result.add((LineString) o);
    }
    return result;
  }

  public static void main(String[] args) throws Exception {
    Case02RoadNetwork app = new Case02RoadNetwork();

    // 三条原始线，其中 (2,0) 是贯穿线与主干道的非端点交点
    Geometry raw = app.reader.read(
        "MULTILINESTRING ((0 0, 10 0), (10 0, 10 10), (0 10, 10 10), (2 -5, 2 5))");

    System.out.println("== 原始线 ==");
    System.out.println("段数 = " + raw.getNumGeometries());
    System.out.println("总长度 = " + raw.getLength());

    // 1. 节点化
    Geometry noded = app.node(raw);
    System.out.println();
    System.out.println("== 节点化后 ==");
    System.out.println("段数 = " + noded.getNumGeometries());
    System.out.println("总长度 = " + noded.getLength());

    // 2. 线合并
    List<LineString> merged = app.merge(noded);
    System.out.println();
    System.out.println("== 线合并后 ==");
    System.out.println("合并线条数 = " + merged.size());
    double total = 0;
    for (LineString line : merged) {
      System.out.println("  " + line);
      total += line.getLength();
    }
    System.out.println("长度合计 = " + total);

    // 3. 欧拉路径判定
    boolean sequenceable;
    LineSequencer sequencer = new LineSequencer();
    sequencer.add(noded);
    sequenceable = sequencer.isSequenceable();
    System.out.println();
    System.out.println("== 欧拉路径 ==");
    System.out.println("isSequenceable = " + sequenceable);

    if (sequenceable) {
      Geometry sequenced = sequencer.getSequencedLineStrings();
      System.out.println("isSequenced = " + LineSequencer.isSequenced(sequenced));
      System.out.println(sequenced);
    }
    else {
      System.out.println("存在超过 2 个奇度节点，无法一次走完不重复");
    }
  }
}

package com.example.cases;

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

import org.locationtech.jts.coverage.CoverageGapFinder;
import org.locationtech.jts.coverage.CoverageUnion;
import org.locationtech.jts.coverage.CoverageValidator;
import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.geom.GeometryFactory;
import org.locationtech.jts.geom.Polygon;
import org.locationtech.jts.geom.util.GeometryFixer;
import org.locationtech.jts.io.WKTReader;
import org.locationtech.jts.operation.valid.IsValidOp;
import org.locationtech.jts.operation.valid.TopologyValidationError;

/**
 * 案例三：宗地质检。
 * <p>
 * 逐宗有效性校验 → 修复 → 覆盖校验 → 缝隙检测 → 覆盖合并。
 */
public class Case03ParcelQC {

  private final GeometryFactory gf = new GeometryFactory();
  private final WKTReader reader = new WKTReader(gf);

  /** 单宗地校验 + 修复。返回修复后的几何（已确保有效）。 */
  public Geometry checkAndFix(String id, Geometry parcel) {
    if (parcel.isValid()) {
      System.out.println("  [" + id + "] valid");
      return parcel;
    }

    IsValidOp op = new IsValidOp(parcel);
    TopologyValidationError err = op.getValidationError();
    System.out.println("  [" + id + "] INVALID: "
        + (err == null ? "unknown" : err.getMessage() + " @ " + err.getCoordinate()));

    Geometry fixed = GeometryFixer.fix(parcel);
    System.out.println("  [" + id + "] fixed  -> " + fixed.getGeometryType()
        + ", area=" + fixed.getArea() + ", valid=" + fixed.isValid());
    return fixed;
  }

  /** 构造一个 10x10 的正方形格网单元。 */
  private static String cell(double x, double y, double size) {
    return "POLYGON ((" + x + " " + y + ", " + (x + size) + " " + y
        + ", " + (x + size) + " " + (y + size) + ", " + x + " " + (y + size)
        + ", " + x + " " + y + "))";
  }

  public static void main(String[] args) throws Exception {
    Case03ParcelQC app = new Case03ParcelQC();

    // ---- 1. 单宗地校验与修复 ----
    System.out.println("== 1. 单宗地校验 ==");
    List<Geometry> parcels = new ArrayList<Geometry>();
    parcels.add(app.checkAndFix("P1-蝴蝶结",
        app.reader.read("POLYGON ((0 0, 10 10, 10 0, 0 10, 0 0))")));
    parcels.add(app.checkAndFix("P2-洞在壳外",
        app.reader.read("POLYGON ((20 0, 30 0, 30 10, 20 10, 20 0),"
            + " (40 0, 50 0, 50 10, 40 10, 40 0))")));

    // ---- 2. 覆盖校验：两个重叠宗地 ----
    System.out.println();
    System.out.println("== 2. 覆盖校验（重叠） ==");
    Geometry[] overlap = new Geometry[] {
        app.reader.read(cell(0, 0, 10)),
        app.reader.read(cell(8, 0, 10))
    };
    System.out.println("CoverageValidator.isValid = " + CoverageValidator.isValid(overlap));

    Geometry[] errors = CoverageValidator.validate(overlap);
    for (int i = 0; i < errors.length; i++) {
      if (!errors[i].isEmpty()) {
        System.out.println("  parcel[" + i + "] 非法边段数 = " + errors[i].getNumGeometries());
        System.out.println("    " + errors[i]);
      }
    }

    // ---- 3. 覆盖校验：严丝合缝的两个宗地 ----
    System.out.println();
    System.out.println("== 3. 覆盖校验（相邻） ==");
    Geometry[] adjacent = new Geometry[] {
        app.reader.read(cell(0, 0, 10)),
        app.reader.read(cell(10, 0, 10))
    };
    System.out.println("CoverageValidator.isValid = " + CoverageValidator.isValid(adjacent));

    // ---- 4. 缝隙检测：九宫格去掉中心，围出一个 10x10 的洞 ----
    System.out.println();
    System.out.println("== 4. 内部孔洞检测 ==");
    List<Geometry> ringCells = new ArrayList<Geometry>();
    for (int x = 0; x < 3; x++) {
      for (int y = 0; y < 3; y++) {
        if (x == 1 && y == 1) continue;      // 挖掉中心
        ringCells.add(app.reader.read(cell(x * 10, y * 10, 10)));
      }
    }
    Geometry[] ring = ringCells.toArray(new Geometry[0]);

    System.out.println("覆盖单元数 = " + ring.length
        + ", coverageValid = " + CoverageValidator.isValid(ring));

    Geometry merged = CoverageUnion.union(ring);
    System.out.println("合并后类型 = " + merged.getGeometryType()
        + ", 面积 = " + merged.getArea()
        + ", 洞数 = " + ((Polygon) merged).getNumInteriorRing());

    // 阈值 20 > 洞宽 10：检出；阈值 5 < 洞宽 10：不检出
    System.out.println("findGaps(20.0) = " + CoverageGapFinder.findGaps(ring, 20.0));
    System.out.println("findGaps(5.0)  = " + CoverageGapFinder.findGaps(ring, 5.0));

    // ---- 5. 覆盖合并 vs 普通 union ----
    System.out.println();
    System.out.println("== 5. 覆盖合并 ==");
    Geometry coverageUnion = CoverageUnion.union(adjacent);
    Geometry plainUnion = adjacent[0].union(adjacent[1]);
    System.out.println("CoverageUnion 面积 = " + coverageUnion.getArea());
    System.out.println("普通 union 面积    = " + plainUnion.getArea());
    System.out.println("面积一致 = "
        + (Math.abs(coverageUnion.getArea() - plainUnion.getArea()) < 1e-9));
  }
}

package com.example.cases;

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.buffer.BufferOp;
import org.locationtech.jts.operation.buffer.BufferParameters;
import org.locationtech.jts.operation.union.UnaryUnionOp;

/**
 * 案例四：缓冲区分区与覆盖率统计。
 * <p>
 * 演示 quadrantSegments 对精度的影响、覆盖合并、覆盖率计算。
 */
public class Case04BufferZones {

  private final GeometryFactory gf = new GeometryFactory();
  private final WKTReader reader = new WKTReader(gf);

  /** 用指定象限段数给点做缓冲。 */
  public Geometry bufferPoint(Geometry pt, double radius, int quadSegs) {
    BufferParameters params = new BufferParameters();
    params.setQuadrantSegments(quadSegs);
    params.setEndCapStyle(BufferParameters.CAP_ROUND);
    return BufferOp.bufferOp(pt, radius, params);
  }

  public static void main(String[] args) throws Exception {
    Case04BufferZones app = new Case04BufferZones();

    Geometry storeA = app.reader.read("POINT (0 0)");
    Geometry storeB = app.reader.read("POINT (60 0)");

    // ---- 1. 不同精度下的圆面积 ----
    System.out.println("== 1. quadrantSegments 对面积的影响 ==");
    double radius = 50.0;
    double exact = Math.PI * radius * radius;
    System.out.printf("理论圆面积 = %.3f%n", exact);
    for (int q : new int[] {4, 8, 16, 64}) {
      Geometry circle = app.bufferPoint(storeA, radius, q);
      System.out.printf("quadSegs=%2d  顶点数=%3d  面积=%.3f  相对误差=%.3f%%%n",
          q, circle.getNumPoints(), circle.getArea(),
          (circle.getArea() - exact) / exact * 100);
    }

    // ---- 2. 两个店的覆盖与合并 ----
    System.out.println();
    System.out.println("== 2. 覆盖合并 ==");
    Geometry zA = app.bufferPoint(storeA, radius, 8);
    Geometry zB = app.bufferPoint(storeB, radius, 8);
    System.out.printf("店A覆盖面积 = %.3f%n", zA.getArea());
    System.out.printf("店B覆盖面积 = %.3f%n", zB.getArea());

    Geometry union = UnaryUnionOp.union(Arrays.asList(zA, zB));
    System.out.printf("合并后面积 = %.3f（若 == 2 倍单店面积说明不重叠）%n", union.getArea());
    System.out.printf("重叠面积   = %.3f%n", zA.getArea() + zB.getArea() - union.getArea());

    // ---- 3. 覆盖率 ----
    System.out.println();
    System.out.println("== 3. 研究区域覆盖率 ==");
    Geometry studyArea = app.reader.read("POLYGON ((0 -100, 200 -100, 200 100, 0 100, 0 -100))");
    Geometry covered = studyArea.intersection(union);
    System.out.printf("研究区面积 = %.3f%n", studyArea.getArea());
    System.out.printf("被覆盖面积 = %.3f%n", covered.getArea());
    System.out.printf("覆盖率     = %.3f%%%n", covered.getArea() / studyArea.getArea() * 100);

    // ---- 4. 等时圈近似：多边形腐蚀/膨胀 ----
    System.out.println();
    System.out.println("== 4. 成本距离近似（腐蚀/膨胀） ==");
    Geometry zone = app.reader.read("POLYGON ((0 0, 100 0, 100 100, 0 100, 0 0))");
    System.out.printf("原始面积           = %.3f%n", zone.getArea());
    System.out.printf("向内腐蚀 10 后面积 = %.3f%n", zone.buffer(-10).getArea());
    System.out.printf("向外膨胀 10 后面积 = %.3f%n", zone.buffer(10).getArea());

    // ---- 5. 点集批量缓冲（设施覆盖率） ----
    System.out.println();
    System.out.println("== 5. 多点批量缓冲 ==");
    List<Geometry> points = new ArrayList<Geometry>();
    for (int i = 0; i < 5; i++) {
      points.add(app.gf.createPoint(new Coordinate(i * 30, 0)));
    }
    List<Geometry> buffers = new ArrayList<Geometry>();
    for (Geometry p : points) {
      buffers.add(app.bufferPoint(p, 50, 8));
    }
    Geometry merged = UnaryUnionOp.union(buffers);
    System.out.println("5 个点合并后类型 = " + merged.getGeometryType()
        + ", 面积 = " + String.format("%.3f", merged.getArea()));
  }
}

package com.example.cases;

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.linearref.LengthIndexedLine;
import org.locationtech.jts.simplify.DouglasPeuckerSimplifier;
import org.locationtech.jts.simplify.TopologyPreservingSimplifier;
import org.locationtech.jts.simplify.VWSimplifier;

/**
 * 案例五：轨迹抽稀与里程统计。
 */
public class Case05TrackSimplify {

  private final GeometryFactory gf = new GeometryFactory();
  private final WKTReader reader = new WKTReader(gf);

  /** 抽稀并报告里程误差。 */
  public void report(String name, Geometry original, Geometry simplified) {
    double len0 = original.getLength();
    double len1 = simplified.getLength();
    double err = len0 == 0 ? 0 : (len0 - len1) / len0 * 100;
    System.out.printf("%-28s 点数 %2d -> %2d   里程 %.6f -> %.6f   损失 %.4f%%%n",
        name, original.getNumPoints(), simplified.getNumPoints(), len0, len1, err);
  }

  public static void main(String[] args) throws Exception {
    Case05TrackSimplify app = new Case05TrackSimplify();
    Geometry track = app.reader.read(
        "LINESTRING (0 0, 10 0, 20 0.5, 30 1, 40 1, 50 0.5, 60 0)");

    // ---- 1. 原始里程 ----
    System.out.println("== 1. 原始轨迹 ==");
    System.out.printf("点数 = %d, 里程 = %.6f%n",
        track.getNumPoints(), track.getLength());

    // ---- 2. 三种化简器对比 ----
    System.out.println();
    System.out.println("== 2. 化简器对比（容差 = 1.0） ==");
    double tol = 1.0;
    app.report("DouglasPeucker", track, DouglasPeuckerSimplifier.simplify(track, tol));
    app.report("TopologyPreserving", track, TopologyPreservingSimplifier.simplify(track, tol));
    app.report("Visvalingam-Whyatt", track, VWSimplifier.simplify(track, tol));

    System.out.println();
    System.out.println("== 3. 不同容差下的 DP 化简 ==");
    for (double t : new double[] {0.1, 0.5, 1.0, 2.0}) {
      app.report("DP tol=" + t, track, DouglasPeuckerSimplifier.simplify(track, t));
    }

    // ---- 4. 线性参考：按里程截取 ----
    System.out.println();
    System.out.println("== 4. 线性参考 ==");
    LengthIndexedLine lil = new LengthIndexedLine(track);
    System.out.printf("起止里程 = [%.6f, %.6f]%n", lil.getStartIndex(), lil.getEndIndex());

    Coordinate query = new Coordinate(30, 5);
    System.out.printf("点 (30,5) 投影里程 = %.6f%n", lil.project(query));

    Geometry seg = lil.extractLine(20, 40);      // 第 20~40 里程段
    System.out.println("20~40 里程段 = " + seg);
    System.out.printf("该段长度 = %.6f（应 ≈ 20）%n", seg.getLength());

    Coordinate mid = lil.extractPoint(25);
    System.out.printf("里程 25 处坐标 = (%.6f, %.6f)%n", mid.x, mid.y);

    // ---- 5. 抽稀后的里程校验 ----
    System.out.println();
    System.out.println("== 5. 抽稀后按里程截取（注意里程已变） ==");
    Geometry simplified = DouglasPeuckerSimplifier.simplify(track, 1.0);
    LengthIndexedLine lil2 = new LengthIndexedLine(simplified);
    System.out.printf("抽稀后总里程 = %.6f%n", lil2.getEndIndex());
    Geometry seg2 = lil2.extractLine(20, 40);
    System.out.println("抽稀后 20~40 段 = " + seg2);
    System.out.printf("长度 = %.6f%n", seg2.getLength());
  }
}

package com.example.cases;

import java.io.File;
import java.io.FileInputStream;
import java.io.InputStreamReader;
import java.io.Reader;
import java.nio.charset.StandardCharsets;
import java.util.ArrayList;
import java.util.HashMap;
import java.util.HashSet;
import java.util.List;
import java.util.Map;
import java.util.Set;
import java.util.stream.Collectors;
import java.util.stream.IntStream;

import org.locationtech.jts.geom.Coordinate;
import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.geom.GeometryFactory;
import org.locationtech.jts.geom.Point;
import org.locationtech.jts.geom.Polygon;
import org.locationtech.jts.index.strtree.STRtree;
import org.locationtech.jts.io.geojson.GeoJsonReader;

/**
 * 案例六：GeoJSON 批量 ETL + 空间连接。
 * <p>
 * 演示流式解析、STRtree 索引、envelope 粗筛 + contains 精判、并行分片。
 */
public class Case06GeoJsonETL {

  /** 每个线程一个 Reader（GeoJsonReader 非线程安全）。 */
  private static final ThreadLocal<GeoJsonReader> READERS =
      new ThreadLocal<GeoJsonReader>() {
        @Override
        protected GeoJsonReader initialValue() {
          return new GeoJsonReader();
        }
      };

  /** 行政区：名称 + 几何，作为 STRtree 的条目。 */
  private static class Zone {
    final String name;
    final Geometry geom;

    Zone(String name, Geometry geom) {
      this.name = name;
      this.geom = geom;
    }
  }

  private final GeometryFactory gf = new GeometryFactory();
  private final STRtree zoneIndex = new STRtree();

  /** 构建行政区索引：只在启动时做一次。 */
  public void loadZones(List<Geometry> zones, List<String> names) {
    for (int i = 0; i < zones.size(); i++) {
      zoneIndex.insert(zones.get(i).getEnvelopeInternal(), new Zone(names.get(i), zones.get(i)));
    }
    zoneIndex.build();
  }

  /**
   * 把一个点挂到所在行政区。返回 null 表示没有命中（数据问题）。
   * <p>
   * envelope 由 STRtree 负责粗筛，这里只做 contains 精判；
   * 索引条目直接带了几何，避免再按名字反查数组。
   */
  public String locate(Geometry point) {
    @SuppressWarnings("unchecked")
    List<Zone> candidates = zoneIndex.query(point.getEnvelopeInternal());
    for (Zone zone : candidates) {
      if (zone.geom.contains(point)) {
        return zone.name;
      }
    }
    return null;
  }

  /** 从 GeoJSON 文件流式读 Feature 的几何（示例只取第一个几何）。 */
  public List<Geometry> readFeatures(File file) throws Exception {
    List<Geometry> geoms = new ArrayList<Geometry>();
    try (Reader r = new InputStreamReader(new FileInputStream(file), StandardCharsets.UTF_8)) {
      GeoJsonReader reader = READERS.get();
      // GeoJsonReader.read(Reader) 一次读一个完整 GeoJSON 对象；
      // 生产环境建议用流式 JSON 解析器（如 Jackson）拆分 Feature 后逐个 read(String)。
      Geometry g = reader.read(r);
      for (int i = 0; i < g.getNumGeometries(); i++) {
        geoms.add(g.getGeometryN(i));
      }
    }
    return geoms;
  }

  public static void main(String[] args) throws Exception {
    Case06GeoJsonETL app = new Case06GeoJsonETL();

    // ---- 构造 4 个行政区（2x2 网格） ----
    List<Geometry> zones = new ArrayList<Geometry>();
    List<String> names = new ArrayList<String>();
    for (int x = 0; x < 2; x++) {
      for (int y = 0; y < 2; y++) {
        double minX = x * 10, minY = y * 10;
        Polygon poly = app.gf.toGeometry(
            new org.locationtech.jts.geom.Envelope(minX, minX + 10, minY, minY + 10));
        zones.add(poly);
        names.add("Z" + x + y);
      }
    }
    app.loadZones(zones, names);

    // ---- 构造 100 个测试点 ----
    List<Point> points = new ArrayList<Point>();
    for (int i = 0; i < 100; i++) {
      points.add(app.gf.createPoint(new Coordinate(i % 20, i / 5)));
    }
    System.out.println("点数 = " + points.size());

    // ---- 串行空间连接 ----
    long t0 = System.nanoTime();
    Map<String, Integer> counts = new HashMap<String, Integer>();
    int miss = 0;
    for (Point p : points) {
      String zone = app.locate(p);
      if (zone == null) {
        miss++;
      } else {
        Integer c = counts.get(zone);
        counts.put(zone, c == null ? 1 : c + 1);
      }
    }
    long t1 = System.nanoTime();
    System.out.printf("串行耗时 = %.2f ms%n", (t1 - t0) / 1e6);
    System.out.println("命中分布 = " + counts);
    System.out.println("未命中 = " + miss);

    // ---- 并行空间连接（按索引分片，结果保持顺序） ----
    long t2 = System.nanoTime();
    List<String> results = IntStream.range(0, points.size())
        .parallel()
        .mapToObj(i -> app.locate(points.get(i)))
        .collect(Collectors.toList());
    long t3 = System.nanoTime();
    System.out.printf("并行耗时 = %.2f ms%n", (t3 - t2) / 1e6);

    // ---- 一致性校验：并行结果必须与串行一致 ----
    List<String> serial = new ArrayList<String>();
    for (Point p : points) serial.add(app.locate(p));
    System.out.println("并行与串行结果一致 = " + serial.equals(results));

    // ---- 数据质量问题：重复点/空点 ----
    Set<String> unique = new HashSet<String>();
    int dup = 0;
    for (Point p : points) {
      if (!unique.add(p.toText())) dup++;
    }
    System.out.println("重复点 = " + dup);
  }
}

// 用 PreparedGeometry 加速“多对一”查询
PreparedGeometryFactory pgf = new PreparedGeometryFactory();
List<PreparedGeometry> prepared = new ArrayList<PreparedGeometry>();
for (Geometry zone : zones) {
  prepared.add(pgf.create(zone));
}
// 判定：prepared.get(i).contains(point)