Aislo/B07_DesignDetail/B07_DesignDetail_Router_Support_Basin.py

"""B07 유역도 배경 자료 — GeoJSON 읽기·좌표계 환산·도곽 클리핑.

지원 모듈이 700줄을 넘어 떼어냈다(2026-09-04). 도면 목록·확정 저장은 본체에 남는다.
"""

import json
import logging
import math
import re
from functools import lru_cache
from pathlib import Path
from typing import Any

from pyproj import Transformer

from B04_PreProcess.B04_PreProcess_Router_Watershed import CONTOUR_FILE, STREAM_FILE
from B07_DesignDetail.B07_DesignDetail_Engine_Cad_Basin import map_area_mm
from B07_DesignDetail.B07_DesignDetail_Engine_Cad_Plan import (
    plan_area_mm,
    plan_chunks,
    plan_drawing_label,
)
from common_util.common_util_drainage_pipes import detail_basins_path
from config.config_system import DRAWING_SCALE_BASIN, DRAWING_SCALE_PLAN

logger = logging.getLogger(__name__)


# 도면 id 상수·빈 도면 목록은 `B07_DesignDetail_Router_Support` 한 곳이 정본이다.
# 이 모듈에 있던 같은 이름의 사본은 아무도 읽지 않으면서 값만 어긋나 지웠다(2026-09-04).


def _geojson_payload(path: Path) -> dict[str, Any]:
    """GeoJSON 전체를 읽는다. 파일이 없거나 깨졌으면 빈 dict."""
    if not path.is_file():
        return {}
    try:
        payload = json.loads(path.read_text(encoding="utf-8"))
    except (OSError, ValueError):
        logger.warning("B07 유역도: GeoJSON을 읽지 못했습니다 — %s", path)
        return {}
    return payload if isinstance(payload, dict) else {}


def _geojson_features(path: Path) -> list[dict[str, Any]]:
    """GeoJSON 피처 목록만 읽는다. 파일이 없거나 깨졌으면 빈 목록."""
    features = _geojson_payload(path).get("features")
    return features if isinstance(features, list) else []


def _basins_crs(context: Any, payload: dict[str, Any]) -> str:
    """저장본을 미터로 되돌릴 좌표계.

    저장할 때 쓴 좌표계를 파일이 갖고 있으면 그 값이다. 없으면 좌표계 기록 이전에 저장된
    파일이라 노선 CSV의 EPSG 라벨로 쓰였다 — 그 라벨로 되돌려야 왕복이 맞는다
    (2026-09-01). 라벨을 못 읽으면 현재 사업지 좌표계로 둔다.
    """
    from common_util.common_util_route_geometry import (
        find_planned_route_file,
        read_planned_route,
    )

    stored = payload.get("crs_input")
    if isinstance(stored, str) and stored:
        return stored
    route_file = find_planned_route_file(context.project_root / "B03_FileInput" / "input")
    planned = read_planned_route(route_file) if route_file else None
    if planned is not None and planned.epsg:
        logger.warning(
            "B07 유역도: 좌표계 기록이 없는 옛 저장본 — 노선 CSV 라벨 EPSG:%s로 되돌립니다. "
            "B04에서 유역을 다시 확정하면 현재 좌표계로 새로 남습니다.",
            planned.epsg,
        )
        return f"EPSG:{planned.epsg}"
    return context.crs


def _geometry_lines(geometry: Any) -> list[list[tuple[float, float]]]:
    """LineString·MultiLineString·Polygon을 점열 목록으로 편다."""
    if not isinstance(geometry, dict):
        return []
    kind = geometry.get("type")
    coordinates = geometry.get("coordinates")
    if not isinstance(coordinates, list) or not coordinates:
        return []
    if kind == "LineString":
        return [[(float(p[0]), float(p[1])) for p in coordinates if isinstance(p, list)]]
    if kind in ("MultiLineString", "Polygon"):
        return [
            [(float(p[0]), float(p[1])) for p in part if isinstance(p, list)]
            for part in coordinates
            if isinstance(part, list)
        ]
    return []


