"""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_Landuse import LANDUSE_LABEL
from B07_DesignDetail.B07_DesignDetail_Engine_Cad_Lidar import LIDAR_LABEL, hillshade_png
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·MultiPolygon을 점열 목록으로 편다."""
    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)
        ]
    # 연속지적도·행정구역은 MultiPolygon 이다 — 폴리곤마다 고리를 모두 편다(2026-09-04).
    if kind == "MultiPolygon":
        return [
            [(float(p[0]), float(p[1])) for p in ring if isinstance(p, list)]
            for polygon in coordinates
            if isinstance(polygon, list)
            for ring in polygon
            if isinstance(ring, 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)]


LANDUSE_ID = re.compile(r"^landuse(?:_(\d+))?$")

# B04 가 내려받아 저장하는 지적·행정구역 GeoJSON (전부 WGS84).
PARCEL_FILE = "연속지적도_bounds.geojson"
EMD_FILE = "행정구역_읍면동_bounds.geojson"
SGG_FILE = "행정구역_시군구_bounds.geojson"


def _clip_rings(
    path: Path, crs: str, box: tuple[float, float, float, float]
) -> list[list[tuple[float, float]]]:
    """행정구역 경계를 사업지 좌표계로 돌려 도곽 범위로 절취한다(속성은 안 씀)."""
    rings: list[list[tuple[float, float]]] = []
    for ring in _metric_lines(path, crs):
        rings.extend(clip_line_to_box(ring, box))
    return rings


def _contains(ring: list[tuple[float, float]], point: tuple[float, float]) -> bool:
    """점이 고리 안에 드는지 (반직선 교차 판정)."""
    x, y = point
    inside = False
    for index in range(len(ring)):
        x1, y1 = ring[index - 1]
        x2, y2 = ring[index]
        if (y1 > y) != (y2 > y) and x < (x2 - x1) * (y - y1) / ((y2 - y1) or 1e-12) + x1:
            inside = not inside
    return inside


def _clip_parcels(
    path: Path, crs: str, box: tuple[float, float, float, float]
) -> list[dict[str, Any]]:
    """연속지적도를 사업지 좌표계로 돌려 도곽 안 필지만 남긴다 (지번 표기용 속성 포함).

    필지는 지번을 적어야 하므로 경계선만 자르는 `_metric_lines` 캐시를 쓰지 못한다 —
    피처와 속성을 짝지어 읽는다. 도곽 밖 필지는 여기서 버려 도면이 무거워지지 않게 한다.

    두 가지를 함께 낸다.
      - `ring` : 도곽으로 자른 경계선(밖으로 나가는 부분은 버린다). 자르지 않으면
        산지 대필지 하나가 도면을 10 km 로 벌린다(2026-09-04 실측: 콘텐츠 8,368 mm).
      - `label_at` : 도곽을 **통째로 감싸는** 필지의 지번 자리. 임야 대필지 안에 노선이
        들어앉으면 경계선이 도곽 안에 하나도 없어 지번이 사라진다(2026-09-04 실측:
        용화_LAS 노선이 「산77-1임 일월면 용화리」 한 필지 안에 통째로 들어감).
    """
    if not path.is_file():
        return []
    transformer = Transformer.from_crs("EPSG:4326", crs, always_xy=True)
    min_x, min_y, max_x, max_y = box
    center = ((min_x + max_x) / 2.0, (min_y + max_y) / 2.0)
    parcels: list[dict[str, Any]] = []
    for feature in _geojson_features(path):
        properties = feature.get("properties") or {}
        for ring in _geometry_lines(feature.get("geometry")):
            converted = [
                (float(x), float(y))
                for x, y in (transformer.transform(point[0], point[1]) for point in ring)
            ]
            if len(converted) < 3:
                continue
            if max(x for x, _y in converted) < min_x or min(x for x, _y in converted) > max_x:
                continue
            if max(y for _x, y in converted) < min_y or min(y for _x, y in converted) > max_y:
                continue
            parts = [part for part in clip_line_to_box(converted, box) if len(part) >= 2]
            for part in parts:
                parcels.append({"ring": part, "props": properties})
            if not parts and _contains(converted, center):
                parcels.append({"ring": [], "props": properties, "label_at": center})
    return parcels


