1571 lines
55 KiB
Python
1571 lines
55 KiB
Python
#!/usr/bin/env python3
|
||
# -*- coding: utf-8 -*-
|
||
"""
|
||
CM93 -> MySQL Dataset ID=1 实体空间 Geometry 只读审计
|
||
Version: v3.0
|
||
|
||
目的：
|
||
1. 审计已经真正写入 MySQL 的 Point / Line / Area / Sounding / Coverage geometry。
|
||
2. 不再把 raw_geotype=16 的正常 META parent 当作 geometry 错误。
|
||
3. 检查单个 geometry 的“内在有效性”，不要求当前日本局部 Dataset 与未导入海域连续。
|
||
4. 不修改数据库，不修 geometry，不做 PBF。
|
||
|
||
特别边界规则：
|
||
当前 CM93JAPAN 只是 CM93 的日本沿海局部子集，并非完整 CM93。
|
||
因此：
|
||
- 不做跨 Dataset 外边界连续性 FAIL 判定；
|
||
- 不因 Line 在 Cell/Dataset 边缘终止而判坏；
|
||
- 不因缺少未导入邻区对象而判坏；
|
||
- 对靠近 Cell/Dataset 边界的异常只增加 boundary_context 提示，不自动归因于 importer。
|
||
|
||
依赖：
|
||
pip install pymysql shapely
|
||
|
||
示例：
|
||
python cm93_mysql_spatial_audit_v3.py ^
|
||
--host 127.0.0.1 ^
|
||
--port 3306 ^
|
||
--user root ^
|
||
--database cm93_chart ^
|
||
--dataset-id 1
|
||
|
||
密码缺省时交互输入。
|
||
|
||
输出：
|
||
- spatial_audit_summary.txt
|
||
- spatial_audit_report.json
|
||
- feature_geometry_mapping_issues.csv
|
||
- intrinsic_geometry_issues.csv
|
||
- intrinsic_issue_category_counts.csv
|
||
- invalid_geometry_reason_counts.csv
|
||
- geometry_stats_by_kind.csv
|
||
- issue_object_class_counts.csv
|
||
- issue_cell_counts.csv
|
||
- issue_samples.csv
|
||
- sounding_feature_bbox_issues.csv
|
||
- coverage_geometry_issues.csv
|
||
并自动生成同名 zip。
|
||
|
||
安全：
|
||
- SQL 仅包含 SELECT / SHOW / SET TRANSACTION READ ONLY
|
||
- 无 INSERT / UPDATE / DELETE / DDL
|
||
"""
|
||
|
||
from __future__ import annotations
|
||
|
||
import argparse
|
||
import csv
|
||
import getpass
|
||
import json
|
||
import math
|
||
import os
|
||
import sys
|
||
import traceback
|
||
import zipfile
|
||
from collections import Counter, defaultdict
|
||
from dataclasses import dataclass, asdict
|
||
from datetime import datetime
|
||
from pathlib import Path
|
||
from typing import Any, Dict, Iterable, Iterator, List, Optional, Sequence, Tuple
|
||
|
||
try:
|
||
import pymysql
|
||
from pymysql.cursors import DictCursor
|
||
except ImportError:
|
||
print("缺少 pymysql。请执行：pip install pymysql shapely")
|
||
sys.exit(2)
|
||
|
||
try:
|
||
import shapely
|
||
try:
|
||
from shapely import from_wkb as _from_wkb
|
||
def load_wkb(data: bytes):
|
||
return _from_wkb(data)
|
||
except ImportError:
|
||
from shapely import wkb as _wkb
|
||
def load_wkb(data: bytes):
|
||
return _wkb.loads(data)
|
||
try:
|
||
from shapely.validation import explain_validity
|
||
except Exception:
|
||
def explain_validity(g):
|
||
return "invalid"
|
||
except ImportError:
|
||
print("缺少 shapely。请执行：pip install pymysql shapely")
|
||
sys.exit(2)
|
||
|
||
|
||
EXPECTED_DATASET1 = {
|
||
"cells": 784,
|
||
"features": 322465,
|
||
"attributes": 405153,
|
||
"sounding_points": 275700,
|
||
}
|
||
|
||
EXPECTED_GEOMETRY_ISSUE_INTERPRETATION = {
|
||
"raw_geotype_16_expected_parent_no_geometry": 4406,
|
||
}
|
||
|
||
|
||
def qi(name: str) -> str:
|
||
return "`" + name.replace("`", "``") + "`"
|
||
|
||
|
||
def finite(v: Optional[float]) -> bool:
|
||
return v is not None and math.isfinite(float(v))
|
||
|
||
|
||
def safe_float(v: Any) -> Optional[float]:
|
||
if v is None:
|
||
return None
|
||
try:
|
||
x = float(v)
|
||
return x if math.isfinite(x) else None
|
||
except Exception:
|
||
return None
|
||
|
||
|
||
def bbox_tuple(row: Dict[str, Any], prefix: str = "") -> Optional[Tuple[float, float, float, float]]:
|
||
vals = [
|
||
safe_float(row.get(prefix + "bbox_min_lon")),
|
||
safe_float(row.get(prefix + "bbox_min_lat")),
|
||
safe_float(row.get(prefix + "bbox_max_lon")),
|
||
safe_float(row.get(prefix + "bbox_max_lat")),
|
||
]
|
||
if any(v is None for v in vals):
|
||
return None
|
||
return vals[0], vals[1], vals[2], vals[3]
|
||
|
||
|
||
def bbox_close(a: Tuple[float, float, float, float],
|
||
b: Tuple[float, float, float, float],
|
||
tol: float) -> bool:
|
||
return all(abs(float(x) - float(y)) <= tol for x, y in zip(a, b))
|
||
|
||
|
||
def touches_bbox_edge(inner: Optional[Tuple[float, float, float, float]],
|
||
outer: Optional[Tuple[float, float, float, float]],
|
||
tol: float) -> bool:
|
||
if not inner or not outer:
|
||
return False
|
||
return (
|
||
abs(inner[0] - outer[0]) <= tol or
|
||
abs(inner[1] - outer[1]) <= tol or
|
||
abs(inner[2] - outer[2]) <= tol or
|
||
abs(inner[3] - outer[3]) <= tol
|
||
)
|
||
|
||
|
||
def bbox_in_lonlat_range(b: Optional[Tuple[float, float, float, float]]) -> bool:
|
||
if not b:
|
||
return False
|
||
minx, miny, maxx, maxy = b
|
||
return (
|
||
-180.0 <= minx <= 180.0 and
|
||
-180.0 <= maxx <= 180.0 and
|
||
-90.0 <= miny <= 90.0 and
|
||
-90.0 <= maxy <= 90.0 and
|
||
minx <= maxx and miny <= maxy
|
||
)
|
||
|
||
|
||
def write_csv(path: Path, rows: Iterable[Dict[str, Any]], fieldnames: Optional[List[str]] = None):
|
||
rows = list(rows)
|
||
path.parent.mkdir(parents=True, exist_ok=True)
|
||
if fieldnames is None:
|
||
seen = set()
|
||
fieldnames = []
|
||
for r in rows:
|
||
for k in r.keys():
|
||
if k not in seen:
|
||
seen.add(k)
|
||
fieldnames.append(k)
|
||
with path.open("w", encoding="utf-8-sig", newline="") as f:
|
||
w = csv.DictWriter(f, fieldnames=fieldnames or [], extrasaction="ignore")
|
||
if fieldnames:
|
||
w.writeheader()
|
||
for r in rows:
|
||
w.writerow({k: "" if r.get(k) is None else r.get(k) for k in fieldnames})
|
||
|
||
|
||
def scalar(cur, sql: str, params: Sequence[Any] = ()) -> Any:
|
||
cur.execute(sql, params)
|
||
row = cur.fetchone()
|
||
if not row:
|
||
return None
|
||
return next(iter(row.values()))
|
||
|
||
|
||
def rows(cur, sql: str, params: Sequence[Any] = ()) -> List[Dict[str, Any]]:
|
||
cur.execute(sql, params)
|
||
return list(cur.fetchall())
|
||
|
||
|
||
def begin_read_only(cur, warnings: List[str]):
|
||
try:
|
||
cur.execute("SET SESSION TRANSACTION READ ONLY")
|
||
except Exception as e:
|
||
warnings.append(f"SET SESSION TRANSACTION READ ONLY 不支持/失败：{e}")
|
||
try:
|
||
cur.execute("START TRANSACTION READ ONLY")
|
||
except Exception as e:
|
||
warnings.append(f"START TRANSACTION READ ONLY 不支持/失败：{e}")
|
||
try:
|
||
cur.execute("START TRANSACTION")
|
||
except Exception:
|
||
pass
|
||
|
||
|
||
def probe_spatial_functions(cur) -> Dict[str, Optional[str]]:
|
||
"""
|
||
探测把 MySQL/MariaDB geometry 转成标准 WKB 的函数和 SRID 函数。
|
||
不使用 ST_IsValid，实际 validity 统一交给 Shapely，避免 DB 方言差异。
|
||
"""
|
||
result = {"wkb_fn": None, "srid_fn": None}
|
||
|
||
for fn in ["ST_AsBinary", "ST_AsWKB", "AsBinary", "AsWKB"]:
|
||
try:
|
||
cur.execute(
|
||
f"SELECT {fn}(geom) AS x "
|
||
f"FROM chart_point_geometry LIMIT 1"
|
||
)
|
||
cur.fetchone()
|
||
result["wkb_fn"] = fn
|
||
break
|
||
except Exception:
|
||
continue
|
||
|
||
for fn in ["ST_SRID", "SRID"]:
|
||
try:
|
||
cur.execute(
|
||
f"SELECT {fn}(geom) AS x "
|
||
f"FROM chart_point_geometry LIMIT 1"
|
||
)
|
||
cur.fetchone()
|
||
result["srid_fn"] = fn
|
||
break
|
||
except Exception:
|
||
continue
|
||
|
||
return result
|
||
|
||
|
||
def classify_invalid_reason(reason: str) -> str:
|
||
s = (reason or "").lower()
|
||
if "self-intersection" in s:
|
||
if "ring self-intersection" in s:
|
||
return "RING_SELF_INTERSECTION"
|
||
return "SELF_INTERSECTION"
|
||
if "too few points" in s:
|
||
return "TOO_FEW_POINTS"
|
||
if "hole lies outside shell" in s:
|
||
return "HOLE_OUTSIDE_SHELL"
|
||
if "nested shells" in s:
|
||
return "NESTED_SHELLS"
|
||
if "nested holes" in s:
|
||
return "NESTED_HOLES"
|
||
if "disconnected interior" in s:
|
||
return "DISCONNECTED_INTERIOR"
|
||
if "invalid coordinate" in s or "nan" in s or "inf" in s:
|
||
return "INVALID_COORDINATE"
|
||
if "ring" in s and "closed" in s:
|
||
return "RING_NOT_CLOSED"
|
||
return "OTHER_INVALID_GEOMETRY"
|
||
|
||
|
||
def actual_type_ok(kind: str, geom_type: str) -> bool:
|
||
if kind == "POINT":
|
||
return geom_type == "Point"
|
||
if kind == "LINE":
|
||
return geom_type in {"LineString", "MultiLineString"}
|
||
if kind == "AREA":
|
||
return geom_type in {"Polygon", "MultiPolygon"}
|
||
if kind == "SOUNDING":
|
||
return geom_type == "Point"
|
||
if kind == "COVERAGE":
|
||
return geom_type in {"Polygon", "MultiPolygon"}
|
||
return True
|
||
|
||
|
||
def intrinsic_check(
|
||
*,
|
||
kind: str,
|
||
feature_id: Optional[int],
|
||
object_class: str,
|
||
cell_id: Optional[int],
|
||
cell_name: str,
|
||
source_record_index: Optional[int],
|
||
wkb_data: Any,
|
||
srid: Any,
|
||
normalized_srid: Any,
|
||
feature_bbox: Optional[Tuple[float, float, float, float]],
|
||
cell_bbox: Optional[Tuple[float, float, float, float]],
|
||
dataset_bbox: Optional[Tuple[float, float, float, float]],
|
||
bbox_tolerance: float,
|
||
boundary_tolerance: float,
|
||
extra_id: Optional[int] = None,
|
||
) -> Tuple[Dict[str, Any], List[Dict[str, Any]]]:
|
||
"""
|
||
返回：
|
||
stat: 单 geometry 统计
|
||
issues: 0..N 条 issue
|
||
"""
|
||
base = {
|
||
"kind": kind,
|
||
"feature_id": feature_id,
|
||
"extra_id": extra_id,
|
||
"object_class": object_class,
|
||
"cell_id": cell_id,
|
||
"cell_name": cell_name,
|
||
"source_record_index": source_record_index,
|
||
"srid": srid,
|
||
"normalized_srid": normalized_srid,
|
||
}
|
||
|
||
stat = dict(base)
|
||
stat.update({
|
||
"parsed": 0,
|
||
"empty": 0,
|
||
"valid": None,
|
||
"geom_type": None,
|
||
"simple": None,
|
||
"bounds": None,
|
||
"touches_cell_boundary": 0,
|
||
"touches_dataset_bbox": 0,
|
||
})
|
||
issues: List[Dict[str, Any]] = []
|
||
|
||
if wkb_data is None:
|
||
x = dict(base)
|
||
x.update({
|
||
"issue_category": "WKB_MISSING",
|
||
"issue_detail": "DB geometry exists but WKB conversion returned NULL",
|
||
"invalid_reason": "",
|
||
"boundary_context": "",
|
||
})
|
||
issues.append(x)
|
||
return stat, issues
|
||
|
||
if isinstance(wkb_data, memoryview):
|
||
wkb_data = wkb_data.tobytes()
|
||
if isinstance(wkb_data, bytearray):
|
||
wkb_data = bytes(wkb_data)
|
||
|
||
try:
|
||
geom = load_wkb(wkb_data)
|
||
stat["parsed"] = 1
|
||
except Exception as e:
|
||
x = dict(base)
|
||
x.update({
|
||
"issue_category": "WKB_PARSE_ERROR",
|
||
"issue_detail": f"{type(e).__name__}: {e}",
|
||
"invalid_reason": "",
|
||
"boundary_context": "",
|
||
})
|
||
issues.append(x)
|
||
return stat, issues
|
||
|
||
geom_type = getattr(geom, "geom_type", None)
|
||
stat["geom_type"] = geom_type
|
||
|
||
try:
|
||
is_empty = bool(geom.is_empty)
|
||
except Exception:
|
||
is_empty = False
|
||
stat["empty"] = int(is_empty)
|
||
|
||
if is_empty:
|
||
x = dict(base)
|
||
x.update({
|
||
"issue_category": "EMPTY_GEOMETRY",
|
||
"issue_detail": f"{kind} geometry is EMPTY",
|
||
"invalid_reason": "",
|
||
"actual_geom_type": geom_type,
|
||
"boundary_context": "",
|
||
})
|
||
issues.append(x)
|
||
return stat, issues
|
||
|
||
try:
|
||
b = tuple(float(v) for v in geom.bounds)
|
||
if len(b) == 4 and all(math.isfinite(v) for v in b):
|
||
stat["bounds"] = b
|
||
else:
|
||
b = None
|
||
except Exception:
|
||
b = None
|
||
|
||
stat["touches_cell_boundary"] = int(touches_bbox_edge(b, cell_bbox, boundary_tolerance))
|
||
stat["touches_dataset_bbox"] = int(touches_bbox_edge(b, dataset_bbox, boundary_tolerance))
|
||
|
||
boundary_bits = []
|
||
if stat["touches_cell_boundary"]:
|
||
boundary_bits.append("CELL_EDGE")
|
||
if stat["touches_dataset_bbox"]:
|
||
boundary_bits.append("DATASET_BBOX_EDGE")
|
||
boundary_context = "|".join(boundary_bits)
|
||
|
||
if not actual_type_ok(kind, geom_type or ""):
|
||
x = dict(base)
|
||
x.update({
|
||
"issue_category": "GEOMETRY_TYPE_MISMATCH",
|
||
"issue_detail": f"expected={kind}, actual={geom_type}",
|
||
"invalid_reason": "",
|
||
"actual_geom_type": geom_type,
|
||
"boundary_context": boundary_context,
|
||
})
|
||
issues.append(x)
|
||
|
||
# 坐标范围：CM93 normalized SRID 当前应为经纬度数据。
|
||
if b is None:
|
||
x = dict(base)
|
||
x.update({
|
||
"issue_category": "NON_FINITE_OR_UNREADABLE_BOUNDS",
|
||
"issue_detail": "geometry bounds unavailable/non-finite",
|
||
"invalid_reason": "",
|
||
"actual_geom_type": geom_type,
|
||
"boundary_context": boundary_context,
|
||
})
|
||
issues.append(x)
|
||
elif not bbox_in_lonlat_range(b):
|
||
x = dict(base)
|
||
x.update({
|
||
"issue_category": "OUT_OF_RANGE_COORDINATE",
|
||
"issue_detail": f"bounds={b}",
|
||
"invalid_reason": "",
|
||
"actual_geom_type": geom_type,
|
||
"boundary_context": boundary_context,
|
||
})
|
||
issues.append(x)
|
||
|
||
# SRID 一致性。无法探测 SRID 函数时 srid=None，不判 FAIL。
|
||
if srid is not None and normalized_srid is not None:
|
||
try:
|
||
if int(srid) != int(normalized_srid):
|
||
x = dict(base)
|
||
x.update({
|
||
"issue_category": "SRID_MISMATCH",
|
||
"issue_detail": f"geom_srid={srid}, feature_normalized_srid={normalized_srid}",
|
||
"invalid_reason": "",
|
||
"actual_geom_type": geom_type,
|
||
"boundary_context": boundary_context,
|
||
})
|
||
issues.append(x)
|
||
except Exception:
|
||
pass
|
||
|
||
# Feature bbox 对账。Sounding 在后面按 feature 聚合单独核对。
|
||
if kind in {"POINT", "LINE", "AREA"}:
|
||
if feature_bbox is None:
|
||
x = dict(base)
|
||
x.update({
|
||
"issue_category": "FEATURE_BBOX_MISSING",
|
||
"issue_detail": f"{kind} has geometry but feature bbox is NULL/incomplete",
|
||
"invalid_reason": "",
|
||
"actual_geom_type": geom_type,
|
||
"geometry_bounds": b,
|
||
"boundary_context": boundary_context,
|
||
})
|
||
issues.append(x)
|
||
elif b is not None and not bbox_close(b, feature_bbox, bbox_tolerance):
|
||
x = dict(base)
|
||
x.update({
|
||
"issue_category": "FEATURE_BBOX_MISMATCH",
|
||
"issue_detail": f"feature_bbox={feature_bbox}, geometry_bounds={b}, tol={bbox_tolerance}",
|
||
"invalid_reason": "",
|
||
"actual_geom_type": geom_type,
|
||
"geometry_bounds": b,
|
||
"feature_bbox": feature_bbox,
|
||
"boundary_context": boundary_context,
|
||
})
|
||
issues.append(x)
|
||
|
||
# Shapely validity：
|
||
# - Area/Coverage 的 invalid 是核心检查；
|
||
# - Point/Line 也记录 validity；
|
||
# - Line 非 simple 不等同 invalid，只作为统计，不生成 issue。
|
||
try:
|
||
valid = bool(geom.is_valid)
|
||
stat["valid"] = int(valid)
|
||
except Exception:
|
||
valid = True
|
||
stat["valid"] = None
|
||
|
||
try:
|
||
stat["simple"] = int(bool(geom.is_simple))
|
||
except Exception:
|
||
stat["simple"] = None
|
||
|
||
if not valid:
|
||
try:
|
||
reason = explain_validity(geom)
|
||
except Exception as e:
|
||
reason = f"validity explanation unavailable: {e}"
|
||
reason_class = classify_invalid_reason(reason)
|
||
x = dict(base)
|
||
x.update({
|
||
"issue_category": reason_class,
|
||
"issue_detail": reason,
|
||
"invalid_reason": reason,
|
||
"actual_geom_type": geom_type,
|
||
"geometry_bounds": b,
|
||
"boundary_context": boundary_context,
|
||
})
|
||
issues.append(x)
|
||
|
||
return stat, issues
|
||
|
||
|
||
def get_dataset_bbox(cur, dataset_id: int) -> Optional[Tuple[float, float, float, float]]:
|
||
cur.execute("""
|
||
SELECT
|
||
MIN(bbox_min_lon) AS minx,
|
||
MIN(bbox_min_lat) AS miny,
|
||
MAX(bbox_max_lon) AS maxx,
|
||
MAX(bbox_max_lat) AS maxy
|
||
FROM chart_cell
|
||
WHERE dataset_id=%s
|
||
AND bbox_min_lon IS NOT NULL
|
||
AND bbox_min_lat IS NOT NULL
|
||
AND bbox_max_lon IS NOT NULL
|
||
AND bbox_max_lat IS NOT NULL
|
||
""", (dataset_id,))
|
||
r = cur.fetchone() or {}
|
||
vals = [safe_float(r.get("minx")), safe_float(r.get("miny")),
|
||
safe_float(r.get("maxx")), safe_float(r.get("maxy"))]
|
||
if any(v is None for v in vals):
|
||
return None
|
||
return vals[0], vals[1], vals[2], vals[3]
|
||
|
||
|
||
def audit_mapping(cur, dataset_id: int) -> Tuple[List[Dict[str, Any]], Dict[str, Any]]:
|
||
"""
|
||
检查 feature.geometry_type 与实际存储表是否一致。
|
||
"""
|
||
mapping_issues: List[Dict[str, Any]] = []
|
||
stats: Dict[str, Any] = {}
|
||
|
||
# Feature type counts
|
||
ft = rows(cur, """
|
||
SELECT geometry_type, COUNT(*) AS n
|
||
FROM chart_feature
|
||
WHERE dataset_id=%s
|
||
GROUP BY geometry_type
|
||
ORDER BY geometry_type
|
||
""", (dataset_id,))
|
||
stats["feature_type_counts"] = ft
|
||
|
||
# 实体表 counts（只统计 Dataset 1）
|
||
for name, table in [
|
||
("point_rows", "chart_point_geometry"),
|
||
("line_rows", "chart_line_geometry"),
|
||
("area_rows", "chart_area_geometry"),
|
||
]:
|
||
stats[name] = int(scalar(cur, f"""
|
||
SELECT COUNT(*)
|
||
FROM {table} g
|
||
JOIN chart_feature f ON f.id=g.feature_id
|
||
WHERE f.dataset_id=%s
|
||
""", (dataset_id,)) or 0)
|
||
|
||
stats["sounding_rows"] = int(scalar(cur, """
|
||
SELECT COUNT(*)
|
||
FROM chart_sounding_point s
|
||
JOIN chart_feature f ON f.id=s.feature_id
|
||
WHERE f.dataset_id=%s
|
||
""", (dataset_id,)) or 0)
|
||
|
||
# 缺 expected geometry
|
||
checks = [
|
||
("POINT_MISSING_GEOMETRY", "POINT", "chart_point_geometry"),
|
||
("LINE_MISSING_GEOMETRY", "LINE", "chart_line_geometry"),
|
||
("AREA_MISSING_GEOMETRY", "AREA", "chart_area_geometry"),
|
||
]
|
||
for issue_cat, gtype, table in checks:
|
||
rs = rows(cur, f"""
|
||
SELECT
|
||
f.id AS feature_id, f.object_class, f.cell_id, c.cell_name,
|
||
f.source_record_index, f.raw_geotype, f.geometry_type
|
||
FROM chart_feature f
|
||
JOIN chart_cell c ON c.id=f.cell_id
|
||
LEFT JOIN {table} g ON g.feature_id=f.id
|
||
WHERE f.dataset_id=%s
|
||
AND f.geometry_type=%s
|
||
AND g.feature_id IS NULL
|
||
ORDER BY f.id
|
||
""", (dataset_id, gtype))
|
||
for r in rs:
|
||
r["issue_category"] = issue_cat
|
||
mapping_issues.append(r)
|
||
|
||
# SOUNDING feature 没有任何 sounding point
|
||
rs = rows(cur, """
|
||
SELECT
|
||
f.id AS feature_id, f.object_class, f.cell_id, c.cell_name,
|
||
f.source_record_index, f.raw_geotype, f.geometry_type
|
||
FROM chart_feature f
|
||
JOIN chart_cell c ON c.id=f.cell_id
|
||
LEFT JOIN chart_sounding_point s ON s.feature_id=f.id
|
||
WHERE f.dataset_id=%s
|
||
AND f.geometry_type='SOUNDING'
|
||
GROUP BY f.id, f.object_class, f.cell_id, c.cell_name,
|
||
f.source_record_index, f.raw_geotype, f.geometry_type
|
||
HAVING COUNT(s.id)=0
|
||
ORDER BY f.id
|
||
""", (dataset_id,))
|
||
for r in rs:
|
||
r["issue_category"] = "SOUNDING_MISSING_POINTS"
|
||
mapping_issues.append(r)
|
||
|
||
# geometry 表串类型
|
||
for issue_cat, table, expected in [
|
||
("POINT_TABLE_WRONG_FEATURE_TYPE", "chart_point_geometry", "POINT"),
|
||
("LINE_TABLE_WRONG_FEATURE_TYPE", "chart_line_geometry", "LINE"),
|
||
("AREA_TABLE_WRONG_FEATURE_TYPE", "chart_area_geometry", "AREA"),
|
||
]:
|
||
rs = rows(cur, f"""
|
||
SELECT
|
||
f.id AS feature_id, f.object_class, f.cell_id, c.cell_name,
|
||
f.source_record_index, f.raw_geotype, f.geometry_type
|
||
FROM {table} g
|
||
JOIN chart_feature f ON f.id=g.feature_id
|
||
JOIN chart_cell c ON c.id=f.cell_id
|
||
WHERE f.dataset_id=%s
|
||
AND f.geometry_type<>%s
|
||
ORDER BY f.id
|
||
""", (dataset_id, expected))
|
||
for r in rs:
|
||
r["issue_category"] = issue_cat
|
||
mapping_issues.append(r)
|
||
|
||
# sounding point 串到非 SOUNDING feature
|
||
rs = rows(cur, """
|
||
SELECT DISTINCT
|
||
f.id AS feature_id, f.object_class, f.cell_id, c.cell_name,
|
||
f.source_record_index, f.raw_geotype, f.geometry_type
|
||
FROM chart_sounding_point s
|
||
JOIN chart_feature f ON f.id=s.feature_id
|
||
JOIN chart_cell c ON c.id=f.cell_id
|
||
WHERE f.dataset_id=%s
|
||
AND f.geometry_type<>'SOUNDING'
|
||
ORDER BY f.id
|
||
""", (dataset_id,))
|
||
for r in rs:
|
||
r["issue_category"] = "SOUNDING_POINT_WRONG_FEATURE_TYPE"
|
||
mapping_issues.append(r)
|
||
|
||
# META/UNKNOWN 不应偷偷挂实际 geometry
|
||
rs = rows(cur, """
|
||
SELECT
|
||
f.id AS feature_id, f.object_class, f.cell_id, c.cell_name,
|
||
f.source_record_index, f.raw_geotype, f.geometry_type,
|
||
CASE WHEN pg.feature_id IS NOT NULL THEN 1 ELSE 0 END AS has_point,
|
||
CASE WHEN lg.feature_id IS NOT NULL THEN 1 ELSE 0 END AS has_line,
|
||
CASE WHEN ag.feature_id IS NOT NULL THEN 1 ELSE 0 END AS has_area,
|
||
COUNT(sp.id) AS sounding_point_count
|
||
FROM chart_feature f
|
||
JOIN chart_cell c ON c.id=f.cell_id
|
||
LEFT JOIN chart_point_geometry pg ON pg.feature_id=f.id
|
||
LEFT JOIN chart_line_geometry lg ON lg.feature_id=f.id
|
||
LEFT JOIN chart_area_geometry ag ON ag.feature_id=f.id
|
||
LEFT JOIN chart_sounding_point sp ON sp.feature_id=f.id
|
||
WHERE f.dataset_id=%s
|
||
AND f.geometry_type IN ('META','UNKNOWN')
|
||
GROUP BY f.id, f.object_class, f.cell_id, c.cell_name,
|
||
f.source_record_index, f.raw_geotype, f.geometry_type,
|
||
pg.feature_id, lg.feature_id, ag.feature_id
|
||
HAVING has_point=1 OR has_line=1 OR has_area=1 OR sounding_point_count>0
|
||
ORDER BY f.id
|
||
""", (dataset_id,))
|
||
for r in rs:
|
||
r["issue_category"] = "META_OR_UNKNOWN_HAS_UNEXPECTED_GEOMETRY"
|
||
mapping_issues.append(r)
|
||
|
||
stats["mapping_issue_count"] = len(mapping_issues)
|
||
stats["mapping_issue_categories"] = dict(Counter(r["issue_category"] for r in mapping_issues))
|
||
return mapping_issues, stats
|
||
|
||
|
||
def scan_feature_geometry_table(
|
||
cur,
|
||
*,
|
||
table: str,
|
||
kind: str,
|
||
dataset_id: int,
|
||
wkb_fn: str,
|
||
srid_fn: Optional[str],
|
||
dataset_bbox: Optional[Tuple[float, float, float, float]],
|
||
chunk_size: int,
|
||
bbox_tolerance: float,
|
||
boundary_tolerance: float,
|
||
) -> Tuple[Dict[str, Any], List[Dict[str, Any]]]:
|
||
last_id = 0
|
||
stats = Counter()
|
||
geom_types = Counter()
|
||
invalid_reasons = Counter()
|
||
issues: List[Dict[str, Any]] = []
|
||
|
||
srid_expr = f"{srid_fn}(g.geom)" if srid_fn else "NULL"
|
||
|
||
while True:
|
||
sql = f"""
|
||
SELECT
|
||
g.feature_id,
|
||
{wkb_fn}(g.geom) AS geom_wkb,
|
||
{srid_expr} AS geom_srid,
|
||
f.object_class,
|
||
f.cell_id,
|
||
c.cell_name,
|
||
f.source_record_index,
|
||
f.normalized_srid,
|
||
f.bbox_min_lon, f.bbox_min_lat, f.bbox_max_lon, f.bbox_max_lat,
|
||
c.bbox_min_lon AS cell_bbox_min_lon,
|
||
c.bbox_min_lat AS cell_bbox_min_lat,
|
||
c.bbox_max_lon AS cell_bbox_max_lon,
|
||
c.bbox_max_lat AS cell_bbox_max_lat
|
||
FROM {table} g
|
||
JOIN chart_feature f ON f.id=g.feature_id
|
||
JOIN chart_cell c ON c.id=f.cell_id
|
||
WHERE f.dataset_id=%s
|
||
AND g.feature_id>%s
|
||
ORDER BY g.feature_id
|
||
LIMIT %s
|
||
"""
|
||
cur.execute(sql, (dataset_id, last_id, chunk_size))
|
||
batch = list(cur.fetchall())
|
||
if not batch:
|
||
break
|
||
|
||
for r in batch:
|
||
last_id = int(r["feature_id"])
|
||
fb = bbox_tuple(r)
|
||
cb = bbox_tuple(r, "cell_")
|
||
|
||
stat, found = intrinsic_check(
|
||
kind=kind,
|
||
feature_id=r["feature_id"],
|
||
object_class=r.get("object_class") or "",
|
||
cell_id=r.get("cell_id"),
|
||
cell_name=r.get("cell_name") or "",
|
||
source_record_index=r.get("source_record_index"),
|
||
wkb_data=r.get("geom_wkb"),
|
||
srid=r.get("geom_srid"),
|
||
normalized_srid=r.get("normalized_srid"),
|
||
feature_bbox=fb,
|
||
cell_bbox=cb,
|
||
dataset_bbox=dataset_bbox,
|
||
bbox_tolerance=bbox_tolerance,
|
||
boundary_tolerance=boundary_tolerance,
|
||
)
|
||
|
||
stats["rows"] += 1
|
||
stats["parsed"] += stat["parsed"]
|
||
stats["empty"] += stat["empty"]
|
||
if stat["valid"] == 1:
|
||
stats["valid"] += 1
|
||
elif stat["valid"] == 0:
|
||
stats["invalid"] += 1
|
||
if kind == "LINE" and stat["simple"] == 0:
|
||
stats["non_simple_line_info"] += 1
|
||
if stat["touches_cell_boundary"]:
|
||
stats["touch_cell_boundary_info"] += 1
|
||
if stat["touches_dataset_bbox"]:
|
||
stats["touch_dataset_bbox_info"] += 1
|
||
geom_types[stat.get("geom_type") or "UNPARSED"] += 1
|
||
|
||
for x in found:
|
||
issues.append(x)
|
||
invalid_reasons[x["issue_category"]] += 1
|
||
|
||
print(
|
||
f" {kind:<8} scanned={stats['rows']:,} "
|
||
f"issues={len(issues):,} last_feature_id={last_id}"
|
||
)
|
||
|
||
return {
|
||
"kind": kind,
|
||
"table": table,
|
||
"rows": stats["rows"],
|
||
"parsed": stats["parsed"],
|
||
"empty": stats["empty"],
|
||
"valid": stats["valid"],
|
||
"invalid": stats["invalid"],
|
||
"non_simple_line_info": stats["non_simple_line_info"],
|
||
"touch_cell_boundary_info": stats["touch_cell_boundary_info"],
|
||
"touch_dataset_bbox_info": stats["touch_dataset_bbox_info"],
|
||
"actual_geom_types": dict(geom_types),
|
||
"issue_count": len(issues),
|
||
"issue_categories": dict(invalid_reasons),
|
||
}, issues
|
||
|
||
|
||
def scan_soundings(
|
||
cur,
|
||
*,
|
||
dataset_id: int,
|
||
wkb_fn: str,
|
||
srid_fn: Optional[str],
|
||
dataset_bbox: Optional[Tuple[float, float, float, float]],
|
||
chunk_size: int,
|
||
bbox_tolerance: float,
|
||
boundary_tolerance: float,
|
||
) -> Tuple[Dict[str, Any], List[Dict[str, Any]], List[Dict[str, Any]]]:
|
||
last_id = 0
|
||
stats = Counter()
|
||
geom_types = Counter()
|
||
issues: List[Dict[str, Any]] = []
|
||
|
||
# feature_id -> [minx,miny,maxx,maxy,count]
|
||
feature_bounds: Dict[int, List[float]] = {}
|
||
|
||
srid_expr = f"{srid_fn}(s.geom)" if srid_fn else "NULL"
|
||
|
||
while True:
|
||
sql = f"""
|
||
SELECT
|
||
s.id AS sounding_id,
|
||
s.feature_id,
|
||
s.depth,
|
||
{wkb_fn}(s.geom) AS geom_wkb,
|
||
{srid_expr} AS geom_srid,
|
||
f.object_class,
|
||
f.cell_id,
|
||
c.cell_name,
|
||
f.source_record_index,
|
||
f.normalized_srid,
|
||
f.bbox_min_lon, f.bbox_min_lat, f.bbox_max_lon, f.bbox_max_lat,
|
||
c.bbox_min_lon AS cell_bbox_min_lon,
|
||
c.bbox_min_lat AS cell_bbox_min_lat,
|
||
c.bbox_max_lon AS cell_bbox_max_lon,
|
||
c.bbox_max_lat AS cell_bbox_max_lat
|
||
FROM chart_sounding_point s
|
||
JOIN chart_feature f ON f.id=s.feature_id
|
||
JOIN chart_cell c ON c.id=f.cell_id
|
||
WHERE f.dataset_id=%s
|
||
AND s.id>%s
|
||
ORDER BY s.id
|
||
LIMIT %s
|
||
"""
|
||
cur.execute(sql, (dataset_id, last_id, chunk_size))
|
||
batch = list(cur.fetchall())
|
||
if not batch:
|
||
break
|
||
|
||
for r in batch:
|
||
last_id = int(r["sounding_id"])
|
||
fb = bbox_tuple(r)
|
||
cb = bbox_tuple(r, "cell_")
|
||
|
||
stat, found = intrinsic_check(
|
||
kind="SOUNDING",
|
||
feature_id=r["feature_id"],
|
||
extra_id=r["sounding_id"],
|
||
object_class=r.get("object_class") or "",
|
||
cell_id=r.get("cell_id"),
|
||
cell_name=r.get("cell_name") or "",
|
||
source_record_index=r.get("source_record_index"),
|
||
wkb_data=r.get("geom_wkb"),
|
||
srid=r.get("geom_srid"),
|
||
normalized_srid=r.get("normalized_srid"),
|
||
feature_bbox=fb,
|
||
cell_bbox=cb,
|
||
dataset_bbox=dataset_bbox,
|
||
bbox_tolerance=bbox_tolerance,
|
||
boundary_tolerance=boundary_tolerance,
|
||
)
|
||
|
||
stats["rows"] += 1
|
||
stats["parsed"] += stat["parsed"]
|
||
stats["empty"] += stat["empty"]
|
||
if stat["valid"] == 1:
|
||
stats["valid"] += 1
|
||
elif stat["valid"] == 0:
|
||
stats["invalid"] += 1
|
||
geom_types[stat.get("geom_type") or "UNPARSED"] += 1
|
||
|
||
# depth 必须 finite
|
||
depth = safe_float(r.get("depth"))
|
||
if depth is None:
|
||
issues.append({
|
||
"kind": "SOUNDING",
|
||
"feature_id": r["feature_id"],
|
||
"extra_id": r["sounding_id"],
|
||
"object_class": r.get("object_class") or "",
|
||
"cell_id": r.get("cell_id"),
|
||
"cell_name": r.get("cell_name") or "",
|
||
"source_record_index": r.get("source_record_index"),
|
||
"issue_category": "INVALID_SOUNDING_DEPTH",
|
||
"issue_detail": f"depth={r.get('depth')}",
|
||
"boundary_context": "",
|
||
})
|
||
|
||
issues.extend(found)
|
||
|
||
b = stat.get("bounds")
|
||
if b and len(b) == 4:
|
||
fid = int(r["feature_id"])
|
||
if fid not in feature_bounds:
|
||
feature_bounds[fid] = [b[0], b[1], b[2], b[3], 1]
|
||
else:
|
||
a = feature_bounds[fid]
|
||
a[0] = min(a[0], b[0])
|
||
a[1] = min(a[1], b[1])
|
||
a[2] = max(a[2], b[2])
|
||
a[3] = max(a[3], b[3])
|
||
a[4] += 1
|
||
|
||
print(
|
||
f" SOUNDING scanned={stats['rows']:,} "
|
||
f"issues={len(issues):,} last_sounding_id={last_id}"
|
||
)
|
||
|
||
# 按 feature 聚合 sounding bbox 与 chart_feature bbox 对账。
|
||
bbox_issues: List[Dict[str, Any]] = []
|
||
feature_ids = list(feature_bounds.keys())
|
||
# 704 个 feature 左右，可直接分批 IN
|
||
for pos in range(0, len(feature_ids), 500):
|
||
ids = feature_ids[pos:pos+500]
|
||
placeholders = ",".join(["%s"] * len(ids))
|
||
rs = rows(cur, f"""
|
||
SELECT
|
||
f.id AS feature_id, f.object_class, f.cell_id, c.cell_name,
|
||
f.source_record_index,
|
||
f.bbox_min_lon, f.bbox_min_lat, f.bbox_max_lon, f.bbox_max_lat,
|
||
c.bbox_min_lon AS cell_bbox_min_lon,
|
||
c.bbox_min_lat AS cell_bbox_min_lat,
|
||
c.bbox_max_lon AS cell_bbox_max_lon,
|
||
c.bbox_max_lat AS cell_bbox_max_lat
|
||
FROM chart_feature f
|
||
JOIN chart_cell c ON c.id=f.cell_id
|
||
WHERE f.id IN ({placeholders})
|
||
""", ids)
|
||
for r in rs:
|
||
fid = int(r["feature_id"])
|
||
a = feature_bounds[fid]
|
||
actual = (a[0], a[1], a[2], a[3])
|
||
fb = bbox_tuple(r)
|
||
cb = bbox_tuple(r, "cell_")
|
||
context = []
|
||
if touches_bbox_edge(actual, cb, boundary_tolerance):
|
||
context.append("CELL_EDGE")
|
||
if touches_bbox_edge(actual, dataset_bbox, boundary_tolerance):
|
||
context.append("DATASET_BBOX_EDGE")
|
||
if fb is None:
|
||
bbox_issues.append({
|
||
"feature_id": fid,
|
||
"object_class": r.get("object_class"),
|
||
"cell_id": r.get("cell_id"),
|
||
"cell_name": r.get("cell_name"),
|
||
"source_record_index": r.get("source_record_index"),
|
||
"issue_category": "SOUNDING_FEATURE_BBOX_MISSING",
|
||
"point_count": int(a[4]),
|
||
"actual_bbox": actual,
|
||
"feature_bbox": None,
|
||
"boundary_context": "|".join(context),
|
||
})
|
||
elif not bbox_close(actual, fb, bbox_tolerance):
|
||
bbox_issues.append({
|
||
"feature_id": fid,
|
||
"object_class": r.get("object_class"),
|
||
"cell_id": r.get("cell_id"),
|
||
"cell_name": r.get("cell_name"),
|
||
"source_record_index": r.get("source_record_index"),
|
||
"issue_category": "SOUNDING_FEATURE_BBOX_MISMATCH",
|
||
"point_count": int(a[4]),
|
||
"actual_bbox": actual,
|
||
"feature_bbox": fb,
|
||
"boundary_context": "|".join(context),
|
||
})
|
||
|
||
return {
|
||
"kind": "SOUNDING",
|
||
"table": "chart_sounding_point",
|
||
"rows": stats["rows"],
|
||
"parsed": stats["parsed"],
|
||
"empty": stats["empty"],
|
||
"valid": stats["valid"],
|
||
"invalid": stats["invalid"],
|
||
"actual_geom_types": dict(geom_types),
|
||
"issue_count_excluding_feature_bbox": len(issues),
|
||
"feature_bbox_issue_count": len(bbox_issues),
|
||
}, issues, bbox_issues
|
||
|
||
|
||
def scan_coverage(
|
||
cur,
|
||
*,
|
||
dataset_id: int,
|
||
wkb_fn: str,
|
||
srid_fn: Optional[str],
|
||
dataset_bbox: Optional[Tuple[float, float, float, float]],
|
||
bbox_tolerance: float,
|
||
boundary_tolerance: float,
|
||
) -> Tuple[Dict[str, Any], List[Dict[str, Any]]]:
|
||
srid_expr = f"{srid_fn}(v.geom)" if srid_fn else "NULL"
|
||
cur.execute(f"""
|
||
SELECT
|
||
v.id AS coverage_id,
|
||
v.cell_id,
|
||
c.cell_name,
|
||
v.source_record_index,
|
||
{wkb_fn}(v.geom) AS geom_wkb,
|
||
{srid_expr} AS geom_srid,
|
||
c.bbox_min_lon AS cell_bbox_min_lon,
|
||
c.bbox_min_lat AS cell_bbox_min_lat,
|
||
c.bbox_max_lon AS cell_bbox_max_lon,
|
||
c.bbox_max_lat AS cell_bbox_max_lat
|
||
FROM chart_coverage v
|
||
JOIN chart_cell c ON c.id=v.cell_id
|
||
WHERE v.dataset_id=%s
|
||
ORDER BY v.id
|
||
""", (dataset_id,))
|
||
rs = list(cur.fetchall())
|
||
|
||
stats = Counter()
|
||
geom_types = Counter()
|
||
issues: List[Dict[str, Any]] = []
|
||
for r in rs:
|
||
cb = bbox_tuple(r, "cell_")
|
||
stat, found = intrinsic_check(
|
||
kind="COVERAGE",
|
||
feature_id=None,
|
||
extra_id=r["coverage_id"],
|
||
object_class="COVERAGE",
|
||
cell_id=r.get("cell_id"),
|
||
cell_name=r.get("cell_name") or "",
|
||
source_record_index=r.get("source_record_index"),
|
||
wkb_data=r.get("geom_wkb"),
|
||
srid=r.get("geom_srid"),
|
||
normalized_srid=None,
|
||
feature_bbox=None,
|
||
cell_bbox=cb,
|
||
dataset_bbox=dataset_bbox,
|
||
bbox_tolerance=bbox_tolerance,
|
||
boundary_tolerance=boundary_tolerance,
|
||
)
|
||
stats["rows"] += 1
|
||
stats["parsed"] += stat["parsed"]
|
||
stats["empty"] += stat["empty"]
|
||
if stat["valid"] == 1:
|
||
stats["valid"] += 1
|
||
elif stat["valid"] == 0:
|
||
stats["invalid"] += 1
|
||
geom_types[stat.get("geom_type") or "UNPARSED"] += 1
|
||
issues.extend(found)
|
||
|
||
return {
|
||
"kind": "COVERAGE",
|
||
"table": "chart_coverage",
|
||
"rows": stats["rows"],
|
||
"parsed": stats["parsed"],
|
||
"empty": stats["empty"],
|
||
"valid": stats["valid"],
|
||
"invalid": stats["invalid"],
|
||
"actual_geom_types": dict(geom_types),
|
||
"issue_count": len(issues),
|
||
}, issues
|
||
|
||
|
||
def zip_directory(folder: Path) -> Path:
|
||
zpath = folder.with_suffix(".zip")
|
||
with zipfile.ZipFile(zpath, "w", compression=zipfile.ZIP_DEFLATED) as z:
|
||
for p in sorted(folder.rglob("*")):
|
||
if p.is_file():
|
||
z.write(p, p.relative_to(folder.parent))
|
||
return zpath
|
||
|
||
|
||
def parse_args():
|
||
p = argparse.ArgumentParser(
|
||
description="CM93 MySQL 实体空间 geometry 只读审计 v3"
|
||
)
|
||
p.add_argument("--host", default=os.getenv("MYSQL_HOST", "127.0.0.1"))
|
||
p.add_argument("--port", type=int, default=int(os.getenv("MYSQL_PORT", "3306")))
|
||
p.add_argument("--user", default=os.getenv("MYSQL_USER", "root"))
|
||
p.add_argument("--password", default=os.getenv("MYSQL_PASSWORD"))
|
||
p.add_argument("--database", default=os.getenv("MYSQL_DATABASE", "cm93_chart"))
|
||
p.add_argument("--dataset-id", type=int, default=1)
|
||
p.add_argument("--chunk-size", type=int, default=2000)
|
||
p.add_argument(
|
||
"--bbox-tolerance",
|
||
type=float,
|
||
default=1e-6,
|
||
help="Feature bbox 与实际 geometry bounds 对账容差（degree），默认 1e-6"
|
||
)
|
||
p.add_argument(
|
||
"--boundary-tolerance",
|
||
type=float,
|
||
default=1e-5,
|
||
help="仅用于标注 geometry 是否碰 Cell/Dataset bbox 边缘，不用于 PASS/FAIL"
|
||
)
|
||
p.add_argument("--output", default=None)
|
||
p.add_argument(
|
||
"--skip-coverage",
|
||
action="store_true",
|
||
help="跳过 chart_coverage geometry 审计"
|
||
)
|
||
return p.parse_args()
|
||
|
||
|
||
def main():
|
||
args = parse_args()
|
||
if not args.password:
|
||
args.password = getpass.getpass("MySQL password: ")
|
||
|
||
stamp = datetime.now().strftime("%Y%m%d_%H%M%S")
|
||
outdir = Path(
|
||
args.output or
|
||
f"cm93_spatial_audit_v3_dataset_{args.dataset_id}_{stamp}"
|
||
)
|
||
outdir.mkdir(parents=True, exist_ok=True)
|
||
|
||
report: Dict[str, Any] = {
|
||
"version": "v3.0",
|
||
"generated_at": datetime.now().isoformat(timespec="seconds"),
|
||
"dataset_id": args.dataset_id,
|
||
"database": args.database,
|
||
"read_only": True,
|
||
"pbf": "OUT_OF_SCOPE",
|
||
"dataset_scope_note": (
|
||
"当前 CM93JAPAN 只是日本沿海局部子集，并非完整 CM93。"
|
||
"不把跨未导入海域/邻区的不连续性作为 geometry corruption。"
|
||
),
|
||
"warnings": [],
|
||
"server": {},
|
||
"dataset_bbox": None,
|
||
"counts": {},
|
||
"mapping": {},
|
||
"geometry_stats": [],
|
||
"issue_summary": {},
|
||
}
|
||
|
||
conn = None
|
||
try:
|
||
conn = pymysql.connect(
|
||
host=args.host,
|
||
port=args.port,
|
||
user=args.user,
|
||
password=args.password,
|
||
database=args.database,
|
||
charset="utf8mb4",
|
||
cursorclass=DictCursor,
|
||
autocommit=False,
|
||
connect_timeout=20,
|
||
read_timeout=600,
|
||
write_timeout=30,
|
||
)
|
||
with conn.cursor() as cur:
|
||
begin_read_only(cur, report["warnings"])
|
||
|
||
print("[1/8] 读取服务器与 Dataset 基线...")
|
||
report["server"]["version"] = scalar(cur, "SELECT VERSION()")
|
||
try:
|
||
report["server"]["version_comment"] = scalar(cur, "SELECT @@version_comment")
|
||
except Exception:
|
||
report["server"]["version_comment"] = None
|
||
|
||
funcs = probe_spatial_functions(cur)
|
||
report["server"]["spatial_functions"] = funcs
|
||
if not funcs["wkb_fn"]:
|
||
raise RuntimeError(
|
||
"无法找到可用的 geometry -> WKB 函数 "
|
||
"(ST_AsBinary/ST_AsWKB/AsBinary/AsWKB)"
|
||
)
|
||
print(f" DB: {report['server']['version']}")
|
||
print(f" WKB function: {funcs['wkb_fn']}")
|
||
print(f" SRID function: {funcs['srid_fn']}")
|
||
|
||
dataset_bbox = get_dataset_bbox(cur, args.dataset_id)
|
||
report["dataset_bbox"] = dataset_bbox
|
||
print(f" Dataset bbox: {dataset_bbox}")
|
||
|
||
# 基线 counts
|
||
report["counts"]["cells"] = int(scalar(cur,
|
||
"SELECT COUNT(*) FROM chart_cell WHERE dataset_id=%s",
|
||
(args.dataset_id,)) or 0)
|
||
report["counts"]["features"] = int(scalar(cur,
|
||
"SELECT COUNT(*) FROM chart_feature WHERE dataset_id=%s",
|
||
(args.dataset_id,)) or 0)
|
||
report["counts"]["attributes"] = int(scalar(cur, """
|
||
SELECT COUNT(*)
|
||
FROM chart_feature_attribute a
|
||
JOIN chart_feature f ON f.id=a.feature_id
|
||
WHERE f.dataset_id=%s
|
||
""", (args.dataset_id,)) or 0)
|
||
report["counts"]["sounding_points"] = int(scalar(cur, """
|
||
SELECT COUNT(*)
|
||
FROM chart_sounding_point s
|
||
JOIN chart_feature f ON f.id=s.feature_id
|
||
WHERE f.dataset_id=%s
|
||
""", (args.dataset_id,)) or 0)
|
||
|
||
report["counts"]["raw_geotype16_meta"] = int(scalar(cur, """
|
||
SELECT COUNT(*)
|
||
FROM chart_feature
|
||
WHERE dataset_id=%s
|
||
AND raw_geotype=16
|
||
AND geometry_type='META'
|
||
""", (args.dataset_id,)) or 0)
|
||
|
||
report["counts"]["raw_geotype16_geometry_issue"] = int(scalar(cur, """
|
||
SELECT COUNT(*)
|
||
FROM chart_feature
|
||
WHERE dataset_id=%s
|
||
AND raw_geotype=16
|
||
AND geometry_type='META'
|
||
AND geometry_issue IS NOT NULL
|
||
""", (args.dataset_id,)) or 0)
|
||
|
||
print("[2/8] 检查 Feature geometry_type ↔ 实际 geometry 表映射...")
|
||
mapping_issues, mapping_stats = audit_mapping(cur, args.dataset_id)
|
||
report["mapping"] = mapping_stats
|
||
write_csv(outdir / "feature_geometry_mapping_issues.csv", mapping_issues)
|
||
print(f" mapping issues: {len(mapping_issues):,}")
|
||
|
||
all_intrinsic_issues: List[Dict[str, Any]] = []
|
||
geometry_stats: List[Dict[str, Any]] = []
|
||
|
||
print("[3/8] 扫描 Point geometry...")
|
||
st, iss = scan_feature_geometry_table(
|
||
cur,
|
||
table="chart_point_geometry",
|
||
kind="POINT",
|
||
dataset_id=args.dataset_id,
|
||
wkb_fn=funcs["wkb_fn"],
|
||
srid_fn=funcs["srid_fn"],
|
||
dataset_bbox=dataset_bbox,
|
||
chunk_size=args.chunk_size,
|
||
bbox_tolerance=args.bbox_tolerance,
|
||
boundary_tolerance=args.boundary_tolerance,
|
||
)
|
||
geometry_stats.append(st)
|
||
all_intrinsic_issues.extend(iss)
|
||
|
||
print("[4/8] 扫描 Line geometry...")
|
||
st, iss = scan_feature_geometry_table(
|
||
cur,
|
||
table="chart_line_geometry",
|
||
kind="LINE",
|
||
dataset_id=args.dataset_id,
|
||
wkb_fn=funcs["wkb_fn"],
|
||
srid_fn=funcs["srid_fn"],
|
||
dataset_bbox=dataset_bbox,
|
||
chunk_size=args.chunk_size,
|
||
bbox_tolerance=args.bbox_tolerance,
|
||
boundary_tolerance=args.boundary_tolerance,
|
||
)
|
||
geometry_stats.append(st)
|
||
all_intrinsic_issues.extend(iss)
|
||
|
||
print("[5/8] 扫描 Area geometry...")
|
||
st, iss = scan_feature_geometry_table(
|
||
cur,
|
||
table="chart_area_geometry",
|
||
kind="AREA",
|
||
dataset_id=args.dataset_id,
|
||
wkb_fn=funcs["wkb_fn"],
|
||
srid_fn=funcs["srid_fn"],
|
||
dataset_bbox=dataset_bbox,
|
||
chunk_size=args.chunk_size,
|
||
bbox_tolerance=args.bbox_tolerance,
|
||
boundary_tolerance=args.boundary_tolerance,
|
||
)
|
||
geometry_stats.append(st)
|
||
all_intrinsic_issues.extend(iss)
|
||
|
||
print("[6/8] 扫描 Sounding points，并按 Feature 对账 bbox...")
|
||
st, iss, sounding_bbox_issues = scan_soundings(
|
||
cur,
|
||
dataset_id=args.dataset_id,
|
||
wkb_fn=funcs["wkb_fn"],
|
||
srid_fn=funcs["srid_fn"],
|
||
dataset_bbox=dataset_bbox,
|
||
chunk_size=max(args.chunk_size, 5000),
|
||
bbox_tolerance=args.bbox_tolerance,
|
||
boundary_tolerance=args.boundary_tolerance,
|
||
)
|
||
geometry_stats.append(st)
|
||
all_intrinsic_issues.extend(iss)
|
||
write_csv(outdir / "sounding_feature_bbox_issues.csv", sounding_bbox_issues)
|
||
|
||
coverage_issues: List[Dict[str, Any]] = []
|
||
if not args.skip_coverage:
|
||
print("[7/8] 扫描 Coverage geometry...")
|
||
st, coverage_issues = scan_coverage(
|
||
cur,
|
||
dataset_id=args.dataset_id,
|
||
wkb_fn=funcs["wkb_fn"],
|
||
srid_fn=funcs["srid_fn"],
|
||
dataset_bbox=dataset_bbox,
|
||
bbox_tolerance=args.bbox_tolerance,
|
||
boundary_tolerance=args.boundary_tolerance,
|
||
)
|
||
geometry_stats.append(st)
|
||
write_csv(outdir / "coverage_geometry_issues.csv", coverage_issues)
|
||
else:
|
||
print("[7/8] Coverage skipped by option.")
|
||
|
||
print("[8/8] 汇总输出...")
|
||
|
||
report["geometry_stats"] = geometry_stats
|
||
|
||
# 主 intrinsic issue 输出（Feature Point/Line/Area + sounding point）
|
||
write_csv(outdir / "intrinsic_geometry_issues.csv", all_intrinsic_issues)
|
||
write_csv(outdir / "geometry_stats_by_kind.csv", geometry_stats)
|
||
|
||
cat_counts = Counter(
|
||
x.get("issue_category", "UNKNOWN")
|
||
for x in all_intrinsic_issues
|
||
)
|
||
cat_rows = [
|
||
{"issue_category": k, "count": v}
|
||
for k, v in cat_counts.most_common()
|
||
]
|
||
write_csv(
|
||
outdir / "intrinsic_issue_category_counts.csv",
|
||
cat_rows,
|
||
["issue_category", "count"]
|
||
)
|
||
|
||
invalid_reason_counts = Counter(
|
||
x.get("invalid_reason")
|
||
for x in all_intrinsic_issues
|
||
if x.get("invalid_reason")
|
||
)
|
||
write_csv(
|
||
outdir / "invalid_geometry_reason_counts.csv",
|
||
[
|
||
{"reason": k, "count": v}
|
||
for k, v in invalid_reason_counts.most_common()
|
||
],
|
||
["reason", "count"]
|
||
)
|
||
|
||
obj_counts = Counter(
|
||
(x.get("object_class") or "", x.get("issue_category") or "")
|
||
for x in all_intrinsic_issues
|
||
)
|
||
write_csv(
|
||
outdir / "issue_object_class_counts.csv",
|
||
[
|
||
{"object_class": k[0], "issue_category": k[1], "count": v}
|
||
for k, v in obj_counts.most_common()
|
||
],
|
||
["object_class", "issue_category", "count"]
|
||
)
|
||
|
||
cell_counts = Counter(
|
||
(x.get("cell_name") or "", x.get("issue_category") or "")
|
||
for x in all_intrinsic_issues
|
||
)
|
||
write_csv(
|
||
outdir / "issue_cell_counts.csv",
|
||
[
|
||
{"cell_name": k[0], "issue_category": k[1], "count": v}
|
||
for k, v in cell_counts.most_common()
|
||
],
|
||
["cell_name", "issue_category", "count"]
|
||
)
|
||
|
||
# 每类最多 10 条样本
|
||
samples = []
|
||
sample_seen = Counter()
|
||
for x in all_intrinsic_issues:
|
||
cat = x.get("issue_category", "UNKNOWN")
|
||
if sample_seen[cat] < 10:
|
||
samples.append(x)
|
||
sample_seen[cat] += 1
|
||
write_csv(outdir / "issue_samples.csv", samples)
|
||
|
||
report["issue_summary"] = {
|
||
"mapping_issue_count": len(mapping_issues),
|
||
"intrinsic_issue_count": len(all_intrinsic_issues),
|
||
"sounding_feature_bbox_issue_count": len(sounding_bbox_issues),
|
||
"coverage_issue_count": len(coverage_issues),
|
||
"intrinsic_by_category": dict(cat_counts),
|
||
"issues_with_cell_boundary_context": sum(
|
||
1 for x in all_intrinsic_issues
|
||
if "CELL_EDGE" in (x.get("boundary_context") or "")
|
||
),
|
||
"issues_with_dataset_bbox_context": sum(
|
||
1 for x in all_intrinsic_issues
|
||
if "DATASET_BBOX_EDGE" in (x.get("boundary_context") or "")
|
||
),
|
||
}
|
||
|
||
# 基线匹配
|
||
acceptance = {}
|
||
for k, exp in EXPECTED_DATASET1.items():
|
||
actual = report["counts"].get(k)
|
||
acceptance[k] = {
|
||
"actual": actual,
|
||
"expected": exp,
|
||
"pass": actual == exp,
|
||
}
|
||
|
||
acceptance["raw_geotype16_expected_parent_count"] = {
|
||
"actual": report["counts"]["raw_geotype16_meta"],
|
||
"expected": 4406,
|
||
"pass": report["counts"]["raw_geotype16_meta"] == 4406,
|
||
"interpretation": "EXPECTED_PARENT_NO_GEOMETRY, not intrinsic geometry corruption",
|
||
}
|
||
|
||
acceptance["mapping"] = {
|
||
"issue_count": len(mapping_issues),
|
||
"pass": len(mapping_issues) == 0,
|
||
}
|
||
|
||
# 注意：intrinsic issue 非 0 时不自动称“Importer FAIL”。
|
||
# 需要结合 category、对象类型和 boundary_context 分析。
|
||
acceptance["intrinsic_geometry"] = {
|
||
"issue_count": len(all_intrinsic_issues),
|
||
"pass_if_zero": len(all_intrinsic_issues) == 0,
|
||
"note": (
|
||
"非零时需要人工分类；局部 Dataset 边界不连续本身不属于 FAIL。"
|
||
),
|
||
}
|
||
report["acceptance"] = acceptance
|
||
|
||
# JSON
|
||
(outdir / "spatial_audit_report.json").write_text(
|
||
json.dumps(report, ensure_ascii=False, indent=2, default=str),
|
||
encoding="utf-8"
|
||
)
|
||
|
||
# TXT summary
|
||
lines = []
|
||
lines.append("CM93 -> MySQL 实体空间 Geometry 只读审计 v3")
|
||
lines.append("=" * 78)
|
||
lines.append(f"Generated: {report['generated_at']}")
|
||
lines.append(f"Database: {args.database}")
|
||
lines.append(f"Dataset ID: {args.dataset_id}")
|
||
lines.append(f"DB version: {report['server']['version']}")
|
||
lines.append(f"WKB fn: {funcs['wkb_fn']}")
|
||
lines.append(f"SRID fn: {funcs['srid_fn']}")
|
||
lines.append(f"Dataset bbox: {dataset_bbox}")
|
||
lines.append("")
|
||
lines.append("DATASET SCOPE")
|
||
lines.append(
|
||
" 当前 CM93JAPAN 是日本沿海局部子集，不是完整 CM93。"
|
||
)
|
||
lines.append(
|
||
" 跨未导入海域/邻区的几何或对象连续性不作为本轮 FAIL 条件。"
|
||
)
|
||
lines.append(
|
||
" boundary_context 仅用于后续解释，不会自动把记录判成 importer bug。"
|
||
)
|
||
lines.append("")
|
||
lines.append("BASELINE COUNTS")
|
||
for k in ["cells", "features", "attributes", "sounding_points",
|
||
"raw_geotype16_meta", "raw_geotype16_geometry_issue"]:
|
||
lines.append(f" {k:36s} {report['counts'].get(k)}")
|
||
lines.append("")
|
||
lines.append("FEATURE ↔ GEOMETRY TABLE MAPPING")
|
||
lines.append(f" mapping issues: {len(mapping_issues)}")
|
||
for k, v in sorted(mapping_stats.get("mapping_issue_categories", {}).items()):
|
||
lines.append(f" {k:40s} {v}")
|
||
lines.append("")
|
||
lines.append("INTRINSIC GEOMETRY STATS")
|
||
for s in geometry_stats:
|
||
lines.append(
|
||
f" {s['kind']:<9} rows={s.get('rows',0):>9} "
|
||
f"parsed={s.get('parsed',0):>9} "
|
||
f"empty={s.get('empty',0):>6} "
|
||
f"valid={s.get('valid',0):>9} "
|
||
f"invalid={s.get('invalid',0):>6} "
|
||
f"issues={s.get('issue_count', s.get('issue_count_excluding_feature_bbox',0)):>6}"
|
||
)
|
||
lines.append(
|
||
f" actual types={s.get('actual_geom_types')}"
|
||
)
|
||
if s["kind"] == "LINE":
|
||
lines.append(
|
||
f" non_simple_line_info={s.get('non_simple_line_info',0)} "
|
||
f"(informational; not automatically invalid)"
|
||
)
|
||
lines.append("")
|
||
lines.append("INTRINSIC ISSUE CATEGORIES")
|
||
if cat_counts:
|
||
for k, v in cat_counts.most_common():
|
||
lines.append(f" {k:36s} {v}")
|
||
else:
|
||
lines.append(" none")
|
||
lines.append("")
|
||
lines.append("SOUNDING FEATURE BBOX")
|
||
lines.append(f" bbox issues: {len(sounding_bbox_issues)}")
|
||
lines.append("")
|
||
lines.append("COVERAGE")
|
||
lines.append(f" coverage issues: {len(coverage_issues)}")
|
||
lines.append("")
|
||
lines.append("BOUNDARY CONTEXT (INFORMATIONAL)")
|
||
lines.append(
|
||
f" intrinsic issues touching cell bbox edge: "
|
||
f"{report['issue_summary']['issues_with_cell_boundary_context']}"
|
||
)
|
||
lines.append(
|
||
f" intrinsic issues touching dataset bbox edge: "
|
||
f"{report['issue_summary']['issues_with_dataset_bbox_context']}"
|
||
)
|
||
lines.append("")
|
||
lines.append("IMPORTANT INTERPRETATION")
|
||
lines.append(
|
||
" raw_geotype=16 / META 的 4406 条正常 parent-no-geometry "
|
||
"不进入本轮 intrinsic geometry issue。"
|
||
)
|
||
lines.append(
|
||
" Line 在 Cell/Dataset 边缘结束、缺少未导入邻区对象，不属于本轮 geometry corruption。"
|
||
)
|
||
lines.append(
|
||
" Polygon 自相交、WKB 无法解析、坐标越界、类型串表等属于单 geometry 内在问题，"
|
||
"需要单独分析。"
|
||
)
|
||
lines.append("")
|
||
lines.append("PBF: OUT OF SCOPE / NOT STARTED")
|
||
|
||
(outdir / "spatial_audit_summary.txt").write_text(
|
||
"\n".join(lines),
|
||
encoding="utf-8"
|
||
)
|
||
|
||
conn.rollback()
|
||
|
||
zpath = zip_directory(outdir)
|
||
print("")
|
||
print("审计完成。")
|
||
print(f"输出目录: {outdir.resolve()}")
|
||
print(f"ZIP: {zpath.resolve()}")
|
||
print("")
|
||
print("请把 ZIP 整包发给 ChatGPT 分析。")
|
||
|
||
except KeyboardInterrupt:
|
||
if conn:
|
||
try:
|
||
conn.rollback()
|
||
except Exception:
|
||
pass
|
||
print("\n用户中断。")
|
||
sys.exit(130)
|
||
except Exception:
|
||
if conn:
|
||
try:
|
||
conn.rollback()
|
||
except Exception:
|
||
pass
|
||
tb = traceback.format_exc()
|
||
print(tb)
|
||
try:
|
||
(outdir / "fatal_error.txt").write_text(tb, encoding="utf-8")
|
||
except Exception:
|
||
pass
|
||
sys.exit(1)
|
||
finally:
|
||
if conn:
|
||
try:
|
||
conn.close()
|
||
except Exception:
|
||
pass
|
||
|
||
|
||
if __name__ == "__main__":
|
||
main()
|