diff --git a/B04_PreProcess/B04_PreProcess_Engine_SheetParser.py b/B04_PreProcess/B04_PreProcess_Engine_SheetParser.py index 8d1ef96dd..59bf7027c 100644 --- a/B04_PreProcess/B04_PreProcess_Engine_SheetParser.py +++ b/B04_PreProcess/B04_PreProcess_Engine_SheetParser.py @@ -33,6 +33,7 @@ SHEET_LAYERS: dict[str, tuple[str, str]] = { "유수방향": ("N3P", "E0042326"), "기준점": ("N3P", "H0020000"), "지명": ("N3P", "H0040000"), + "도로중심선": ("N3L", "A0020000"), } # 산출 파일 접두사 (기존 processed/의 국가 GIS geojson과 구분) diff --git a/D03_PlanAnalysis/D03_PlanAnalysis_Layers.py b/D03_PlanAnalysis/D03_PlanAnalysis_Layers.py new file mode 100644 index 000000000..2e0b531a1 --- /dev/null +++ b/D03_PlanAnalysis/D03_PlanAnalysis_Layers.py @@ -0,0 +1,259 @@ +"""D 노선분석 지도 겹 — 도엽 겹 · 지적도 · 급경사 그림 (요구 395-2 · PLAN 569-1). + +- 도엽 겹(하천중심선 · 도로중심선 …) = 지형 받기 병합본 · 없으면 지난 병합 도엽에서 그 겹만 풀어 냄 +- 지적도 = 브이월드 WFS LX맵(`lt_c_landinfobasemap`) · B04 `request_vworld_wfs` 를 1 km 조각으로 · + 한 번 1,000 필지 한도 → 꽉 차면 넷으로 다시 쪼갬 · `pnu` 로 겹친 필지 거름 · + 저장 = `processed/연속지적도_bounds.geojson`(B04 `/geojson?layer=지적도` 자리) + 곁 열쇠 +- 급경사 = 평활 안 한 사업 DEM(`models_plan/dtm_sheet_<방법>.npz`) · E07 `slope_rad` · + 15° · 35° · 40° 단계 PNG · 파생 저장본(열쇠 = 격자 지문) + +계산 자리 = ③ 서버 단독. +""" + +from __future__ import annotations + +import base64 +import io +import json +import logging +import math +import tempfile +from concurrent.futures import ThreadPoolExecutor +from pathlib import Path +from typing import Any + +import numpy as np + +from B04_PreProcess.B04_PreProcess_Engine_GisVector import request_vworld_wfs +from B04_PreProcess.B04_PreProcess_Engine_SheetParser import SHEET_LAYERS, parse_and_merge_sheets +from B04_PreProcess.B04_PreProcess_Engine_SheetStore import get_project_map_sheets_dir +from common_util.common_util_derived_store import ( + derived_key, + file_stamp, + read_derived, + write_derived, +) +from common_util.common_util_json import atomic_write_json +from config.config_system_terrain import SHEET_SURFACE_DEFAULT_METHOD + +logger = logging.getLogger(__name__) + +Box = tuple[float, float, float, float] # x0, y0, x1, y1(작업 좌표계 m) + +CADASTRAL_TYPENAME = "lt_c_landinfobasemap" +CADASTRAL_FILE = "연속지적도_bounds.geojson" +CADASTRAL_KEY = "연속지적도_bounds.key.json" +TILE_M = 1000.0 +WFS_LIMIT = 1000 +MAX_SPLIT = 3 # 1 km → 125 m 까지 쪼갬 +WFS_WORKERS = 4 + +SLOPE_KIND = "plan_slope" +SLOPE_VERSION = 1 +# (아래 °, 위 °, 색 RGBA) — 임도망_편성 완 · 중 · 급경사지 · 타당성평가 35° · 노선_배치 40° 연속 +SLOPE_BANDS = ( + (15.0, 35.0, (250, 204, 21, 110)), + (35.0, 40.0, (249, 115, 22, 140)), + (40.0, 90.0, (220, 38, 38, 160)), +) +MAX_SIDE_PX = 4096 # ponytail: 한 변 넘으면 솎음 · 넓은 사업은 조각 그림으로 + + +def _read_json(path: Path) -> Any: + try: + return json.loads(path.read_text(encoding="utf-8")) + except (OSError, ValueError): + return None + + +# ── 도엽 겹 ───────────────────────────────────────────── + + +def ensure_sheet_layer(project_root: Path, key: str) -> Path: + """`processed/도엽_{key}.geojson` — 없으면 지난 병합 도엽에서 그 겹만 풀어 냄(0 건 = 빈 것).""" + if key not in SHEET_LAYERS: + raise KeyError(key) + processed = project_root / "B04_PreProcess" / "processed" + path = processed / f"도엽_{key}.geojson" + if path.is_file(): + return path + sheets = (_read_json(processed / "sheet_merge.key.json") or {}).get("sheets") or [] + if not sheets: + raise FileNotFoundError("병합한 도엽이 없습니다 — 지형 받기 뒤 보입니다.") + parse_and_merge_sheets(get_project_map_sheets_dir(project_root), sheets, processed, [key]) + if not path.is_file(): + atomic_write_json(path, {"type": "FeatureCollection", "features": []}) + return path + + +# ── 지적도 ────────────────────────────────────────────── + + +def tiles(box: Box, size: float = TILE_M) -> list[Box]: + """상자를 `size` 칸으로 쪼갬(끝 칸은 남은 만큼).""" + x0, y0, x1, y1 = box + nx, ny = max(1, math.ceil((x1 - x0) / size)), max(1, math.ceil((y1 - y0) / size)) + return [ + (x0 + i * size, y0 + j * size, min(x1, x0 + (i + 1) * size), min(y1, y0 + (j + 1) * size)) + for j in range(ny) + for i in range(nx) + ] + + +def quarter(box: Box) -> list[Box]: + x0, y0, x1, y1 = box + xm, ym = (x0 + x1) / 2, (y0 + y1) / 2 + return [(x0, y0, xm, ym), (xm, y0, x1, ym), (x0, ym, xm, y1), (xm, ym, x1, y1)] + + +def _lonlat_bounds(box: Box, epsg: int) -> dict[str, float]: + from common_util.common_util_crs import tm_to_lonlat + + x0, y0, x1, y1 = box + ring = [tm_to_lonlat(x, y, epsg) for x, y in ((x0, y0), (x1, y0), (x1, y1), (x0, y1))] + lons, lats = [p[0] for p in ring], [p[1] for p in ring] + return {"min_lon": min(lons), "min_lat": min(lats), "max_lon": max(lons), "max_lat": max(lats)} + + +def _wfs(box: Box, epsg: int) -> list[dict[str, Any]] | None: + """조각 하나 — B04 WFS 받기(실패 = 파일 없음 · 키 값은 남기지 않음). None = 실패.""" + with tempfile.TemporaryDirectory() as folder: + request_vworld_wfs( + CADASTRAL_TYPENAME, "tile.geojson", _lonlat_bounds(box, epsg), Path(folder) + ) + data = _read_json(Path(folder) / "tile.geojson") + if not isinstance(data, dict) or not isinstance(data.get("features"), list): + return None + return data["features"] + + +def fetch_tile(box: Box, epsg: int, depth: int = 0) -> tuple[list[dict[str, Any]], int]: + """(필지, 실패 조각 수) — 1,000 꽉 차면 넷으로 다시 쪼갬.""" + features = _wfs(box, epsg) + if features is None: + return [], 1 + if len(features) < WFS_LIMIT or depth >= MAX_SPLIT: + return features, 0 + out: list[dict[str, Any]] = [] + failed = 0 + for part in quarter(box): + got, bad = fetch_tile(part, epsg, depth + 1) + out.extend(got) + failed += bad + return out, failed + + +def dedup_parcels(features: list[dict[str, Any]]) -> list[dict[str, Any]]: + """조각 경계에 걸친 필지 — `pnu`(없으면 지물 id) 첫 것만.""" + seen: set[str] = set() + out = [] + for feature in features: + props = feature.get("properties") or {} + key = str(props.get("pnu") or feature.get("id") or json.dumps(feature.get("geometry"))) + if key not in seen: + seen.add(key) + out.append(feature) + return out + + +def _covers(have: Any, box: Box) -> bool: + try: + return have[0] <= box[0] and have[1] <= box[1] and have[2] >= box[2] and have[3] >= box[3] + except (TypeError, IndexError): + return False + + +def ensure_cadastral(processed: Path, box: Box, epsg: int) -> dict[str, Any]: + """지적도 — 열쇠 상자가 덮고 실패 조각 없으면 그대로 · 아니면 상자 통째로 다시 받음.""" + path, key_path = processed / CADASTRAL_FILE, processed / CADASTRAL_KEY + key = _read_json(key_path) or {} + if path.is_file() and not key.get("failed") and _covers(key.get("box"), box): + return {"path": path, "fetched": False, "count": key.get("count"), "failed": 0} + parts = tiles(box) + with ThreadPoolExecutor(max_workers=WFS_WORKERS) as pool: + results = list(pool.map(lambda part: fetch_tile(part, epsg), parts)) + failed = sum(bad for _, bad in results) + if failed == len(parts): + raise RuntimeError(f"지적도를 받지 못했습니다 — 브이월드 WFS 조각 {failed}개 실패") + features = dedup_parcels([f for got, _ in results for f in got]) + processed.mkdir(parents=True, exist_ok=True) + atomic_write_json(path, {"type": "FeatureCollection", "features": features}) + atomic_write_json( + key_path, + { + "box": [round(v, 3) for v in box], + "epsg": int(epsg), + "count": len(features), + "failed": failed, + }, + ) + logger.info("D 지적도 받기: 조각 %d · 실패 %d · 필지 %d", len(parts), failed, len(features)) + return {"path": path, "fetched": True, "count": len(features), "failed": failed} + + +# ── 급경사 ────────────────────────────────────────────── + + +def slope_bands_png(z: np.ndarray, valid: np.ndarray, cell_m: float) -> dict[str, Any]: + """격자(행 = 남 → 북) → 단계 PNG(위 = 북 · 15° 미만 · 빈 칸 투명) + 단계 칸 비율.""" + from PIL import Image + + from E07_StabArea.E07_StabArea_Engine import slope_rad + + stride = max(1, math.ceil(max(z.shape) / MAX_SIDE_PX)) + grid = np.where(valid, np.asarray(z, float), np.nan)[::stride, ::stride] + deg = np.degrees(slope_rad(grid, cell_m * stride)) + known = np.isfinite(deg) + index = np.zeros(deg.shape, np.uint8) + for i, (low, high, _) in enumerate(SLOPE_BANDS, start=1): + index[known & (deg >= low) & (deg < high)] = i + total = max(1, int(known.sum())) + image = Image.fromarray(np.flipud(index), mode="P") + palette = [0, 0, 0] + [c for _, _, rgba in SLOPE_BANDS for c in rgba[:3]] + image.putpalette(palette + [0] * (768 - len(palette))) + alpha = bytes([0] + [rgba[3] for _, _, rgba in SLOPE_BANDS]) + buffer = io.BytesIO() + image.save(buffer, format="PNG", optimize=True, transparency=alpha) + return { + "width": image.width, + "height": image.height, + "image": "data:image/png;base64," + base64.b64encode(buffer.getvalue()).decode("ascii"), + "bands": [ + { + "min_deg": low, + "max_deg": high, + "color": "#{:02x}{:02x}{:02x}".format(*rgba[:3]), + "ratio": round(int((index == i).sum()) / total, 4), + } + for i, (low, high, rgba) in enumerate(SLOPE_BANDS, start=1) + ], + } + + +def slope_layer(project_root: Path, models_dir: Path, epsg: int) -> dict[str, Any]: + """급경사 겹 — 지문(격자 파일 · 단계 · 판) 같으면 저장본 그대로.""" + from common_util.common_util_crs import tm_to_lonlat + + npz = models_dir / f"dtm_sheet_{SHEET_SURFACE_DEFAULT_METHOD}.npz" + stamp = file_stamp(npz) + if stamp is None: + raise FileNotFoundError("사업 지형면이 없습니다 — 지형 받기 뒤 보입니다.") + key = derived_key( + SLOPE_KIND, SLOPE_VERSION, npz.name, stamp, int(epsg), SLOPE_BANDS, MAX_SIDE_PX + ) + saved = read_derived(project_root, SLOPE_KIND, key) + if saved is not None: + return saved["value"] + with np.load(npz) as dem: + x, y = np.asarray(dem["x"], float), np.asarray(dem["y"], float) + cell = float(dem["resolution"][0]) if "resolution" in dem.files else float(x[1] - x[0]) + value = slope_bands_png(dem["z"], np.asarray(dem["valid_mask"], bool), cell) + half = cell / 2 + box = (x[0] - half, y[0] - half, x[-1] + half, y[-1] + half) + value["box_m"] = [float(v) for v in box] + value["box"] = [ + list(tm_to_lonlat(px, py, epsg)) + for px, py in ((box[0], box[1]), (box[2], box[1]), (box[2], box[3]), (box[0], box[3])) + ] + write_derived(project_root, SLOPE_KIND, key, value) + return value diff --git a/D03_PlanAnalysis/D03_PlanAnalysis_Router.py b/D03_PlanAnalysis/D03_PlanAnalysis_Router.py index c4a4be961..1de51d747 100644 --- a/D03_PlanAnalysis/D03_PlanAnalysis_Router.py +++ b/D03_PlanAnalysis/D03_PlanAnalysis_Router.py @@ -3,6 +3,7 @@ - `GET /api/projects/{id}/plan/terrain/extent` — 받기 전 범위(상자 · 도엽 · 칸 수 · 받은 범위) - `POST …/plan/terrain` — 사업 전체 지형 한 벌 받기(뒤에서 · 누를 때만) · 끝나면 `surface_models` 행 - `GET …/plan/progress` — 일의 진행률 +- `GET …/plan/layers/sheet?key=` · `…/cadastral` · `…/raster?kind=slope` — 지도 겹(PLAN 569-1) 범위 = 전체 시점 · 종점(D02 정본 `load_base`) + 그은 점(`extra`) 상자 + 여유. 노선 정본(손 그은 선)은 @@ -18,13 +19,15 @@ from pathlib import Path from typing import Any from uuid import UUID -from fastapi import APIRouter, HTTPException +from fastapi import APIRouter, HTTPException, Request, Response from pydantic import BaseModel from B03_FileInput.B03_FileInput_Repository import get_project_storage_relative_path +from common_util.common_util_http_cache import cached_file_response from common_util.common_util_storage import resolve_stored_project_path from config.config_db import get_db_pool from D02_PlanBase import D02_PlanBase_Router as base +from D03_PlanAnalysis import D03_PlanAnalysis_Layers as layers from D03_PlanAnalysis import D03_PlanAnalysis_Terrain as terrain logger = logging.getLogger(__name__) @@ -153,3 +156,53 @@ async def plan_progress(project_id: UUID) -> dict[str, Any]: "status": "success", **{kind: dict(_JOBS.get((str(project_id), kind)) or idle) for kind in KINDS}, } + + +async def _layer_root(project_id: UUID) -> tuple[Path, int, terrain.Box]: + """(실경로, 좌표계, 받은 범위) — 지형 받기 전이면 404.""" + root, epsg, _ = await _plan_ends(project_id) + box = terrain.received_box(root) + if box is None: + raise HTTPException(status_code=404, detail="지형 받기 뒤 보입니다.") + return root, epsg, box + + +@router.get("/{project_id}/plan/layers/sheet", response_model=None) +async def plan_layer_sheet(project_id: UUID, key: str, request: Request) -> Response: + """도엽 겹 GeoJSON(하천중심선 · 도로중심선 …) — 없으면 지난 병합 도엽에서 풀어 냄.""" + from B04_PreProcess.B04_PreProcess_Router_GIS import _display_copy + + root, _, _ = await _layer_root(project_id) + try: + path = await asyncio.to_thread(layers.ensure_sheet_layer, root, key) + except KeyError as error: + raise HTTPException(status_code=400, detail=f"없는 도엽 겹: {key}") from error + except FileNotFoundError as error: + raise HTTPException(status_code=404, detail=str(error)) from error + display = await asyncio.to_thread(_display_copy, path) + return cached_file_response(request, display, "application/geo+json") + + +@router.get("/{project_id}/plan/layers/cadastral", response_model=None) +async def plan_layer_cadastral(project_id: UUID, request: Request) -> Response: + """지적도 GeoJSON(브이월드 LX맵 · 참고용) — 받은 범위를 덮지 않으면 지금 받음.""" + root, epsg, box = await _layer_root(project_id) + processed = root / terrain.STAGE / "processed" + try: + got = await asyncio.to_thread(layers.ensure_cadastral, processed, box, epsg) + except RuntimeError as error: + raise HTTPException(status_code=502, detail=str(error)) from error + return cached_file_response(request, got["path"], "application/geo+json") + + +@router.get("/{project_id}/plan/layers/raster", response_model=None) +async def plan_layer_raster(project_id: UUID, kind: str = "slope") -> dict[str, Any]: + """급경사 단계 그림(15° · 35° · 40°) — 격자 지문 같으면 저장본.""" + if kind != "slope": + raise HTTPException(status_code=400, detail=f"없는 그림 겹: {kind}") + root, epsg, _ = await _layer_root(project_id) + try: + value = await asyncio.to_thread(layers.slope_layer, root, terrain.models_dir_of(root), epsg) + except FileNotFoundError as error: + raise HTTPException(status_code=404, detail=str(error)) from error + return {"status": "success", "kind": kind, **value} diff --git a/D03_PlanAnalysis/D03_PlanAnalysis_Terrain.py b/D03_PlanAnalysis/D03_PlanAnalysis_Terrain.py index 9e909f37b..d59576bf3 100644 --- a/D03_PlanAnalysis/D03_PlanAnalysis_Terrain.py +++ b/D03_PlanAnalysis/D03_PlanAnalysis_Terrain.py @@ -6,6 +6,7 @@ - 도엽 = 범위 상자에 닿는 도엽 모두 → `ensure_sheets`(없는 것만 · prune 없음) - 병합 = 지난번 병합 도엽 ∪ 이번 도엽(옛 도엽 유지 · 열쇠 같으면 건너뜀) - 위성 = 범위 도엽 도곽 ∪ 옛 위성 범위(덮으면 건너뜀) +- 지적도 = 범위 상자 · 브이월드 WFS 1 km 조각(`D03_PlanAnalysis_Layers` · 실패해도 계속) - DEM = 한 벌 `models_plan/`(연도는 D03 나눔 뒤에 정해짐) · 받은 범위 = `plan_terrain.json` - 작업 좌표계 = 지점 epsg · `B03_FileInput/input/prj/plan.prj` 한 장(B04 · B05 좌표계 창구가 읽음) - 받기 전 범위 보기 = `plan_extent`(내려받기 없음) · 받기는 사용자가 누를 때만 @@ -42,6 +43,7 @@ from config.config_system_terrain import ( SHEET_SURFACE_DEFAULT_METHOD, SHEET_SURFACE_MARGIN_M, ) +from D03_PlanAnalysis import D03_PlanAnalysis_Layers as layers logger = logging.getLogger(__name__) @@ -254,6 +256,12 @@ def ensure_plan_terrain( merged = _merge_sheets(processed, store, list(got["available"])) report(60, "위성 지도 받는 중") satellite = _satellite(processed, prj, [(p[1], p[0]) for p in box_ring(box, epsg)], epsg) + report(62, "지적도 받는 중") + try: + cadastral = layers.ensure_cadastral(processed, box, epsg)["count"] + except Exception as exc: # noqa: BLE001 — 지적도는 겹 · 지형 받기를 막지 않음(켤 때 다시) + logger.warning("D 지적도 받기 실패: %s", exc) + cadastral = None report(65, "사업 지형면 만드는 중") models_dir = models_dir_of(project_root) @@ -287,6 +295,7 @@ def ensure_plan_terrain( "merged": merged, }, "satellite": satellite, + "cadastral": cadastral, "box_m": list(box), "models": models, } diff --git a/resources/tester/test_530_1_plan_terrain.py b/resources/tester/test_530_1_plan_terrain.py index a019f61ec..082c77c00 100644 --- a/resources/tester/test_530_1_plan_terrain.py +++ b/resources/tester/test_530_1_plan_terrain.py @@ -54,6 +54,7 @@ def fake(tmp_path, monkeypatch): monkeypatch.setattr(terrain, "parse_and_merge_sheets", merge) monkeypatch.setattr(terrain, "download_vworld_satellite_map", satellite) monkeypatch.setattr(terrain, "build_sheet_surface_model", dem) + monkeypatch.setattr(terrain.layers, "ensure_cadastral", lambda *a: {"count": 0}) # 569-1 return tmp_path, seen diff --git a/resources/tester/test_569_1_plan_layers.py b/resources/tester/test_569_1_plan_layers.py new file mode 100644 index 000000000..d6a057a32 --- /dev/null +++ b/resources/tester/test_569_1_plan_layers.py @@ -0,0 +1,181 @@ +"""D 노선분석 지도 겹(PLAN 569-1) — 도엽 겹 풀어 냄 · 지적도 1 km 조각 · 1,000 꽉 차면 넷 · +겹친 필지 거름 · 급경사 15° · 35° · 40° 단계 · 격자 지문 같으면 다시 안 셈. + +./venv/Scripts/python.exe -m pytest resources/tester/test_569_1_plan_layers.py -q +""" + +import base64 +import io +import itertools +import json +import os +from pathlib import Path + +import numpy as np +import pytest +from PIL import Image + +from B04_PreProcess.B04_PreProcess_Engine_SheetParser import SHEET_LAYERS +from D03_PlanAnalysis import D03_PlanAnalysis_Layers as layers + +EPSG = 5186 +BOX = (230000.0, 450000.0, 234700.0, 452500.0) # 4.7 × 2.5 km + + +# ── 도엽 겹 ── + + +def test_road_layer_listed() -> None: + assert SHEET_LAYERS["도로중심선"] == ("N3L", "A0020000") + assert SHEET_LAYERS["하천중심선"] == ("N3L", "E0020000") + + +def test_sheet_layer_unpacks_once(tmp_path, monkeypatch) -> None: + processed = tmp_path / "B04_PreProcess" / "processed" + processed.mkdir(parents=True) + calls = [] + + def merge(store, sheets, out, keys): + calls.append((list(sheets), keys)) # 지물 0 건 = 파일 안 씀 + + monkeypatch.setattr(layers, "parse_and_merge_sheets", merge) + with pytest.raises(FileNotFoundError): + layers.ensure_sheet_layer(tmp_path, "도로중심선") + with pytest.raises(KeyError): + layers.ensure_sheet_layer(tmp_path, "없는겹") + (processed / "sheet_merge.key.json").write_text( + json.dumps({"sheets": ["36712088", "36712089"], "parser": 2}), "utf-8" + ) + path = layers.ensure_sheet_layer(tmp_path, "도로중심선") + assert calls == [(["36712088", "36712089"], ["도로중심선"])] + assert json.loads(path.read_text("utf-8"))["features"] == [] + layers.ensure_sheet_layer(tmp_path, "도로중심선") + assert len(calls) == 1 # 있으면 그대로 + + +# ── 지적도 ── + + +def test_tiles_cover_box() -> None: + parts = layers.tiles(BOX) + assert len(parts) == 15 + assert sum((b[2] - b[0]) * (b[3] - b[1]) for b in parts) == pytest.approx(4700 * 2500) + assert max(b[2] for b in parts) == BOX[2] and max(b[3] for b in parts) == BOX[3] + assert len(layers.quarter(parts[0])) == 4 + + +class FakeWfs: + """조각 너비가 `full_deg` 보다 넓으면 1,000 필지(꽉 참) · 아니면 셋(하나는 이웃과 같은 pnu) · + 처음 `fail` 번은 실패(B04 꼴 — 파일 없음).""" + + def __init__(self, full_deg: float = 1.0, fail: int = 0) -> None: + self.full_deg, self.fail, self._count = full_deg, fail, itertools.count(1) + self.calls = 0 + + def __call__(self, typename, filename, bounds, output_dir) -> None: + assert typename == "lt_c_landinfobasemap" + self.calls = next(self._count) # 조각을 여러 줄기로 받음 — count 는 한 번에 하나 + if self.calls <= self.fail: + return + width = bounds["max_lon"] - bounds["min_lon"] + tag = f"{bounds['min_lon']:.5f},{bounds['min_lat']:.5f}" + n = 1000 if width > self.full_deg else 3 + feats = [{"id": f"{tag}/{i}", "properties": {"pnu": f"{tag}/{i}"}} for i in range(n)] + feats[-1]["properties"]["pnu"] = "SHARED" + body = {"type": "FeatureCollection", "features": feats} + (Path(output_dir) / filename).write_text(json.dumps(body), "utf-8") + + +def test_full_tile_splits_into_four(monkeypatch) -> None: + fake = FakeWfs(full_deg=0.008) # 1 km 칸 ≈ 0.011° 꽉 참 · 500 m 칸 = 셋 + monkeypatch.setattr(layers, "request_vworld_wfs", fake) + got, failed = layers.fetch_tile((230000.0, 450000.0, 231000.0, 451000.0), EPSG) + assert fake.calls == 5 and failed == 0 and len(got) == 12 + + +def test_dedup_by_pnu() -> None: + feats = [{"properties": {"pnu": "1"}}, {"properties": {"pnu": "1"}}, {"id": "x"}] + assert len(layers.dedup_parcels(feats)) == 2 + + +def test_cadastral_fetch_once_then_grow(tmp_path, monkeypatch) -> None: + fake = FakeWfs() + monkeypatch.setattr(layers, "request_vworld_wfs", fake) + got = layers.ensure_cadastral(tmp_path, BOX, EPSG) + assert got["fetched"] and fake.calls == 15 and got["failed"] == 0 + assert got["count"] == 15 * 2 + 1 # 조각마다 둘 + 같은 pnu 하나 + saved = json.loads((tmp_path / layers.CADASTRAL_FILE).read_text("utf-8")) + assert len(saved["features"]) == got["count"] + again = layers.ensure_cadastral(tmp_path, (231000.0, 450500.0, 233000.0, 452000.0), EPSG) + assert not again["fetched"] and fake.calls == 15 # 덮으면 다시 안 받음 + layers.ensure_cadastral(tmp_path, (*BOX[:2], BOX[2] + 500, BOX[3]), EPSG) + assert fake.calls == 15 + 18 # 넓어지면 상자 통째로 다시(5.2 km = 6 × 3 조각) + + +def test_cadastral_failures(tmp_path, monkeypatch) -> None: + monkeypatch.setattr(layers, "request_vworld_wfs", FakeWfs(fail=99)) + with pytest.raises(RuntimeError): + layers.ensure_cadastral(tmp_path, BOX, EPSG) + assert not (tmp_path / layers.CADASTRAL_FILE).exists() + monkeypatch.setattr(layers, "request_vworld_wfs", FakeWfs(fail=1)) + first = layers.ensure_cadastral(tmp_path, BOX, EPSG) + assert first["failed"] == 1 and (tmp_path / layers.CADASTRAL_FILE).exists() + again = FakeWfs() + monkeypatch.setattr(layers, "request_vworld_wfs", again) + assert layers.ensure_cadastral(tmp_path, BOX, EPSG)["fetched"] and again.calls == 15 + + +# ── 급경사 ── + + +def _plane(deg_south: float, deg_north: float, n: int = 40, cell: float = 2.0) -> np.ndarray: + """남쪽 반 = deg_south · 북쪽 반 = deg_north 경사(동쪽으로 오름).""" + x = np.arange(n) * cell + z = np.empty((n, n)) + z[: n // 2] = x * np.tan(np.radians(deg_south)) + z[n // 2 :] = x * np.tan(np.radians(deg_north)) + return z + + +def _indices(value) -> np.ndarray: + raw = base64.b64decode(value["image"].split(",", 1)[1]) + return np.asarray(Image.open(io.BytesIO(raw))) + + +def test_slope_bands_and_north_up() -> None: + z = _plane(10.0, 37.0) + value = layers.slope_bands_png(z, np.ones(z.shape, bool), 2.0) + ratios = [b["ratio"] for b in value["bands"]] + assert [b["min_deg"] for b in value["bands"]] == [15.0, 35.0, 40.0] + assert ratios[0] < 0.1 and 0.4 < ratios[1] < 0.6 and ratios[2] < 0.1 # 경계 줄만 + index = _indices(value) + assert value["width"] == value["height"] == 40 + assert (index[5, 5:35] == 2).all() and (index[-5, 5:35] == 0).all() # 위 = 북 = 37° + steep = layers.slope_bands_png(_plane(45.0, 20.0), np.ones((40, 40), bool), 2.0) + assert steep["bands"][2]["ratio"] > 0.4 and steep["bands"][0]["ratio"] > 0.4 + hole = np.ones(z.shape, bool) + hole[25:30, 10:20] = False + assert (_indices(layers.slope_bands_png(z, hole, 2.0))[-29:-24, 10:20] == 0).all() + + +def test_slope_layer_fingerprint(tmp_path, monkeypatch) -> None: + models = tmp_path / "B04_PreProcess" / "models_plan" + models.mkdir(parents=True) + npz = models / "dtm_sheet_laplace.npz" + z = _plane(37.0, 37.0) + x = 200000.0 + np.arange(40) * 2.0 + y = 450000.0 + np.arange(40) * 2.0 + np.savez(npz, x=x, y=y, z=z, valid_mask=np.ones(z.shape, bool), resolution=[2.0]) + counted = [] + real = layers.slope_bands_png + monkeypatch.setattr(layers, "slope_bands_png", lambda *a: counted.append(1) or real(*a)) + first = layers.slope_layer(tmp_path, models, EPSG) + assert first["box_m"] == [199999.0, 449999.0, 200079.0, 450079.0] and len(first["box"]) == 4 + assert layers.slope_layer(tmp_path, models, EPSG) == first and len(counted) == 1 + np.savez(npz, x=x, y=y, z=z * 0.5, valid_mask=np.ones(z.shape, bool), resolution=[2.0]) + os.utime(npz, ns=(1, 1)) + layers.slope_layer(tmp_path, models, EPSG) + assert len(counted) == 2 # 격자가 바뀌면 다시 셈 + (npz).unlink() + with pytest.raises(FileNotFoundError): + layers.slope_layer(tmp_path, models, EPSG)