Files
Aislo/D03_PlanAnalysis/D03_PlanAnalysis_Terrain.py
T
eomsangdonandClaude Opus 5.5 0c463bcc64 feat(D 예정노선): 지형 받기 서버 — 지적 뺌 · 실패 조각만 · 위치 바꿈 멀다 · 띠 · reset · 칸 상한 · 끊김 · 단계 글 · 임도망 ±5 km (PLAN 595 · 요구 463)
- 지적도를 지형 일에서 뺌(겹 GET 이 받음) · 지적 실패 조각만 다시(열쇠 failed_tiles)
- 받기 몸 center · reset — 5 km 넘음 far · 다른 TM 띠 band 물음 · reset = 옛 범위 · 병합 도엽 · 위성 버림
- 받을 칸 500 만 넘음 413 · extent max_cells · lon/lat
- 진행 기록 파일 — 서버 다시 켜짐 = interrupted(percent < 100 · message 없음) · 받은 적 없음 percent 100
- 진행 단계 글 n/5 · 받은 범위에 epsg · center 기록
- 임도망 자르기 = 기준점 ±5 km ∪ 받은 지형 · 지형 전에도 올림 · bounds · distance_m · 밖뿐 = 거리 글
- 시험 530_1 · 530_2 · 569_1 · 584_3 덧

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-30 14:45:26 +09:00

346 lines
14 KiB
Python
Raw Blame History

