Files
CMPlanning/scripts/build_solids.py
T
조선두andClaude Haiku 5.5 ff2c057aa7 feat: solid model of main line 4+745~5+480 (STEP + GLB), per-pattern selection, viewer rework (Task 8)
- scripts/build_solids.py: B-rep solids from standard section 03-001 (rev2, approved 2026-10-21)
  - lining per 지보패턴 run (independently selectable), excavation band per run
  - road surface (계획고, 280 mm) with per-station superelevation; base/invert ASSUMED
  - side drains (box minus pipe), H/I boxes = 공동구 (assumed); 시설한계 50 mm as reference wires only
  - staged run (4 GB RAM / 180 s per call); STEP and GLB from the same B-rep
- model/tunnel_solids.glb, model/tunnel_solids_meta.json, data/verify_solids.json
- viewer/index.html: solids, pattern and component toggles, click readout, no development text
- STEP (84.6 MB) and model/parts/ are git-ignored; the STEP is generated locally by the script

Co-Authored-By: Claude Haiku 5.5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01A5Cun3RHQKna9qxE6ybZ4A
2026-10-21 15:32:54 +00:00

328 lines
19 KiB
Python

#!/usr/bin/env python3
"""Task 8: B-rep solids of the main line (본선) 4+745..5+480 -> STEP (Rhino/Autodesk) + GLB (web viewer).
Section geometry: standard section C0060303-001 (drawing 03-001), revision rev2 approved 2026-10-21.
- lining (라이닝): arch R6982 (inner) / R7282 (outer) about centre C; haunch arcs R3712/R3910 (inner),
R4454/R4699 (outer) run below road level as drawn; inner walls and inner chord (road level) as drawn.
Lining is symmetric (no superelevation). H/I boxes removed from lining (2D).
- excavation band (굴착 여유대): R7402 / R4574 / R4819 minus lining outer, above road level. Symmetric.
- road surface / 계획고 (노면, 포장 280 mm): cross slope S per station (superelevation). Only this part tilts.
- base / invert (기층·인버트, ASSUMED: not in drawing): region between lining outer and inner chord minus
pavement and side drains. Its top follows the tilted pavement bottom.
- side drains (측구): box 300 x 836 mm minus pipe r150 (drawing circles).
- H/I boxes (공동구, ASSUMED by user decision 2026-10-21).
- 시설한계 50 mm: REFERENCE polylines only (not solids), as decided by the user.
Coordinates: local section x = lateral along the LEFT normal of the alignment (sign convention flagged),
z = height above road centre. Superelevation: z += (S/100)*x (sign flagged).
Outputs: model/tunnel_solids.step, model/tunnel_solids.glb, model/tunnel_solids_meta.json, data/verify_solids.json
Usage: python3 scripts/build_solids.py "$PWD"
"""
import csv, json, math, os, sys, tempfile
import numpy as np
from shapely.geometry import Polygon, box
from build123d import Vector, Wire, Solid, Compound, Color, export_step, export_stl, import_step
import trimesh
ROOT = sys.argv[1]
MODE = sys.argv[2] if len(sys.argv) > 2 else 'all' # parts_a | parts_b | assemble (4 GB RAM: run in stages)
PARTS = f"{ROOT}/model/parts"
RANGE = (int(sys.argv[3]), int(sys.argv[4])) if len(sys.argv) > 4 else (0, 10**6) # run index range for parts_a
os.makedirs(PARTS, exist_ok=True)
C_CSV = f"{ROOT}/data/alignment_plan_001_centre.csv"
P_CSV = f"{ROOT}/data/profile_gl_002_draft.csv"
PAT_CSV = f"{ROOT}/data/pattern_stations_dangjin.csv"
S0, S1 = 4745, 5480
STEP_M = 10.0 # section spacing (m): ruled loft between sections (STEP size vs curvature)
CHUNK_M = 100
PIPE_SEG = 24
CLEAR_OFFSET = 0.05
CLEAR_LINES = 6
CLEAR_RIB_W = 0.02
# ---------- drawing geometry (mm) -> local m ----------
CX, CY, Y0 = 14126.29, 15305.02, 15672.0
def P(x, y): return ((x - CX) / 1000.0, (y - Y0) / 1000.0)
def arc(cx, cy, r, a0, a1, n=120):
return [(cx + r * math.cos(math.radians(a0 + (a1 - a0) * i / n)),
cy + r * math.sin(math.radians(a0 + (a1 - a0) * i / n))) for i in range(n + 1)]
def loc(pts): return [P(x, y) for x, y in pts]
WALL_L_INNER = [(8771, 15751), (8771, 15826), (8651, 16001), (8621, 16411), (8471, 16411),
(8171, 16411), (8045, 16411), (8034, 16446)]
WALL_R_INNER = [(19481, 15567), (19481, 15642), (19601, 15817), (19631, 16227), (19781, 16227),
(20181, 16227), (20251, 16239)]
INNER_MM = (WALL_L_INNER + arc(11621.3, 17406.9, 3712, 195.6, 140) + arc(CX, CY, 6982, 140, 40)
+ arc(16479.6, 17279.7, 3910, 40, -15.6) + WALL_R_INNER[::-1])
OUTER_MM = (arc(11959.9, 17122.8, 4454, 204.7, 140) + arc(CX, CY, 7282, 140, 40)
+ arc(16105.0, 16965.3, 4699, 40, -23.7))
EX_MM = (arc(11959.9, 17122.8, 4574, 204.0, 140) + arc(CX, CY, 7402, 140, 40)
+ arc(16105.0, 16965.3, 4819, 40, -23.0))
inner_poly = Polygon(loc(INNER_MM)).buffer(0)
outer_poly = Polygon(loc(OUTER_MM)).buffer(0)
ex_poly = Polygon(loc(EX_MM)).buffer(0)
def _fill_slivers(g, tol=1e-3):
if g.geom_type != 'Polygon': return g
return Polygon(g.exterior, [h for h in g.interiors if Polygon(h).area > tol])
ring = _fill_slivers(outer_poly.difference(inner_poly)) # 37 mm wall/haunch join leaves 1e-5 m2 slivers
CHORD_A, CHORD_B = P(8771, 15751), P(19481, 15567)
_k = (CHORD_B[1] - CHORD_A[1]) / (CHORD_B[0] - CHORD_A[0])
def chord_z(x): return CHORD_A[1] + _k * (x - CHORD_A[0])
below_chord = Polygon([(-9, chord_z(-9)), (9, chord_z(9)), (9, -5), (-9, -5)])
above_chord = Polygon([(-9, chord_z(-9)), (9, chord_z(9)), (9, 9), (-9, 9)])
XL = P(9526, 15764)[0]; XR = P(18726, 15580)[0]
T_FIX = -0.02
PAV_T = 0.28
def pav_corners(t):
return [(XL, t * XL), (XR, t * XR), (XR, t * XR - PAV_T), (XL, t * XL - PAV_T)]
pav_fix = Polygon(pav_corners(T_FIX))
def rect_local(x0_mm, x1_mm, y0_mm, y1_mm):
(a, b), (c, d) = P(x0_mm, y0_mm), P(x1_mm, y1_mm)
return box(min(a, c), min(b, d), max(a, c), max(b, d))
DRAIN_L = rect_local(9101, 9401, 14864, 15700)
DRAIN_R = rect_local(18851, 19151, 14680, 15520)
PIPE_L = (P(9251.3, 15064)[0], P(9251.3, 15064)[1], 0.150)
PIPE_R = (P(19001.3, 14880)[0], P(19001.3, 14880)[1], 0.150)
HBOX = rect_local(8171, 8471, 16051, 16411)
IBOX = rect_local(19781, 20181, 15867, 16227)
lin_full = _fill_slivers(ring.intersection(above_chord))
assert lin_full.geom_type == 'Polygon' and len(lin_full.interiors) == 0, f"lining topology: {lin_full.geom_type}"
LINING_POLY = lin_full.simplify(0.01, preserve_topology=True) # 10 mm: removes 30 um edges (loft invalid) and cuts ~700 -> ~60 section vertices (STEP size)
EXCAV_POLY = ex_poly.difference(outer_poly).intersection(box(-9, 0.0, 9, 9)).simplify(0.01, preserve_topology=True)
assert EXCAV_POLY.geom_type == 'Polygon' and len(EXCAV_POLY.interiors) == 0, f"excavation topology: {EXCAV_POLY.geom_type}"
base_fix = ring.intersection(below_chord).difference(pav_fix).difference(DRAIN_L).difference(DRAIN_R)
assert base_fix.geom_type == 'Polygon' and len(base_fix.interiors) == 0, f"base topology: {base_fix.geom_type}"
BASE_PTS_FIX = np.array(base_fix.exterior.coords)[:-1]
def _find(pt):
d = np.hypot(BASE_PTS_FIX[:, 0] - pt[0], BASE_PTS_FIX[:, 1] - pt[1])
i = int(np.argmin(d)); assert d[i] < 1e-6, f"pavement corner not on base boundary: {d[i]}"
return i
IDX_BL = _find((XL, T_FIX * XL - PAV_T)); IDX_BR = _find((XR, T_FIX * XR - PAV_T))
def base_pts(t):
pts = BASE_PTS_FIX.copy()
pts[IDX_BL, 1] = t * XL - PAV_T
pts[IDX_BR, 1] = t * XR - PAV_T
return pts
LIN_PTS = np.array(LINING_POLY.exterior.coords)[:-1]
EXC_PTS = np.array(EXCAV_POLY.exterior.coords)[:-1]
def PIPE_PTS(cx, cz, r):
return np.array([(cx + r * math.cos(2 * math.pi * k / PIPE_SEG), cz + r * math.sin(2 * math.pi * k / PIPE_SEG)) for k in range(PIPE_SEG)])
DRAIN_L_PTS = np.array(DRAIN_L.exterior.coords)[:-1]
DRAIN_R_PTS = np.array(DRAIN_R.exterior.coords)[:-1]
HBOX_PTS = np.array(HBOX.exterior.coords)[:-1]
IBOX_PTS = np.array(IBOX.exterior.coords)[:-1]
# ---------- alignment, profile ----------
cn = np.array([[int(r["station_m"]), float(r["centre_x"]), float(r["centre_y"])]
for r in csv.DictReader(open(C_CSV, encoding="utf-8"))])
dg_d = {}
for r in csv.DictReader(open(P_CSV, encoding="utf-8")):
if float(r["across_m"]) >= -25:
dg_d[int(r["station_anchor"])] = float(r["value"])
dg_s = np.array(sorted(dg_d)); dg_v = np.array([dg_d[k] for k in dg_s])
cn_s = cn[:, 0]
def axis(s):
x = np.interp(s, cn_s, cn[:, 1]); y = np.interp(s, cn_s, cn[:, 2])
i = min(max(np.searchsorted(cn_s, s) - 1, 0), len(cn_s) - 2)
dx, dy = cn[i + 1, 1] - cn[i, 1], cn[i + 1, 2] - cn[i, 2]
L = math.hypot(dx, dy)
return x, y, dx / L, dy / L
ORX, ORY = axis(S0)[0], axis(S0)[1]
def slope_pct(s):
if s < 5105: return -2.0
if s < 5165: return -2.0 + 4.0 * (s - 5105) / 60.0
return 2.0
def world(s, x, z):
ax, ay, tx, ty = axis(s)
nx, ny = -ty, tx
zr = float(np.interp(s, dg_s, dg_v))
return (ax - ORX + x * nx, ay - ORY + x * ny, zr + z)
# ---------- B-rep helpers ----------
def loft(stations, pts2_fn):
wires = []
for s in stations:
pts3 = [world(s, float(x), float(z)) for x, z in pts2_fn(s)]
wires.append(Wire.make_polygon([Vector(*p) for p in pts3], close=True))
sol = Solid.make_loft(wires, ruled=True)
if not sol.is_valid: # end caps of lofts with near-degenerate edges; ShapeFix repairs them
sol = sol.fix()
return sol
def stations_of(a, b, step=STEP_M):
n = max(1, int(round((b - a) / step)))
return [float(v) for v in np.linspace(a, b, n + 1)]
def chunks():
out = []; a = S0
while a < S1:
b = min(a + CHUNK_M, S1); out.append((a, b)); a = b
return out
PAL = {"A": (42, 157, 143), "B-1": (231, 111, 81), "B-2": (244, 162, 97), "C-1": (233, 196, 106),
"C-2": (138, 177, 125), "D-1": (38, 70, 83), "D-2": (100, 130, 160), "P-2": (190, 110, 140),
"P-3": (150, 90, 60), "RP-1": (80, 170, 200), "WL-1": (200, 200, 120), "WL-2": (120, 90, 170),
"개착부": (170, 170, 170)}
COL = {"base": (141, 110, 99), "pav": (96, 125, 139), "drain": (30, 136, 229), "hi": (253, 216, 53),
"clear": (255, 111, 0)}
runs = []
for r in csv.DictReader(open(PAT_CSV, encoding="utf-8")):
a, b = int(float(r["sta_from"])), int(float(r["sta_to"])); pat = r["lining_pattern"]
a2, b2 = max(a, S0), min(b, S1)
if b2 > a2: runs.append((a2, b2, pat))
runs.sort()
assert runs[0][0] == S0 and runs[-1][1] == S1, f"pattern coverage {runs[0][0]}..{runs[-1][1]}"
for (a1, b1, _), (a2, _, _) in zip(runs, runs[1:]):
assert b1 == a2, f"pattern gap/overlap at {b1} vs {a2}"
solids, meta_c, clear_wires = [], [], []
def add(sol, cid, label, group, color, pattern=None, a=None, b=None):
if not sol.is_valid:
sol = sol.fix() # booleans on lofts can leave cap faces invalid
assert sol.is_valid, f"invalid solid {cid}"
sol.label = cid; sol.color = Color(*[c / 255 for c in color])
solids.append(sol)
meta_c.append({"id": cid, "label": label, "group": group, "pattern": pattern,
"from_sta": a, "to_sta": b, "volume_m3": round(sol.volume, 4), "color": list(color)})
for i, (a, b, pat) in enumerate(runs):
if MODE != 'parts_a' or not (RANGE[0] <= i < RANGE[1]): continue
print('lining run', i + 1, pat, a, b, flush=True)
lin_run = loft(stations_of(a, b), lambda s: LIN_PTS)
lin_run = lin_run - loft(stations_of(a, b), lambda s: HBOX_PTS) - loft(stations_of(a, b), lambda s: IBOX_PTS)
add(lin_run, f"LIN_{i+1:02d}_{pat}", f"라이닝 {pat} {a/1000:.3f}~{b/1000:.3f}",
"라이닝", PAL[pat], pat, a, b)
for i, (a, b, pat) in enumerate(runs):
if MODE != 'parts_a' or not (RANGE[0] <= i < RANGE[1]): continue
add(loft(stations_of(a, b), lambda s: EXC_PTS), f"EXC_{i+1:02d}_{pat}", f"굴착 여유대 {pat} {a/1000:.3f}~{b/1000:.3f}",
"굴착면", (150, 120, 90), pat, a, b)
for (a, b) in chunks():
if MODE != 'parts_b': break
sts = stations_of(a, b)
add(loft(sts, lambda s: np.array(pav_corners(slope_pct(s) / 100.0))), f"PAV_{a}", f"노면 계획고 포장 {a}~{b}",
"노면(계획고)", COL["pav"], None, a, b)
add(loft(sts, lambda s: base_pts(slope_pct(s) / 100.0)), f"BASE_{a}", f"기층·인버트(가정) {a}~{b}",
"기층·인버트(가정)", COL["base"], None, a, b)
for side, dpts, pipe in (("L", DRAIN_L_PTS, PIPE_L), ("R", DRAIN_R_PTS, PIPE_R)):
bx = loft(sts, lambda s, d=dpts: d)
pp = loft(sts, lambda s, p=pipe: PIPE_PTS(p[0], p[1], p[2]))
add(bx - pp, f"DRN{side}_{a}", f"측구 {side} {a}~{b}", "배수(측구)", COL["drain"], None, a, b)
add(loft(sts, lambda s: HBOX_PTS), f"HBOX_{a}", f"공동구 H(가정) {a}~{b}", "공동구(가정)", COL["hi"], None, a, b)
add(loft(sts, lambda s: IBOX_PTS), f"IBOX_{a}", f"공동구 I(가정) {a}~{b}", "공동구(가정)", COL["hi"], None, a, b)
# 시설한계 50 mm reference lines (wires, not solids)
clr = inner_poly.buffer(CLEAR_OFFSET).exterior
clr_pts = [(x, z) for x, z in list(clr.coords)[:-1] if z >= chord_z(x) - 1e-9 and abs(x) <= 5.4]
idxs = np.round(np.linspace(0, len(clr_pts) - 1, CLEAR_LINES)).astype(int)
clear_lines = []
sts_all = stations_of(S0, S1)
for k, ii in enumerate(idxs):
x0, z0 = clr_pts[ii]
pts3 = [world(s, x0, z0) for s in sts_all]
clear_lines.append(pts3)
if MODE != 'parts_b':
continue
w = Wire.make_polygon([Vector(*p) for p in pts3], close=False)
w.label = f"REF_CLEAR_{k+1}"; w.color = Color(1.0, 0.44, 0.0)
clear_wires.append(w)
meta_c.append({"id": f"REF_CLEAR_{k+1}", "label": f"시설한계 50 mm 참조선 {k+1}", "group": "시설한계(참조선)",
"pattern": None, "from_sta": S0, "to_sta": S1, "volume_m3": 0, "color": list(COL["clear"]), "reference": True})
if MODE in ('parts_a', 'parts_b'):
tag = f'lin_exc_{RANGE[0]:02d}' if MODE == 'parts_a' else 'other'
export_step(Compound(children=solids + clear_wires, label=tag), f"{PARTS}/{tag}.step")
json.dump(meta_c, open(f"{PARTS}/{tag}.json", 'w', encoding='utf-8'), ensure_ascii=False, indent=1)
print('parts written', tag, len(solids), 'solids', len(clear_wires), 'ref wires')
sys.exit(0)
# ---- post-processing modes (each call < 180 s, 4 GB RAM) ----
import glob
TAGS = sorted(os.path.basename(x)[:-5] for x in glob.glob(f"{PARTS}/*.json"))
from build123d import Solid as _Solid
if MODE == 'mesh': # python build_solids.py ROOT mesh <tag>
tag = sys.argv[3]
imp = import_step(f"{PARTS}/{tag}.step")
items = [x for x in (imp.children or list(imp.solids())) if isinstance(x, _Solid)]
pj = [m for m in json.load(open(f"{PARTS}/{tag}.json", encoding='utf-8')) if not m.get('reference')]
assert len(items) == len(pj), f"{tag}: {len(items)} solids vs {len(pj)} meta"
sc = trimesh.Scene(); tmp = tempfile.mkdtemp(prefix='stl_')
for k, (sol, mm) in enumerate(zip(items, pj)):
stl = os.path.join(tmp, f"s{k}.stl")
export_stl(sol, stl, tolerance=0.002, angular_tolerance=0.2)
m = trimesh.load(stl, force='mesh')
m.visual.face_colors = np.tile(np.array(list(mm['color']) + [255], dtype=np.uint8), (len(m.faces), 1))
sc.add_geometry(m, geom_name=mm['id'])
sc.export(f"{PARTS}/{tag}.glb")
print('mesh', tag, len(sc.geometry), 'geometries'); sys.exit(0)
if MODE == 'step': # python build_solids.py ROOT step
allitems = []
for tag in TAGS:
imp = import_step(f"{PARTS}/{tag}.step")
allitems += list(imp.children or list(imp.solids()))
os.makedirs(f"{ROOT}/model", exist_ok=True)
export_step(Compound(children=allitems, label="본선 4+745~5+480 표준단면 모델"), f"{ROOT}/model/tunnel_solids.step")
print('step merged', len(allitems), 'items', os.path.getsize(f"{ROOT}/model/tunnel_solids.step"), 'bytes'); sys.exit(0)
if MODE == 'final': # python build_solids.py ROOT final
meta_c = []
for tag in TAGS:
meta_c += json.load(open(f"{PARTS}/{tag}.json", encoding='utf-8'))
scene = trimesh.Scene()
for tag in TAGS:
g = trimesh.load(f"{PARTS}/{tag}.glb")
for name, geom in g.geometry.items():
scene.add_geometry(geom, geom_name=name)
glb_faces = {}
for k, pts3 in enumerate(clear_lines): # 시설한계 reference ribbons (visual only; exact wires in STEP)
P3 = np.array(pts3); V = []; F = []
for i in range(len(P3)):
d = P3[min(i + 1, len(P3) - 1)] - P3[max(i - 1, 0)]
d = d / (np.linalg.norm(d) + 1e-12)
w = np.cross(d, np.array([0, 0, 1.0])); w = w / (np.linalg.norm(w) + 1e-12) * CLEAR_RIB_W
V += [P3[i] - w, P3[i] + w]
for i in range(len(P3) - 1):
a0, b0, c0, d0 = 2 * i, 2 * i + 1, 2 * i + 2, 2 * i + 3
F += [[a0, c0, b0], [b0, c0, d0]]
m = trimesh.Trimesh(vertices=np.array(V), faces=np.array(F), process=False)
m.visual.face_colors = np.tile(np.array([255, 111, 0, 255], dtype=np.uint8), (len(m.faces), 1))
scene.add_geometry(m, geom_name=f"REF_CLEAR_{k+1}")
glb_path = f"{ROOT}/model/tunnel_solids.glb"
scene.export(glb_path)
step_path = f"{ROOT}/model/tunnel_solids.step"
solid_meta = [x for x in meta_c if not x.get('reference')]
meta = {
"model": ["tunnel_solids.step (STEP B-rep, Rhino/Autodesk)", "tunnel_solids.glb (web viewer, same B-rep tessellated)"],
"scope": f"{S0}~{S1} (4+745~5+480)", "units": "m",
"origin_plan_xy": [round(ORX, 3), round(ORY, 3)],
"sections": {"lining": "rev2: 라이닝 R6982/R7282, 하단 호 도면 끝점까지", "excavation": "굴착 여유대 R7402/R4574/R4819 (노면 위)",
"road": "계획고 노면 포장 280 mm, 측점별 편경사 S (노면만)", "base": "기층·인버트: 가정 (도면 미표현)",
"drains": "측구 300x836 mm, 관 r150", "hi": "H/I 박스: 공동구 가정(사용자 결정)", "clearance": "시설한계 50 mm 참조선 (솔리드 아님)"},
"assumptions": ["기층·인버트 형상·두께는 가정 (도면 미표현, 사용자 결정 '가정하여 생성')",
"H/I 박스 = 공동구 가정 (사용자 결정)",
"편경사 부호 z = (S/100)·x, x = 접선 좌측 법선 (부호 확인 필요)",
"라이닝·굴착 여유대는 대칭 (편경사 미적용)",
"측점 간 축선은 선형 보간, 단면은 10 m 간격 Ruled loft (곡선부 근사)",
"라이닝 H/I 박스 제외 부분은 3D 불리언 차감 (공동구 가정)"],
"pattern_runs": [{"id": f"LIN_{i+1:02d}_{p}", "pattern": p, "from": a, "to": b} for i, (a, b, p) in enumerate(runs)],
"components": meta_c,
"open": ["E-1, E-2, RP-2, WL-3 패턴 배치 [미확정] (당진 표 밖)", "NP-5BE-1-WL-4 해석 [확인 필요]",
"04 공동구·배수 도면 상세 미반영 (측구·공동구 치수는 03-001 단면 기준)", "편경사 부호(좌·우) 확인 필요"]
}
json.dump(meta, open(f"{ROOT}/model/tunnel_solids_meta.json", "w", encoding="utf-8"), ensure_ascii=False, indent=2)
ver = {"solid_count": len(solid_meta), "runs": len(runs), "pattern_coverage": f"{runs[0][0]}..{runs[-1][1]} contiguous",
"lining_runs": sum(1 for x in solid_meta if x['group'] == '라이닝'),
"excavation_runs": sum(1 for x in solid_meta if x['group'] == '굴착면'),
"volume_sum_m3": round(sum(x['volume_m3'] for x in solid_meta), 3),
"step_bytes": os.path.getsize(step_path), "glb_bytes": os.path.getsize(glb_path),
"glb_geometries": len(trimesh.load(glb_path).geometry)}
def _overlap(sample):
base = Polygon(base_pts(slope_pct(sample) / 100.0)); pav = Polygon(pav_corners(slope_pct(sample) / 100.0))
return {"lining_vs_excavation": round(LINING_POLY.intersection(EXCAV_POLY).area, 9),
"base_vs_pavement": round(base.intersection(pav).area, 9),
"lining_vs_base": round(LINING_POLY.intersection(base).area, 9),
"pavement_vs_lining": round(pav.intersection(LINING_POLY).area, 9),
"base_vs_drainL": round(base.intersection(DRAIN_L).area, 9),
"base_vs_drainR": round(base.intersection(DRAIN_R).area, 9)}
ver["section_overlap_m2_by_station"] = {str(st): _overlap(st) for st in [4745, 5105, 5165, 5480]}
json.dump(ver, open(f"{ROOT}/data/verify_solids.json", "w", encoding="utf-8"), ensure_ascii=False, indent=2)
print(json.dumps(ver, ensure_ascii=False, indent=1)); sys.exit(0)