def _clip_segment(
    start: tuple[float, float],
    end: tuple[float, float],
    box: tuple[float, float, float, float],
) -> tuple[tuple[float, float], tuple[float, float]] | None:
    """선분에서 상자 안에 드는 부분만 돌려준다(Liang-Barsky). 겹치지 않으면 None."""
    min_x, min_y, max_x, max_y = box
    dx = end[0] - start[0]
    dy = end[1] - start[1]
    t0, t1 = 0.0, 1.0
    for numerator, denominator in (
        (start[0] - min_x, -dx),
        (max_x - start[0], dx),
        (start[1] - min_y, -dy),
        (max_y - start[1], dy),
    ):
        if denominator == 0.0:
            if numerator < 0.0:
                return None  # 경계와 나란한데 바깥이다
            continue
        t = numerator / denominator
        if denominator < 0.0:
            if t > t1:
                return None
            t0 = max(t0, t)
        else:
            if t < t0:
                return None
            t1 = min(t1, t)
    return (
        (start[0] + t0 * dx, start[1] + t0 * dy),
        (start[0] + t1 * dx, start[1] + t1 * dy),
    )


def clip_line_to_box(
    line: list[tuple[float, float]], box: tuple[float, float, float, float]
) -> list[list[tuple[float, float]]]:
    """범위 안에 든 부분만 조각으로 잘라 낸다 — 끝점은 경계 위에 정확히 놓인다.

    도엽 등고선 한 줄은 도엽 끝까지 이어지므로, 걸치면 통째로 남기면 도면이 A1을 넘는다.
    예전 구현은 경계 바깥 점을 하나씩 물고 나와 외곽선이 들쭉날쭉했다(2026-08-31 사용자
    지적). 이제 교차점을 계산해 자르므로 절취면이 곧게 떨어진다.
    """
    runs: list[list[tuple[float, float]]] = []
    current: list[tuple[float, float]] = []
    for start, end in zip(line, line[1:]):
        piece = _clip_segment(start, end, box)
        if piece is None:
            current = []
            continue
        head, tail = piece
        if head == tail:
            continue  # 경계에 점으로만 닿았다
        if current and current[-1] == head:
            current.append(tail)
        else:
            current = [head, tail]
            runs.append(current)
    return [run for run in runs if len(run) >= 2]


# 세부유역이 노선에서 이만큼 넘게 떨어져 있으면 좌표계를 잘못 되돌린 것으로 본다.
# 임도 한 노선이 담는 유역은 길어야 수 km라 오검출 여지가 없다.
_BASIN_MAX_DISTANCE_M = 50_000.0


def _too_far_from_route(
    ring: list[tuple[float, float]], route_xy: list[tuple[float, float]]
) -> bool:
    """되돌린 유역이 노선 근처에 없으면 True — 좌표계를 되찾지 못한 저장본이다."""
    if not ring or not route_xy:
        return False
    cx = sum(x for x, _ in ring) / len(ring)
    cy = sum(y for _, y in ring) / len(ring)
    return min(math.hypot(cx - x, cy - y) for x, y in route_xy) > _BASIN_MAX_DISTANCE_M


@lru_cache(maxsize=8)
def _metric_lines_cached(
    path_str: str, mtime_ns: int, crs: str
) -> tuple[tuple[tuple[float, float], ...], ...]:
    """도엽 GeoJSON 한 벌을 사업지 좌표계(m) 선 목록으로 돌려 **캐시**한다.

    같은 배경을 유역도와 계획평면도가 나눠 쓴다 — 도면마다 다시 읽고 다시 투영하면
    한 장 여는 데 수 초가 걸린다(등고선 수만 점). 파일이 바뀌면 mtime 이 달라져
    캐시가 저절로 갈린다(`_read_template` 와 같은 방식).
    """
    to_metric = Transformer.from_crs("EPSG:4326", crs, always_xy=True)
    lines: list[tuple[tuple[float, float], ...]] = []
    for feature in _geojson_features(Path(path_str)):
        for part in _geometry_lines(feature.get("geometry")):
            converted = tuple(
                (float(x), float(y))
                for x, y in (to_metric.transform(point[0], point[1]) for point in part)
            )
            if len(converted) >= 2:
                lines.append(converted)
    return tuple(lines)


def _metric_lines(path: Path, crs: str) -> list[list[tuple[float, float]]]:
    """캐시된 배경 선을 쓰기 좋은 형태로 낸다. 파일이 없으면 빈 목록."""
    if not path.is_file():
        return []
    cached = _metric_lines_cached(str(path), path.stat().st_mtime_ns, crs)
    return [list(line) for line in cached]