def landuse_source(context: Any, longitudinal: dict[str, Any], drawing_id: str) -> dict[str, Any]:
    """용지도 한 장의 입력(노선·등고선·연속지적도·행정구역)을 사업지 CRS(m)로 모은다.

    축척·도곽·장 나눔은 계획평면도와 같다 — 배경도 같은 창구(`map_background`)를 쓴다.
    """
    match = LANDUSE_ID.fullmatch(drawing_id)
    if not match:
        raise ValueError("올바르지 않은 용지도 ID입니다.")
    chunks = plan_chunks(plan_stations(longitudinal))
    number = int(match.group(1)) if match.group(1) else 1
    chunk = next((item for item in chunks if item["number"] == number), None)
    if chunk is None:
        raise FileNotFoundError("요청한 용지도 장을 찾을 수 없습니다.")
    total = len(chunks)
    start_m, end_m = float(chunk["start_m"]), float(chunk["end_m"])

    route_xy = [
        (vertex.x, vertex.y)
        for vertex in context.vertices
        if total <= 1 or start_m <= float(getattr(vertex, "chainage_m", 0.0) or 0.0) <= end_m
    ]
    background = map_background(
        Path(context.project_root),
        context.crs,
        DRAWING_SCALE_PLAN,
        plan_area_mm(),
        route_xy,
    )
    # 지적·행정 경계는 등고선과 **같은 범위**로 자른다 — 배경보다 넓으면 도면이 A1을 넘는다.
    area_w_mm, area_h_mm = plan_area_mm()
    half_w = area_w_mm / 2.0 * DRAWING_SCALE_PLAN / 1000.0
    half_h = area_h_mm / 2.0 * DRAWING_SCALE_PLAN / 1000.0
    center_x = (min(x for x, _y in route_xy) + max(x for x, _y in route_xy)) / 2.0
    center_y = (min(y for _x, y in route_xy) + max(y for _x, y in route_xy)) / 2.0
    box = (
        min(center_x - half_w, min(x for x, _y in route_xy)),
        min(center_y - half_h, min(y for _x, y in route_xy)),
        max(center_x + half_w, max(x for x, _y in route_xy)),
        max(center_y + half_h, max(y for _x, y in route_xy)),
    )
    sheet_dir = Path(context.project_root) / "B04_PreProcess" / "processed"
    label = LANDUSE_LABEL if total <= 1 else f"{LANDUSE_LABEL} {number}장"
    return {
        "label": label,
        "route_xy": route_xy,
        "contours": background["contours"],
        "parcels": _clip_parcels(sheet_dir / PARCEL_FILE, context.crs, box),
        "emd_rings": _clip_rings(sheet_dir / EMD_FILE, context.crs, box),
        "sgg_rings": _clip_rings(sheet_dir / SGG_FILE, context.crs, box),
    }


LIDAR_ID = re.compile(r"^plan_lidar(?:_(\d+))?$")


def _sheet_box(
    route_xy: list[tuple[float, float]],
) -> tuple[float, float, float, float]:
    """그 장의 도곽 범위(실좌표 m) — 계획평면도·용지도·라이다가 같은 규칙을 쓴다."""
    area_w_mm, area_h_mm = plan_area_mm()
    half_w = area_w_mm / 2.0 * DRAWING_SCALE_PLAN / 1000.0
    half_h = area_h_mm / 2.0 * DRAWING_SCALE_PLAN / 1000.0
    center_x = (min(x for x, _y in route_xy) + max(x for x, _y in route_xy)) / 2.0
    center_y = (min(y for _x, y in route_xy) + max(y for _x, y in route_xy)) / 2.0
    return (
        min(center_x - half_w, min(x for x, _y in route_xy)),
        min(center_y - half_h, min(y for _x, y in route_xy)),
        max(center_x + half_w, max(x for x, _y in route_xy)),
        max(center_y + half_h, max(y for _x, y in route_xy)),
    )


def _chunk_route(
    context: Any, longitudinal: dict[str, Any], number: int
) -> tuple[list[tuple[float, float]], dict[str, Any], int]:
    """장 번호로 그 장의 노선 구간을 잘라 낸다 (계획평면도 장 나눔과 같은 기준)."""
    chunks = plan_chunks(plan_stations(longitudinal))
    chunk = next((item for item in chunks if item["number"] == number), None)
    if chunk is None:
        raise FileNotFoundError("요청한 장을 찾을 수 없습니다.")
    total = len(chunks)
    start_m, end_m = float(chunk["start_m"]), float(chunk["end_m"])
    route_xy = [
        (vertex.x, vertex.y)
        for vertex in context.vertices
        if total <= 1 or start_m <= float(getattr(vertex, "chainage_m", 0.0) or 0.0) <= end_m
    ]
    return route_xy, chunk, total