This file contains ambiguous Unicode characters
This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.
"""D 예정노선 지형 받기 — 사업 전체 한 벌(도엽 · 위성 · DEM) (요구 345 · 389 · 463 · PLAN 530-1 ·
561-2 · 595).
지형은 D 프로젝트 storage 안 `B04_PreProcess/` 배치 그대로(지형 코드 주인 = B04 · 고칠 곳 0):
- 범위 = 전체 시점 · 종점(+ 그은 점) 상자 + 여유 · 받은 범위가 있으면 그 범위와 합침(줄지 않음) ·
새로(`reset` · 위치를 멀리 옮김) = 옛 받은 범위 · 병합 도엽 · 위성 범위를 버리고 새 상자만 ·
칸 상한 `MAX_CELLS`(넘으면 라우터가 거절)
- 도엽 = 범위 상자에 닿는 도엽 모두 → `ensure_sheets`(없는 것만 · prune 없음)
- 병합 = 지난번 병합 도엽 ∪ 이번 도엽(옛 도엽 유지 · 열쇠 같으면 건너뜀)
- 위성 = 범위 도엽 도곽 ∪ 옛 위성 범위(덮으면 건너뜀)
- 지적도 = 지형 일에서 뺌 — 겹을 켤 때 GET 이 받음(`D03_PlanAnalysis_Layers.ensure_cadastral`)
- DEM = 한 벌 `models_plan/`(연도는 D03 나눔 뒤에 정해짐) · 받은 범위 = `plan_terrain.json`
(상자 · 여유 · 좌표계 · 위치 칸 가운데)
- 작업 좌표계 = 지점 epsg · `B03_FileInput/input/prj/plan.prj` 한 장(B04 · B05 좌표계 창구가 읽음)
- 받기 전 범위 보기 = `plan_extent`(내려받기 없음) · 받기는 사용자가 누를 때만
계산 자리 = ③ 서버 단독. DB 행(`surface_models`)은 라우터가 씀.
"""
from __future__ import annotations
import json
import logging
import math
from pathlib import Path
from typing import Any, Callable
import numpy as np
from B04_PreProcess.B04_PreProcess_Engine_MapSheet import (
SHEET5K_SIZE_DEG,
sheet5k_to_bounds,
sheets_for_points,
)
from B04_PreProcess.B04_PreProcess_Engine_SheetParser import PARSER_VERSION, parse_and_merge_sheets
from B04_PreProcess.B04_PreProcess_Engine_SheetStore import (
ensure_sheets,
get_project_map_sheets_dir,
get_sheet_path,
)
from B04_PreProcess.B04_PreProcess_Engine_SheetSurface import build_sheet_surface_model
from B04_PreProcess.B04_PreProcess_Engine_VWorld import download_vworld_satellite_map
from common_util.common_util_json import atomic_write_json
from config.config_system_terrain import (
PLAN_TERRAIN_GRID_M,
SHEET_SURFACE_DEFAULT_METHOD,
SHEET_SURFACE_MARGIN_M,
)
logger = logging.getLogger(__name__)
# 여유 처음 값 = 옛 연도 구간 DEM 폭(탐색 복도 300 m + 횡단 표본 50 m) · 바닥 = B 절취 여유
DEFAULT_MARGIN_M = 350.0
MAX_MARGIN_M = 5000.0
# 처음 받는 범위 = 등록 기준점 둘레 ± 1 km(요구 453 · PLAN 591-2)
SITE_RADIUS_M = 1000.0
# 위치 바꿈이 받은 범위에서 이만큼 넘게 멀면 묻고 새로 · 받을 DEM 칸 상한(요구서_463 3장 4)
FAR_M = 5000.0
MAX_CELLS = 5_000_000
STAGE = Path("B04_PreProcess")
PLAN_MODELS = "models_plan"
TERRAIN_META = "plan_terrain.json"
PRJ_NAME = "plan.prj"
Progress = Callable[[int, str], None]
Box = tuple[float, float, float, float] # x0, y0, x1, y1(작업 좌표계 m)
XY = tuple[float, float]
def models_dir_of(project_root: Path) -> Path:
"""사업 전체 DEM 자리."""
return project_root / STAGE / PLAN_MODELS
def clamp_margin(margin_m: float | None) -> float:
"""여유 — 빈 값 = 처음 값 · B 절취 여유 ~ 5 km."""
if margin_m is None or not math.isfinite(margin_m):
return DEFAULT_MARGIN_M
return min(MAX_MARGIN_M, max(SHEET_SURFACE_MARGIN_M, float(margin_m)))
def received_meta(project_root: Path) -> dict[str, Any]:
"""받은 범위 기록(`box` · `margin_m` · `epsg` · `center`) — 없으면 빈 dict."""
return _read_json(models_dir_of(project_root) / TERRAIN_META) or {}
def received_box(project_root: Path, epsg: int | None = None) -> Box | None:
"""받은 DEM 범위 — 받은 적 없거나 읽지 못하거나 다른 좌표계(`epsg` 줄 때)면 None."""
meta = received_meta(project_root)
if epsg is not None and meta.get("epsg") not in (None, int(epsg)):
return None
try:
x0, y0, x1, y1 = (float(v) for v in meta["box"])
except (KeyError, TypeError, ValueError):
return None
return x0, y0, x1, y1
def site_corners(x: float, y: float) -> list[XY]:
"""기준점 둘레 상자 두 귀 — 처음 여유를 더하면 ± `SITE_RADIUS_M`."""
r = SITE_RADIUS_M - DEFAULT_MARGIN_M
return [(x - r, y - r), (x + r, y + r)]
def wanted_box(points: list[XY], margin_m: float, have: Box | None) -> Box:
"""받을 범위 = 점 상자 + 여유 ∪ 받은 범위."""
xs, ys = [p[0] for p in points], [p[1] for p in points]
box = (min(xs) - margin_m, min(ys) - margin_m, max(xs) + margin_m, max(ys) + margin_m)
if have is None:
return box
return min(box[0], have[0]), min(box[1], have[1]), max(box[2], have[2]), max(box[3], have[3])
def cells_of(box: Box) -> int:
"""DEM 칸 수 어림 — 받기 시간 · 메모리 가늠 · 상한 판정."""
return int((box[2] - box[0]) * (box[3] - box[1]) / PLAN_TERRAIN_GRID_M**2)
def distance_to_box(xy: XY, box: Box) -> float:
"""점 ~ 상자 거리(m · 안이면 0)."""
dx = max(box[0] - xy[0], 0.0, xy[0] - box[2])
dy = max(box[1] - xy[1], 0.0, xy[1] - box[3])
return math.hypot(dx, dy)
def covers(have: Box | None, box: Box) -> bool:
return (
have is not None
and all(have[i] <= box[i] for i in (0, 1))
and all(have[i] >= box[i] for i in (2, 3))
)
def box_ring(box: Box, epsg: int) -> list[list[float]]:
"""상자 네 귀 경위도([경도, 위도] · 왼아래부터 반시계)."""
from common_util.common_util_crs import tm_to_lonlat
x0, y0, x1, y1 = box
return [list(tm_to_lonlat(x, y, epsg)) for x, y in ((x0, y0), (x1, y0), (x1, y1), (x0, y1))]
def sheets_for_box(box: Box, epsg: int) -> list[str]:
"""상자에 닿는 도엽 — 경위도 외접 상자를 도엽 반 칸 간격으로 짚음."""
ring = box_ring(box, epsg)
lons, lats = [p[0] for p in ring], [p[1] for p in ring]
step = SHEET5K_SIZE_DEG / 2
def ticks(lo: float, hi: float) -> list[float]:
return [lo + step * i for i in range(int((hi - lo) / step) + 1)] + [hi]
grid = [
(lat, lon) for lat in ticks(min(lats), max(lats)) for lon in ticks(min(lons), max(lons))
]
return sorted(set(sheets_for_points(grid)))
def plan_extent(
project_root: Path, points: list[XY], epsg: int, margin_m: float | None = None
) -> dict[str, Any]:
"""받기 전 범위 보기 — 내려받기 없음. `points` = 전체 시점 · 종점(+ 그은 점 · 작업 좌표계 m)."""
margin = clamp_margin(margin_m)
have = received_box(project_root, epsg)
box = wanted_box(points, margin, have)
store = get_project_map_sheets_dir(project_root)
sheets = []
for no in sheets_for_box(box, epsg):
lon0, lat0, lon1, lat1 = sheet5k_to_bounds(no)
sheets.append(
{
"no": no,
"ring": [[lon0, lat0], [lon1, lat0], [lon1, lat1], [lon0, lat1]],
"stored": get_sheet_path(store, no) is not None,
}
)
return {
"epsg": int(epsg),
"margin_m": margin,
"box_m": list(box),
"box": box_ring(box, epsg),
"have_m": list(have) if have else None,
"have": box_ring(have, epsg) if have else None,
"covered": covers(have, box),
"sheets": sheets,
"cells": cells_of(box),
"max_cells": MAX_CELLS,
}
def write_plan_prj(project_root: Path, epsg: int) -> Path:
"""작업 좌표계 .prj — 내용이 같으면 그대로."""
from pyproj import CRS
path = project_root / "B03_FileInput" / "input" / "prj" / PRJ_NAME
text = CRS.from_epsg(int(epsg)).to_wkt("WKT1_ESRI")
if not path.is_file() or path.read_text(encoding="utf-8", errors="ignore") != text:
path.parent.mkdir(parents=True, exist_ok=True)
path.write_text(text, encoding="utf-8")
return path
def _read_json(path: Path) -> dict[str, Any] | None:
try:
return json.loads(path.read_text(encoding="utf-8"))
except (OSError, ValueError):
return None
def _merge_sheets(
processed: Path, store: Path, available: list[str], keep_old: bool = True
) -> list[str]:
"""병합 = 지난번 병합 도엽(zip 남은 것 · 새로면 뺌) ∪ 이번 도엽 · 열쇠 같으면 건너뜀."""
key_path = processed / "sheet_merge.key.json"
old = ((_read_json(key_path) or {}).get("sheets") or []) if keep_old else []
sheets = sorted({*available, *(s for s in old if get_sheet_path(store, s) is not None)})
key = {"sheets": sheets, "parser": PARSER_VERSION}
if sheets and (_read_json(key_path) != key or not any(processed.glob("도엽_*.geojson"))):
parse_and_merge_sheets(store, sheets, processed)
atomic_write_json(key_path, key)
return sheets
def _satellite(
processed: Path,
prj: Path,
points: list[tuple[float, float]],
epsg: int,
keep_old: bool = True,
) -> str:
"""위성 = 기준 도엽 도곽 ∪ 옛 범위(새로면 뺌) · 이미 덮으면 건너뜀. 돌려줌 = skip · done ·
failed."""
from pyproj import Transformer
to_metric = Transformer.from_crs("EPSG:4326", f"EPSG:{epsg}", always_xy=True)
xs: list[float] = []
ys: list[float] = []
for sheet in sheets_for_points(points):
lon0, lat0, lon1, lat1 = sheet5k_to_bounds(sheet)
for lon, lat in ((lon0, lat0), (lon0, lat1), (lon1, lat0), (lon1, lat1)):
x, y = to_metric.transform(lon, lat)
xs.append(x)
ys.append(y)
need = {"x": [min(xs), max(xs)], "y": [min(ys), max(ys)]}
meta = _read_json(processed / "vworld_satellite_meta.json") if keep_old else None
try:
old = [float(meta[k]) for k in ("x_min", "x_max", "y_min", "y_max")] if meta else None
except (KeyError, TypeError, ValueError):
old = None
if old and old[0] <= need["x"][0] and old[1] >= need["x"][1]:
if old[2] <= need["y"][0] and old[3] >= need["y"][1]:
return "skip"
if old:
need = {
"x": [min(old[0], need["x"][0]), max(old[1], need["x"][1])],
"y": [min(old[2], need["y"][0]), max(old[3], need["y"][1])],
}
try:
done = download_vworld_satellite_map(prj, need, processed, layer_name="Satellite")
except Exception as exc: # noqa: BLE001 — 위성은 배경 · 지형 받기를 막지 않음
logger.warning("D 위성 지도 받기 실패: %s", exc)
return "failed"
return "done" if done else "failed"
def ensure_plan_terrain(
project_root: Path,
points: list[XY],
epsg: int,
margin_m: float | None = None,
on_progress: Progress | None = None,
reset: bool = False,
center: tuple[float, float] | None = None,
) -> dict[str, Any]:
"""사업 전체 지형 받기 — `points` = 전체 시점 · 종점(+ 그은 점 · 작업 좌표계 m) ·
`reset` = 옛 받은 범위 버리고 새로 · `center` = 위치 칸 가운데(경위도 · 기록만).
돌려줌 = 도엽(받은 · 실패 · 병합) · 위성 상태 · 받은 범위 · 모델 등록용 dict 목록.
"""
def report(percent: int, message: str) -> None:
if on_progress is not None:
on_progress(percent, message)
if not points:
raise ValueError("시점 · 종점이 없습니다.")
processed = project_root / STAGE / "processed"
processed.mkdir(parents=True, exist_ok=True)
prj = write_plan_prj(project_root, epsg)
margin = clamp_margin(margin_m)
have = None if reset else received_box(project_root, epsg)
box = wanted_box(points, margin, have)
report(5, "1/5 수치지형도 도엽 고르는 중")
wanted = sheets_for_box(box, epsg)
store = get_project_map_sheets_dir(project_root)
before = {s for s in wanted if get_sheet_path(store, s) is not None}
report(10, f"2/5 수치지형도 도엽 {len(wanted) - len(before)}매 받는 중")
got = ensure_sheets(store, wanted)
if not got["available"]:
raise ValueError(f"수치지형도 도엽을 하나도 받지 못했습니다: {got['failed']}")
report(45, "3/5 도엽 레이어 병합 중")
merged = _merge_sheets(processed, store, list(got["available"]), keep_old=have is not None)
report(60, "4/5 위성 지도 받는 중")
ring = [(p[1], p[0]) for p in box_ring(box, epsg)]
satellite = _satellite(processed, prj, ring, epsg, keep_old=have is not None)
report(65, f"5/5 사업 지형면 만드는 중 — 칸 {cells_of(box) / 1e4:.0f}만 · 수 분 걸릴 수 있음")
models_dir = models_dir_of(project_root)
models_dir.mkdir(parents=True, exist_ok=True)
# B 절취 = 넘긴 꼭짓점 상자 + `SHEET_SURFACE_MARGIN_M` → 그만큼 안쪽으로 넘김(여유 ≥ 그 값)
inner = SHEET_SURFACE_MARGIN_M
outline = np.array(
[[box[0] + inner, box[1] + inner], [box[2] - inner, box[3] - inner]], dtype=np.float64
)
models = build_sheet_surface_model(
project_root,
processed,
models_dir,
outline,
f"EPSG:{int(epsg)}",
[SHEET_SURFACE_DEFAULT_METHOD],
PLAN_TERRAIN_GRID_M,
)
if not models:
raise ValueError(
"사업 지형면을 만들지 못했습니다 — 범위 안 도엽 등고선 높이 확인"
f"(병합 도엽 {len(merged)}매)"
)
meta: dict[str, Any] = {"box": list(box), "margin_m": margin, "epsg": int(epsg)}
kept = None if reset else received_meta(project_root).get("center")
if center or kept:
meta["center"] = list(center or kept)
atomic_write_json(models_dir / TERRAIN_META, meta)
report(100, "지형 받기 끝")
return {
"sheets": {
"wanted": len(wanted),
"downloaded": sorted(set(got["available"]) - before),
"failed": got["failed"],
"merged": merged,
},
"satellite": satellite,
"box_m": list(box),
"models": models,
}