def map_background(
    project_root: Path,
    crs: str,
    scale: int,
    area_mm: tuple[float, float],
    extent_points: list[tuple[float, float]],
) -> dict[str, list[list[tuple[float, float]]]]:
    """도엽 등고선·세류선을 사업지 좌표계로 읽어 **그 도면의 도곽 크기로 절취**한다.

    유역도·계획평면도·용지도가 같은 창구를 쓴다 — 자료 읽기·좌표 환산은 한 번뿐이고
    (`_metric_lines` 캐시), 도면마다 다른 것은 축척과 도곽 크기뿐이다.

    `extent_points` 는 그 도면의 주제(노선·유역 등) 좌표다. 도곽보다 크면 그쪽을
    우선한다 — 배경만 잘리고 주제는 다 보인다.
    """
    area_w_mm, area_h_mm = area_mm
    half_w_m = area_w_mm / 2.0 * scale / 1000.0
    half_h_m = area_h_mm / 2.0 * scale / 1000.0
    if extent_points:
        center_x = (min(x for x, _y in extent_points) + max(x for x, _y in extent_points)) / 2.0
        center_y = (min(y for _x, y in extent_points) + max(y for _x, y in extent_points)) / 2.0
        box = (
            min(center_x - half_w_m, min(x for x, _y in extent_points)),
            min(center_y - half_h_m, min(y for _x, y in extent_points)),
            max(center_x + half_w_m, max(x for x, _y in extent_points)),
            max(center_y + half_h_m, max(y for _x, y in extent_points)),
        )
    else:
        box = (-math.inf, -math.inf, math.inf, math.inf)

    sheet_dir = Path(project_root) / "B04_PreProcess" / "processed"
    background: dict[str, list[list[tuple[float, float]]]] = {}
    for key, filename in (("contours", CONTOUR_FILE), ("streams", STREAM_FILE)):
        lines: list[list[tuple[float, float]]] = []
        for line in _metric_lines(sheet_dir / filename, crs):
            lines.extend(clip_line_to_box(line, box))
        background[key] = lines
    return background


PLAN_ID = re.compile(r"^(plan_terrain|plan_route|plan_layout)(?:_(\d+))?$")


def plan_stations(longitudinal: dict[str, Any]) -> list[tuple[float, float, float]]:
    """종단 측점의 (누가거리 m, x, y) — 계획평면도 장 나눔의 유일한 기준 자료."""
    stations: list[tuple[float, float, float]] = []
    for station in longitudinal.get("stations") or []:
        if not isinstance(station, dict):
            continue
        chainage = station.get("chainage_m")
        x, y = station.get("center_x"), station.get("center_y")
        if all(isinstance(value, (int, float)) for value in (chainage, x, y)):
            stations.append((float(chainage), float(x), float(y)))
    return stations


def plan_chunk_for(
    longitudinal: dict[str, Any], drawing_id: str
) -> tuple[str, dict[str, Any], int]:
    """도면 id 에서 (주제, 그 장의 구간, 전체 장수)를 찾는다."""
    match = PLAN_ID.fullmatch(drawing_id)
    if not match:
        raise ValueError("올바르지 않은 계획평면도 ID입니다.")
    kind = match.group(1)
    chunks = plan_chunks(plan_stations(longitudinal))
    number = int(match.group(2)) if match.group(2) else 1
    chunk = next((item for item in chunks if item["number"] == number), None)
    if chunk is None:
        raise FileNotFoundError("요청한 계획평면도 장을 찾을 수 없습니다.")
    return kind, chunk, len(chunks)


def plan_source(context: Any, longitudinal: dict[str, Any], drawing_id: str) -> dict[str, Any]:
    """계획평면도 한 장의 입력(노선·측점·등고선·세류선·구조물)을 사업지 CRS(m)로 모은다.

    배경은 유역도와 **같은 창구**(`map_background`)를 쓴다 — 도엽 GeoJSON 읽기·좌표
    환산이 캐시돼 두 도면이 자료를 나눠 쓴다(2026-09-04 사용자 지시).

    장이 여럿이면 그 장의 구간(누가거리)에 드는 노선·구조물만 싣고, 배경도 그 범위로
    절취한다 — 축척 1/1,200 은 고정이므로 안 들어가면 장을 나눈다.
    """
    kind, chunk, total = plan_chunk_for(longitudinal, drawing_id)
    start_m, end_m = float(chunk["start_m"]), float(chunk["end_m"])

    route_xy: list[tuple[float, float]] = []
    for vertex in context.vertices:
        chainage = float(getattr(vertex, "chainage_m", 0.0) or 0.0)
        if total > 1 and not (start_m <= chainage <= end_m):
            continue
        route_xy.append((vertex.x, vertex.y))
    # 측점 눈금은 **종단 측점**을 쓴다 — 노선 정점은 조밀하고 누가거리가 간격의 배수가 아니다.
    stations = [
        station
        for station in plan_stations(longitudinal)
        if total <= 1 or start_m <= station[0] <= end_m
    ]

    structures = [
        structure
        for structure in _plan_structures(context)
        if total <= 1 or start_m <= float(structure.get("chainage_m") or -1.0) <= end_m
    ]
    background = map_background(
        Path(context.project_root),
        context.crs,
        DRAWING_SCALE_PLAN,
        plan_area_mm(),
        route_xy,
    )
    return {
        "kind": kind,
        "label": plan_drawing_label(kind, chunk, total),
        "route_xy": route_xy,
        "stations": stations,
        "contours": background["contours"],
        "streams": background["streams"],
        "structures": structures,
    }


