diff --git a/resources/tester/test_612_precompute.py b/resources/tester/test_612_precompute.py new file mode 100644 index 000000000..32b06e9d0 --- /dev/null +++ b/resources/tester/test_612_precompute.py @@ -0,0 +1,124 @@ +"""PLAN 612-1 — 전국 미리 계산 일괄 도구(`tools/precompute_terrain.py`) 지킴 시험. + +도엽 목록 칸 모음 · 고름 · 산지 차례 · 좌표계 띠 · 이어 하기(성공 건너뜀 · 실패 다시) · +죽은 잠금 치움 — 셈 함수는 가짜로 바꿔 서버 · 전국 파일 없이 돎. +""" + +from __future__ import annotations + +import json +import os +import sys +import time +from pathlib import Path + +import numpy as np +import pytest +import shapely +from pyproj import Transformer + +ROOT = Path(__file__).resolve().parents[2] +sys.path.insert(0, str(ROOT)) + +from B04_PreProcess.B04_PreProcess_Engine_MapSheet import latlon_to_sheet5k # noqa: E402 +from tools import precompute_terrain as tool # noqa: E402 + + +def _line_5179(lon_lat: list[tuple[float, float]]) -> bytes: + to = Transformer.from_crs("EPSG:4326", "EPSG:5179", always_xy=True) + xs, ys = to.transform(*zip(*lon_lat)) + return shapely.to_wkb(shapely.LineString(list(zip(xs, ys)))) + + +def test_aggregate_cells_length_and_relief() -> None: + # 도엽 하나 안 두 줄(높이 100 · 300) + 다른 도엽 한 줄 + a = [(129.0010 + 0.0001 * i, 36.9901) for i in range(200)] # 약 1.8 km + b = [(129.0010 + 0.0001 * i, 36.9951) for i in range(100)] + c = [(129.1001 + 0.0001 * i, 36.9901) for i in range(50)] + cells = tool.aggregate([_line_5179(a), _line_5179(b), _line_5179(c)], np.array([100, 300, 50])) + assert len(cells) == 2 + big = max(cells.values(), key=lambda v: v[0]) + assert big[1:] == [100.0, 300.0] + assert 2_300 < big[0] < 2_800 # 1.78 + 0.89 km — 솎음 현이라 조금 짧음 + + +def test_pick_sheets_prefix_box_relief_order() -> None: + index = { + "36901001": [10.0, 100.0, 400.0], + "36901002": [30.0, 100.0, 200.0], + "36901003": [30.0, 100.0, 500.0], + "36712089": [50.0, 0.0, 20.0], + } + assert tool.pick_sheets(index) == ["36901003", "36901001"] # 높이 차 150 m 이상 · 길이 순 + # 앞자리 고름은 높이 차 거름 안 씀 · 길이 같으면 높이 차 큰 것 먼저 + assert tool.pick_sheets(index, sheets=["36901"]) == ["36901003", "36901002", "36901001"] + assert tool.pick_sheets(index, sheets=["36901"], limit=1) == ["36901003"] + assert tool.pick_sheets(index, sheets=["99999999"]) == [] # 전국 파일 밖 = 안 고름 + assert tool.pick_sheets(index, min_relief=0, box=(127.5, 36.0, 128.0, 36.5)) == ["36712089"] + + +def test_epsg_band_follows_longitude() -> None: + assert tool.epsg_of(latlon_to_sheet5k(37.5, 125.5)) == 5185 + assert tool.epsg_of(latlon_to_sheet5k(37.5, 127.5)) == 5186 + assert tool.epsg_of("36901074") == 5187 + + +@pytest.fixture() +def fake_store(tmp_path, monkeypatch): + monkeypatch.setattr(tool.store, "TERRAIN_ROOT", tmp_path) + monkeypatch.setattr(tool.store, "_gpkg_stamp", lambda: "stamp-1") + calls: list[tuple[str, str]] = [] + fail = {"36901002"} + + def fake_run(step: str, no: str): + calls.append((step, no)) + if no in fail: + return step, no, {"ok": False, "error": "RuntimeError: 받기 실패"} + return step, no, {"ok": True, "s": 0.01, "fresh": True, "empty": False} + + monkeypatch.setattr(tool, "run_step", fake_run) + return tmp_path, calls, fail + + +def test_resume_skips_done_and_retries_failed(fake_store) -> None: + root, calls, fail = fake_store + sheets = ["36901001", "36901002", "36901003"] + progress = tool.load_progress() + tool.run(sheets, ["e07", "dem"], 1, progress, say=lambda _: None) + assert [s for s, _ in calls] == ["dem"] * 3 + ["e07"] * 3 # 차례 고정 + saved = json.loads((root / tool.PROGRESS_NAME).read_text(encoding="utf-8")) + assert saved["steps"]["dem"]["36901001"]["ok"] is True + assert saved["steps"]["dem"]["36901002"]["ok"] is False + + calls.clear() + fail.clear() + tool.run(sheets, ["dem", "e07"], 1, tool.load_progress(), say=lambda _: None) + assert calls == [("dem", "36901002"), ("e07", "36901002")] # 성공은 건너뜀 · 실패만 다시 + assert tool.summary(sheets, tool.load_progress())["dem"]["ok"] == 3 + + +def test_progress_restarts_when_source_changes(fake_store, monkeypatch) -> None: + root, calls, _ = fake_store + tool.run(["36901001"], ["dem"], 1, tool.load_progress(), say=lambda _: None) + monkeypatch.setattr(tool.store, "_gpkg_stamp", lambda: "stamp-2") + assert tool.load_progress()["steps"]["dem"] == {} + + +def test_clean_dead_locks_and_old_tmp(tmp_path) -> None: + folder = tmp_path / "36901001" + folder.mkdir() + dead = folder / ".dem_5187_2m.lock" + dead.write_text("999999999 0", encoding="utf-8") + live = folder / ".zero_acc_5187_2m.lock" + live.write_text(f"{os.getpid()} 0", encoding="utf-8") + old_tmp = folder / ".dem_5187_2m.npz.abc.tmp" + old_tmp.write_bytes(b"x") + stale = time.time() - tool.store.LOCK_STALE_S - 10 + os.utime(old_tmp, (stale, stale)) + new_tmp = folder / ".e07_5187_2m.npz.def.tmp" + new_tmp.write_bytes(b"x") + kept = folder / "dem_5187_2m.npz" + kept.write_bytes(b"x") + assert tool.clean_dead(tmp_path) == 2 + assert not dead.exists() and not old_tmp.exists() + assert live.exists() and new_tmp.exists() and kept.exists() diff --git a/tools/precompute_terrain.py b/tools/precompute_terrain.py new file mode 100644 index 000000000..5dbfc3e6b --- /dev/null +++ b/tools/precompute_terrain.py @@ -0,0 +1,407 @@ +# -*- coding: utf-8 -*- +"""전국 미리 계산 일괄 도구(PLAN 612-1 · 요구 493) — 공용 지형 창고에 도엽마다 미리 셈. + +실행: ./venv/Scripts/python.exe tools/precompute_terrain.py --sheets 36901 --workers 4 + --steps dem,zero,e07 단계(차례 고정 = 지형 → 0차곡 누적 → E07 · 뒤에 tiles) + --sheets 36901,36712089 도엽 번호 또는 앞자리(1:50,000 묶음) · 주면 높이 차 거름 안 씀 + --box lon0,lat0,lon1,lat1 도엽 가운데가 이 경위도 상자 안 + --min-relief 150 도엽 안 높이 차(m) 이상만(산지) · --limit N 앞 N 장만 + --workers 4 일꾼 수 · --status 진행 기록 요약만 · --rebuild-index 도엽 목록 다시 + +- 도엽 목록 = 전국 등고선 파일 꼭짓점이 닿는 도엽(0.025°) — `_precompute_sheets.json`(전국 파일 + 지문이 같으면 다시 안 셈) · 차례 = 등고선 길이 큰 것(산지)부터 · 같으면 높이 차. +- 좌표계 = 도엽 가운데 경도의 TM 띠(`tm_epsg_for_lon` · 5185~5187) · 격자 2 m. +- 셈 = 공용 함수만 부름(창고 지형 · D 0차곡 누적 · E07 도엽) — 저장 열쇠가 같으면 함수가 건너뜀 · + 쓰기 = 함수 쪽 원자 쓰기(끊겨도 반쯤 쓴 파일 없음). +- 이어 하기 = `_precompute_progress.json`(단계 · 도엽마다 성공 · 시간) — 성공한 도엽은 + 다시 안 부름 · 실패는 다음에 다시 · 전국 파일이 바뀌면 기록을 새로 시작. +- 시작 때 죽은 셈의 잠금(주인 프로세스 없음) · 오래된 임시 파일을 치움 · + 낮은 우선순위(BELOW_NORMAL). +""" + +from __future__ import annotations + +import argparse +import json +import math +import os +import sqlite3 +import sys +import time +from concurrent.futures import ProcessPoolExecutor, as_completed +from pathlib import Path +from typing import Any + +ROOT = Path(__file__).resolve().parent.parent +if str(ROOT) not in sys.path: + sys.path.insert(0, str(ROOT)) + +import numpy as np # noqa: E402 + +from B04_PreProcess.B04_PreProcess_Engine_MapSheet import ( # noqa: E402 + SHEET5K_SIZE_DEG, + latlon_to_sheet5k, + neighbors_3x3, + sheet5k_center, +) +from common_util import common_util_terrain_store as store # noqa: E402 +from common_util.common_util_crs import tm_epsg_for_lon # noqa: E402 +from common_util.common_util_json import atomic_write_json # noqa: E402 + +STEPS = ("dem", "zero", "e07") +CELL_M = 2.0 +PROGRESS_NAME = "_precompute_progress.json" +INDEX_NAME = "_precompute_sheets.json" +INDEX_CHUNK = 50_000 # 가닥 묶음(약 0.5 GB 읽기) +INDEX_STRIDE = 8 # 꼭짓점 솎음(도엽 2 km 에 비해 넉넉) +SAVE_EVERY_S = 10.0 +# 단계별 창고 파일 앞머리(크기 셈) +STEP_FILES = {"dem": ("dem_",), "zero": ("zero_acc_",), "e07": ("e07_", "e07si_")} +_ENVELOPE = {0: 0, 1: 32, 2: 48, 3: 48, 4: 64} + + +def below_normal() -> None: + """이 프로세스를 낮은 우선순위로(Windows) — 다른 일 · 서버를 막지 않게.""" + import psutil + + if os.name == "nt": + psutil.Process().nice(psutil.BELOW_NORMAL_PRIORITY_CLASS) + + +def epsg_of(sheet_no: str) -> int: + return tm_epsg_for_lon(sheet5k_center(sheet_no)[1]) + + +# ── 도엽 목록 ──────────────────────────────────────────────────────────────── + + +def aggregate(wkbs: list[bytes], heights: np.ndarray) -> dict[int, list[float]]: + """EPSG:5179 선(WKB) · 높이 → {칸 열쇠: [길이 m, 최저, 최고]} · 칸 = 0.025°(위도 · 경도).""" + import shapely + from pyproj import Transformer + + geoms = shapely.from_wkb(wkbs) + coords, owner = shapely.get_coordinates(geoms, return_index=True) + coords, owner = coords[::INDEX_STRIDE], owner[::INDEX_STRIDE] + if not len(coords): + return {} + lon, lat = Transformer.from_crs("EPSG:5179", "EPSG:4326", always_xy=True).transform( + coords[:, 0], coords[:, 1] + ) + keys = np.floor(lat / SHEET5K_SIZE_DEG).astype(np.int64) * 100_000 + np.floor( + lon / SHEET5K_SIZE_DEG + ).astype(np.int64) + seg = np.zeros(len(coords)) + same = owner[1:] == owner[:-1] + seg[:-1] = np.where(same, np.hypot(*(coords[1:] - coords[:-1]).T), 0.0) + z = np.asarray(heights, dtype=np.float64)[owner] + uniq, inv = np.unique(keys, return_inverse=True) + length = np.bincount(inv, weights=seg, minlength=len(uniq)) + lo = np.full(len(uniq), np.inf) + hi = np.full(len(uniq), -np.inf) + np.minimum.at(lo, inv, z) + np.maximum.at(hi, inv, z) + return {int(k): [float(a), float(b), float(c)] for k, a, b, c in zip(uniq, length, lo, hi)} + + +def _merge(total: dict[int, list[float]], part: dict[int, list[float]]) -> None: + for k, (length, lo, hi) in part.items(): + if k in total: + t = total[k] + t[0] += length + t[1] = min(t[1], lo) + t[2] = max(t[2], hi) + else: + total[k] = [length, lo, hi] + + +def _gpkg_rows(gpkg: Path): + """전국 등고선 가닥 묶음 — (WKB 목록 · 높이) · 읽기 전용.""" + con = sqlite3.connect(f"file:{gpkg.resolve().as_posix()}?mode=ro", uri=True) + try: + last = -1 + while True: + rows = con.execute( + "SELECT fid, geom, CTRLN_HG FROM contours WHERE fid > ? AND geom IS NOT NULL" + " ORDER BY fid LIMIT ?", + (last, INDEX_CHUNK), + ).fetchall() + if not rows: + return + last = rows[-1][0] + wkbs = [g[8 + _ENVELOPE[(g[3] >> 1) & 7] :] for _, g, _ in rows] + yield wkbs, np.array([h if h is not None else np.nan for _, _, h in rows]) + finally: + con.close() + + +def sheet_index(rebuild: bool = False) -> dict[str, list[float]]: + """{도엽: [등고선 길이 km, 최저 m, 최고 m]} — 전국 파일 지문이 같으면 저장본.""" + path = store.TERRAIN_ROOT / INDEX_NAME + stamp = store._gpkg_stamp() + if not rebuild: + try: + saved = json.loads(path.read_text(encoding="utf-8")) + if saved.get("gpkg") == stamp or stamp is None: + return saved["sheets"] + except (OSError, ValueError, KeyError): + pass + if stamp is None: + raise SystemExit(f"전국 등고선 파일 없음: {store.NATIONAL_GPKG}") + started = time.monotonic() + cells: dict[int, list[float]] = {} + for n, (wkbs, heights) in enumerate(_gpkg_rows(store.NATIONAL_GPKG), 1): + _merge(cells, aggregate(wkbs, heights)) + elapsed = time.monotonic() - started + print(f"도엽 목록 {n * INDEX_CHUNK:,} 가닥 · {elapsed:.0f} s", flush=True) + step = SHEET5K_SIZE_DEG + sheets = { + latlon_to_sheet5k((k // 100_000 + 0.5) * step, (k % 100_000 + 0.5) * step): [ + round(v[0] / 1000, 2), + v[1], + v[2], + ] + for k, v in cells.items() + if math.isfinite(v[1]) + } + atomic_write_json(path, {"gpkg": stamp, "sheets": sheets}) + return sheets + + +def pick_sheets( + index: dict[str, list[float]], + *, + sheets: list[str] | None = None, + box: tuple[float, float, float, float] | None = None, + min_relief: float = 150.0, + limit: int | None = None, +) -> list[str]: + """고른 도엽 — 등고선 길이 큰 것부터(같으면 높이 차) · 전국 파일이 안 덮는 도엽은 뺌.""" + if sheets: + chosen = [no for no in index if any(no.startswith(p) for p in sheets)] + else: + chosen = [no for no in index if index[no][2] - index[no][1] >= min_relief] + if box is not None: + lon0, lat0, lon1, lat1 = box + inside = [] + for no in chosen: + lat, lon = sheet5k_center(no) + if lon0 <= lon <= lon1 and lat0 <= lat <= lat1: + inside.append(no) + chosen = inside + chosen.sort(key=lambda no: (-index[no][0], -(index[no][2] - index[no][1]), no)) + return chosen[:limit] if limit else chosen + + +# ── 셈(일꾼) ───────────────────────────────────────────────────────────────── + + +def _mtime(path: Path) -> int | None: + try: + return path.stat().st_mtime_ns + except OSError: + return None + + +def run_step(step: str, sheet_no: str) -> tuple[str, str, dict[str, Any]]: + """도엽 하나 · 단계 하나 — 공용 함수만 부름 · fresh = 이번에 새로 셈.""" + epsg = epsg_of(sheet_no) + started = time.monotonic() + try: + if step == "dem": + fresh = store._stored_ok(sheet_no, epsg, CELL_M) is None + empty = store.ensure_sheet_terrain(sheet_no, epsg, CELL_M) is None + elif step == "zero": + from D03_PlanAnalysis import D03_PlanAnalysis_ZeroOrder as zero + + for no in neighbors_3x3(sheet_no): # 여유 200 m 이웃 지형(누적 열쇠가 나중에 안 바뀌게) + store.ensure_sheet_terrain(no, epsg, CELL_M) + target = zero._acc_names(sheet_no, epsg, CELL_M)[0] + before = _mtime(target) + empty = zero.sheet_accumulation(sheet_no, epsg, CELL_M) is None + fresh = _mtime(target) != before + elif step == "e07": + from E07_StabArea import E07_StabArea_Router as e07_router + from E07_StabArea import E07_StabArea_Sheet as e07 + + paths = e07._paths(sheet_no, epsg, CELL_M) + before = (_mtime(paths["flow"]), _mtime(paths["si"])) + empty = e07.ensure_sheet(sheet_no, epsg, CELL_M, e07_router.default_params()) is None + fresh = (_mtime(paths["flow"]), _mtime(paths["si"])) != before + else: + raise ValueError(f"모르는 단계: {step}") + except Exception as exc: # noqa: BLE001 — 한 장이 죽어도 나머지는 감 · 기록에 까닭 + return step, sheet_no, {"ok": False, "error": f"{type(exc).__name__}: {exc}"[:300]} + seconds = round(time.monotonic() - started, 2) + return step, sheet_no, {"ok": True, "s": seconds, "fresh": fresh, "empty": empty} + + +# ── 진행 기록 · 정리 ───────────────────────────────────────────────────────── + + +def load_progress() -> dict[str, Any]: + path = store.TERRAIN_ROOT / PROGRESS_NAME + stamp = store._gpkg_stamp() + try: + saved = json.loads(path.read_text(encoding="utf-8")) + except (OSError, ValueError): + saved = None + if not saved or saved.get("gpkg") != stamp or saved.get("cell_m") != CELL_M: + saved = {"gpkg": stamp, "cell_m": CELL_M, "steps": {}} + for step in STEPS: + saved["steps"].setdefault(step, {}) + return saved + + +def save_progress(progress: dict[str, Any]) -> None: + atomic_write_json(store.TERRAIN_ROOT / PROGRESS_NAME, progress) + + +def _alive(pid: int) -> bool: + import psutil + + return psutil.pid_exists(pid) + + +def clean_dead(root: Path | None = None) -> int: + """주인이 죽은 잠금 · 오래된 임시 파일을 지움(끊긴 일괄 셈이 남긴 것) → 지운 수.""" + root = root or store.TERRAIN_ROOT + removed = 0 + if not root.is_dir(): + return 0 + now = time.time() + for folder in root.iterdir(): + if not folder.is_dir(): + continue + for path in folder.glob(".*"): + try: + if path.suffix == ".lock": + pid = int(path.read_text(encoding="utf-8").split()[0]) + if _alive(pid): + continue + elif path.suffix != ".tmp" or now - path.stat().st_mtime < store.LOCK_STALE_S: + continue + path.unlink(missing_ok=True) + removed += 1 + except (OSError, ValueError, IndexError): + continue + return removed + + +def run( + sheets: list[str], + steps: list[str], + workers: int, + progress: dict[str, Any], + say=print, +) -> dict[str, Any]: + """단계 차례대로 · 단계 안은 일꾼 여럿 · 성공 기록된 도엽은 건너뜀 → 단계별 벽시계 초.""" + walls: dict[str, float] = {} + for step in [s for s in STEPS if s in steps]: + done = progress["steps"][step] + todo = [no for no in sheets if not done.get(no, {}).get("ok")] + say(f"[{step}] {len(todo)} 장 (건너뜀 {len(sheets) - len(todo)})") + started = time.monotonic() + saved_at = started + count = 0 + + def record(result: tuple[str, str, dict[str, Any]]) -> None: + nonlocal count, saved_at + _, no, info = result + done[no] = info + count += 1 + state = f"{info['s']} s" if info["ok"] else f"실패 {info['error']}" + say(f"[{step}] {count}/{len(todo)} {no} {state}") + if time.monotonic() - saved_at > SAVE_EVERY_S: + save_progress(progress) + saved_at = time.monotonic() + + try: + if workers <= 1: + for no in todo: + record(run_step(step, no)) + else: + with ProcessPoolExecutor(workers, initializer=below_normal) as pool: + futures = [pool.submit(run_step, step, no) for no in todo] + try: + for future in as_completed(futures): + record(future.result()) + except BaseException: + pool.shutdown(wait=False, cancel_futures=True) + raise + finally: + save_progress(progress) + walls[step] = round(time.monotonic() - started, 1) + return walls + + +def summary(sheets: list[str], progress: dict[str, Any]) -> dict[str, Any]: + """단계별 — 성공 · 새로 셈 · 실패 · 새로 셈 장당 초(가운데 · 평균 · 최대) · 장당 크기 MB.""" + out: dict[str, Any] = {} + for step in STEPS: + done = progress["steps"][step] + infos = [done[no] for no in sheets if no in done] + failed = {no: done[no]["error"] for no in sheets if done.get(no, {}).get("ok") is False} + fresh = [i["s"] for i in infos if i.get("ok") and i.get("fresh")] + sizes = [] + for no in sheets: + if not done.get(no, {}).get("ok"): + continue + stem = f"{epsg_of(no)}_{CELL_M:g}m" + folder = store.TERRAIN_ROOT / no + total = 0 + for head in STEP_FILES[step]: + total += sum(p.stat().st_size for p in folder.glob(f"{head}{stem}*")) + if total: + sizes.append(total / 1e6) + out[step] = { + "ok": sum(1 for i in infos if i.get("ok")), + "fresh": len(fresh), + "failed": failed, + "s_median": round(float(np.median(fresh)), 1) if fresh else None, + "s_mean": round(float(np.mean(fresh)), 1) if fresh else None, + "s_max": round(float(np.max(fresh)), 1) if fresh else None, + "mb_mean": round(float(np.mean(sizes)), 2) if sizes else None, + "mb_total": round(float(np.sum(sizes)), 1) if sizes else 0.0, + } + return out + + +def main(argv: list[str] | None = None) -> int: + parser = argparse.ArgumentParser(description="공용 지형 창고 전국 미리 계산(PLAN 612)") + parser.add_argument("--steps", default="dem,zero,e07") + parser.add_argument("--sheets", default="") + parser.add_argument("--box", default="") + parser.add_argument("--min-relief", type=float, default=150.0) + parser.add_argument("--limit", type=int, default=0) + parser.add_argument("--workers", type=int, default=4) + parser.add_argument("--status", action="store_true") + parser.add_argument("--rebuild-index", action="store_true") + args = parser.parse_args(argv) + steps = [s.strip() for s in args.steps.split(",") if s.strip()] + unknown = [s for s in steps if s not in STEPS] + if unknown: + parser.error(f"모르는 단계 {unknown} — 지금 되는 단계 {STEPS}(tiles 는 다음)") + box = tuple(float(v) for v in args.box.split(",")) if args.box else None + if box is not None and len(box) != 4: + parser.error("--box = lon0,lat0,lon1,lat1") + + below_normal() + index = sheet_index(args.rebuild_index) + wanted = [s.strip() for s in args.sheets.split(",") if s.strip()] + sheets = pick_sheets( + index, sheets=wanted or None, box=box, min_relief=args.min_relief, limit=args.limit or None + ) + print(f"도엽 {len(sheets)} 장 (목록 {len(index)} 장) · 단계 {steps} · 일꾼 {args.workers}") + progress = load_progress() + if not args.status: + print(f"치운 잠금 · 임시 파일 {clean_dead()}") + started = time.monotonic() + walls = run(sheets, steps, args.workers, progress) + print(f"벽시계 {walls} · 전체 {time.monotonic() - started:.1f} s") + print(json.dumps(summary(sheets, progress), ensure_ascii=False, indent=1)) + return 0 + + +if __name__ == "__main__": + raise SystemExit(main())