def lidar_source(context: Any, longitudinal: dict[str, Any], drawing_id: str) -> dict[str, Any]:
    """라이다 계획평면도 한 장의 입력(노선 + 지표면 음영기복 그림)을 모은다.

    지표면은 확정 DTM 격자(`dtm_{필터}[_smooth].npz`)를 도곽 범위로 잘라 쓴다 —
    점구름을 그대로 그리면 수천만 점이라 도면 만들기가 느려진다(2026-09-04 사용자 지시).
    """
    match = LIDAR_ID.fullmatch(drawing_id)
    if not match:
        raise ValueError("올바르지 않은 라이다 계획평면도 ID입니다.")
    number = int(match.group(1)) if match.group(1) else 1
    route_xy, chunk, total = _chunk_route(context, longitudinal, number)
    box = _sheet_box(route_xy)
    label = LIDAR_LABEL if total <= 1 else f"{LIDAR_LABEL} {chunk['number']}장"

    shade_image: str | None = None
    shade_box: tuple[float, float, float, float] | None = None
    try:
        shade_image, shade_box = _hillshade_for_box(context, box)
    except (FileNotFoundError, ValueError, OSError) as exc:
        # 지표면이 없어도 노선·도각은 그린다 — 빈 화면보다 낫다.
        logger.warning("B07 라이다 계획평면도: 음영기복을 만들지 못했습니다 — %s", exc)

    return {
        "label": label,
        "route_xy": route_xy,
        "shade_image": shade_image,
        "shade_box": shade_box,
    }


def _hillshade_for_box(
    context: Any, box: tuple[float, float, float, float]
) -> tuple[str, tuple[float, float, float, float]]:
    """확정 DTM 격자를 도곽 범위로 잘라 음영기복 PNG(data URL)와 실제 덮은 범위를 낸다."""
    import numpy as np

    params = getattr(context, "surface_params", None) or {}
    source_filter = str(params.get("source_filter") or "csf")
    smooth = bool(params.get("smooth", True))
    models_dir = Path(context.project_root) / "B04_PreProcess" / "models"
    candidates = [models_dir / f"dtm_{source_filter}_smooth.npz"] if smooth else []
    candidates.append(models_dir / f"dtm_{source_filter}.npz")
    candidates.extend(sorted(models_dir.glob("dtm_*_smooth.npz")))
    candidates.extend(sorted(models_dir.glob("dtm_*.npz")))
    path = next((item for item in candidates if item.is_file()), None)
    if path is None:
        raise FileNotFoundError("확정 지표면 격자(DTM)가 없습니다.")

    with np.load(path, allow_pickle=False) as data:
        grid_x = np.asarray(data["x"], dtype=np.float64)
        grid_y = np.asarray(data["y"], dtype=np.float64)
        grid_z = np.asarray(data["z"], dtype=np.float64)
        valid = np.asarray(data["valid_mask"], dtype=bool)
        resolution = float(np.asarray(data["resolution"]).reshape(-1)[0])

    min_x, min_y, max_x, max_y = box
    columns = np.where((grid_x >= min_x) & (grid_x <= max_x))[0]
    rows = np.where((grid_y >= min_y) & (grid_y <= max_y))[0]
    if columns.size < 2 or rows.size < 2:
        raise ValueError("도곽 안에 지표면 격자가 없습니다.")
    sliced_z = grid_z[rows[0] : rows[-1] + 1, columns[0] : columns[-1] + 1]
    sliced_valid = valid[rows[0] : rows[-1] + 1, columns[0] : columns[-1] + 1]

    data_url, _width, _height = hillshade_png(sliced_z, sliced_valid, resolution)
    return (
        data_url,
        (
            float(grid_x[columns[0]]),
            float(grid_y[rows[0]]),
            float(grid_x[columns[-1]]),
            float(grid_y[rows[-1]]),
        ),
    )


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"],
    }