def _plan_structures(context: Any) -> list[dict[str, Any]]:
    """배치도에 찍을 구조물 — B04 배수시설 정본(`pipe_points.json`)을 그대로 읽는다.

    좌표는 이미 사업지 CRS(m)다(B04가 그렇게 쓴다). 없으면 빈 목록 — 배치도는 배경과
    노선만으로도 열린다.
    """
    path = Path(context.project_root) / "B04_PreProcess" / "drainage" / "edits" / "pipe_points.json"
    points = _geojson_payload(path).get("points")
    if not isinstance(points, list):
        return []
    return [point for point in points if isinstance(point, dict)]


def watershed_source(context: Any) -> dict[str, Any]:
    """유역도 입력(노선·세부유역·등고선·세류선)을 사업지 CRS(m)로 모은다.

    저장본은 전부 WGS84라 여기서 미터 좌표로 되돌린다(B04가 저장할 때와 반대 방향).

    되돌리는 좌표계가 둘이다. **배경(도엽 등고선·세류선)은 사업지 좌표계**로 돌린다 —
    노선(`context.vertices`)이 그 좌표계에 있으므로 같은 자리에 겹쳐야 한다(그 환산은
    `map_background()` 안에 있다). **세부유역은 그 파일을 쓸 때 쓴 좌표계**로 돌린다 —
    좌표계 기록 이전 저장본은 노선 CSV의 EPSG 라벨로 쓰였고, 그 라벨로 되돌려야 원래
    미터 좌표가 나온다(2026-09-01).
    """
    basins_payload = _geojson_payload(detail_basins_path(context.stored_path))
    to_basin_metric = Transformer.from_crs(
        "EPSG:4326", _basins_crs(context, basins_payload), always_xy=True
    )

    def basin_metric(point: tuple[float, float]) -> tuple[float, float]:
        x, y = to_basin_metric.transform(point[0], point[1])
        return (float(x), float(y))

    route_xy = [(vertex.x, vertex.y) for vertex in context.vertices]
    basins: list[dict[str, Any]] = []
    dropped = 0
    for feature in basins_payload.get("features") or []:
        properties = feature.get("properties") or {}
        if properties.get("kind") != "detail_basin":
            continue
        rings = _geometry_lines(feature.get("geometry"))
        if not rings:
            continue
        ring = [basin_metric(point) for point in rings[0]]
        # 저장본 좌표계를 못 되찾으면 유역이 노선에서 수백 km 밖으로 떨어진다. 그대로 두면
        # 도곽이 그 거리까지 벌어져 도면이 통째로 빈 화면이 된다 — 유역만 버리고 배경·노선은
        # 그린다(2026-09-01 다른 PC 보고: 유역도 그림 자체가 없음).
        if _too_far_from_route(ring, route_xy):
            dropped += 1
            continue
        basins.append({"ring": ring, "props": properties})
    if dropped:
        logger.warning(
            "B07 유역도: 노선에서 %.0fkm 넘게 떨어진 세부유역 %d개를 뺐습니다 — "
            "저장본 좌표계를 되찾지 못했습니다. B04에서 유역을 다시 확정하세요.",
            _BASIN_MAX_DISTANCE_M / 1000.0,
            dropped,
        )

    # 배경(등고선·세류선)은 **여러 도엽을 합쳐 받은 뒤 도곽 크기로 절취**한다
    # (2026-08-30 사용자 지시 — 노선이 도엽 경계에 걸릴 수 있어 주변 도엽까지 받아 둔다).
    # 읽기·환산·절취는 `map_background()` 한 곳에 있고 계획평면도·용지도도 같은 것을 쓴다.
    background = map_background(
        Path(context.project_root),
        context.crs,
        DRAWING_SCALE_BASIN,
        map_area_mm(),
        [*route_xy, *(point for basin in basins for point in basin["ring"])],
    )

    return {
        "route_xy": route_xy,
        "basins": basins,
        "contours": background["contours"],
        "streams": background["streams"],
    }