Files
Aislo/D03_PlanAnalysis/D03_PlanAnalysis_Layers.py
eomsangdonandClaude Opus 5.5 8793b4f6b2 feat(D03): 노선분석 바탕 MapLibre 전체 적용 — 겹 스물 · 등고 글 · 지명 GL · 처음 = 켬 · WebGL 못 쓰면 옛 캔버스 (PLAN 586-2 ①②③⑤ · 요구 448)
- GL 겹 스물: 위성 조각 · 음영 · 급경사 · 산사태위험(image) · 산림입지 · 보호구역 · 띠 · 물 면 · 건물(fill + 테두리)
  · 지적 · 등고선 · 도로 · 임도 짐작 · 올린 임도망도 · 세류 · 단층 · 지질구조선 · 0차곡 · 산사태취약 필지(line · 점선 px 그대로)
  · 산사태취약 점 · 지명 점(circle) — 색 · 굵기 · 투명도 · 차례 = 옛 캔버스 그대로
- 음영 = 위성 위 「곱하기 0.75」와 같은 셈의 검은 반투명 그림(밝기 g → 투명도 0.75·(1−g)) · 위성 끔 = 그대로
- 등고 수치 · 지명 = symbol(겹치면 뺌 · 등고 글은 선 따라) · 글꼴 = 화면 Noto Sans KR(OFL · 자체 보관)을 그 자리에서 그림
  (꼴에 글꼴 주소 없음 = 글꼴 조각 파일 없음 · 받기 전에 그 글자 글꼴 미리 받음)
- 카메라 주인 = D03 화면 view · 돌림 하나 — 휠 · 가운데 끌기 · 돌림 · 맞춤 · 단축키 모두 draw() → glCamera
- 처음 = GL(aislo:pref:d03.canvas = "1" 이면 옛 캔버스) · 겹 단추줄 끝 단추로 바꿈 · WebGL 못 열거나 잃음 = 저절로 옛 캔버스 + 알림
- GL 이면 캔버스 바탕 겹은 쉼(등고 주곡선 간격 셈만) · GL 이면 위성 한 장 그림 안 받음
- 산사태위험 그림도 경위도 귀퉁이(lonlat_corners 한 곳) · 묶음 꼴 표시 없는 합친 자료도 GL 이 받게 다시 쌈
- 캔버스 그리기 목록 = 새 D03_PlanAnalysis_UI_LayerPaint.ts(700줄 규칙) · 검증 창구 layers().gl = GL 겹별 그린 지물 수
- 잼(headless Edge 실 GPU · 검증_sub2_D · 겹 18 켬 = 자료 있는 겹 전부): 줌아웃 화면 간격 중앙 58.4 → 4.2 ms · 끊김 27 → 0/30
  · 끌기 62.5 → 4.2 ms · 끊김 37 → 0/40 · 마우스 4.2 ms 끊김 0 그대로 · 그은 선 대 GL 0.019 px(돌림 30° 0.019 · 확대 0.077)
- 시험: 캔버스 그리기를 옮긴 자리로 지킴 시험 넷 고침(581-2 · 582 · 584-1 · 584-3 — 그리기 차례 · 색 값이 LayerPaint 로)

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-29 20:28:13 +09:00

625 lines
27 KiB
Python

"""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 · 파생 저장본(열쇠 = 격자 지문)
- 음영기복 = 같은 DEM · B07 `hillshade_png`(PLAN 581-1)
- 보호 · 규제 구역 = 브이월드 WFS 여덟 종(겹 이름 = `_분석_440-C` 표 1) ·
받은 범위로 잘라 종류 · 이름 · 넓이만 · 저장 = `processed/보호구역_bounds.geojson` + 곁 열쇠 ·
종류 하나 실패 = 그 종류만 빠짐
- 산림입지도 = 브이월드 WFS `lt_c_fsdifrsts`(원자료 산림청 산림입지토양도 · PLAN 581-2) ·
받은 범위로 잘라 토양형 코드 · 이름 · 넓이만 ·
저장 = `processed/산림입지도_bounds.geojson` + 곁 열쇠
- 선 겹 = 브이월드 WFS 단층 `lt_l_gimsfault` · 지질구조선 `lt_l_gimslinea`(PLAN 584-1) ·
받은 범위로 잘라 가닥 · 길이만 · 저장 = `processed/<이름>_bounds.geojson` + 곁 열쇠
- 0차곡 = `D03_PlanAnalysis_ZeroOrder`
- 계류 100 m · 도로 300 m 띠 = 도엽 하천중심선(세류 뺌) · 도로중심선(미분류 · 소로 ·
면리간도로 뺌) 합친 띠 · 받은 범위로 자름 · 파생 저장본(열쇠 = 두 파일 지문)
계산 자리 = ③ 서버 단독.
"""
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"
HILLSHADE_KIND = "plan_hillshade"
HILLSHADE_VERSION = 2 # 빈 칸 투명
RASTER_KINDS = ("slope", "hillshade")
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, typename: str) -> list[dict[str, Any]] | None:
"""조각 하나 — B04 WFS 받기(실패 = 파일 없음 · 키 값은 남기지 않음). None = 실패."""
with tempfile.TemporaryDirectory() as folder:
request_vworld_wfs(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, typename: str = CADASTRAL_TYPENAME
) -> tuple[list[dict[str, Any]], int]:
"""(지물, 실패 조각 수) — 1,000 꽉 차면 넷으로 다시 쪼갬."""
features = _wfs(box, epsg, typename)
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, typename)
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": [float(v) 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 hillshade_layer_png(z: np.ndarray, valid: np.ndarray, cell_m: float) -> dict[str, Any]:
"""음영기복 PNG(B07 그대로 · 빈 칸 흰색 255 → 투명 · 값 칸은 90~250)."""
from PIL import Image
from B07_DesignDetail.B07_DesignDetail_Engine_Cad_Lidar import hillshade_png
image, width, height = hillshade_png(z, valid, cell_m)
gray = Image.open(io.BytesIO(base64.b64decode(image.split(",", 1)[1]))).convert("L")
gray.putalpha(gray.point(lambda v: 0 if v == 255 else 255))
buffer = io.BytesIO()
gray.save(buffer, format="PNG", optimize=True)
data = "data:image/png;base64," + base64.b64encode(buffer.getvalue()).decode("ascii")
return {"width": width, "height": height, "image": data, "bands": []}
def slope_layer(
project_root: Path, models_dir: Path, epsg: int, kind: str = "slope"
) -> dict[str, Any]:
"""급경사(`slope`) · 음영기복(`hillshade`) — 지문(격자 파일 · 단계 · 판) 같으면 저장본."""
npz = models_dir / f"dtm_sheet_{SHEET_SURFACE_DEFAULT_METHOD}.npz"
stamp = file_stamp(npz)
if stamp is None:
raise FileNotFoundError("사업 지형면이 없습니다 — 지형 받기 뒤 보입니다.")
store_kind, version, build = (
(SLOPE_KIND, SLOPE_VERSION, slope_bands_png)
if kind == "slope"
else (HILLSHADE_KIND, HILLSHADE_VERSION, hillshade_layer_png)
)
key = derived_key(store_kind, version, npz.name, stamp, int(epsg), SLOPE_BANDS, MAX_SIDE_PX)
saved = read_derived(project_root, store_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 = build(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"] = lonlat_corners(box, epsg)
write_derived(project_root, store_kind, key, value)
return value
def lonlat_corners(box: Any, epsg: int) -> list[list[float]]:
"""TM 상자 → 경위도 네 귀퉁이(왼쪽 아래부터 돎) — GL 바탕 그림 자리(586)."""
from common_util.common_util_crs import tm_to_lonlat
x0, y0, x1, y1 = (float(v) for v in box)
return [list(tm_to_lonlat(px, py, epsg)) for px, py in ((x0, y0), (x1, y0), (x1, y1), (x0, y1))]
# ── 보호 · 규제 구역 ─────────────────────────────────────
# (열쇠, 이름, 브이월드 WFS 겹 이름) — `_분석_440-C` 표 1 그대로
ZONE_KINDS: tuple[tuple[str, str, tuple[str, ...]], ...] = (
("forest", "산림보호구역", ("lt_c_uf151",)),
("baekdu", "백두대간보호지역", ("lt_c_uf901",)),
("water", "상수원보호구역", ("lt_c_um710",)),
("park", "자연공원", ("lt_c_wgisnpgug", "lt_c_wgisnpdo", "lt_c_wgisnpgun")),
("heritage", "국가유산 지정 · 보호구역", ("lt_c_uo301",)),
("wetland", "습지보호지역", ("lt_c_um901",)),
("wildlife", "야생생물보호구역", ("lt_c_um221",)),
("greenbelt", "개발제한구역", ("lt_c_ud801",)),
)
ZONE_FILE = "보호구역_bounds.geojson"
ZONE_KEY = "보호구역_bounds.key.json"
ZONE_VERSION = 1
ZONE_NAME_PROPS = ("alias", "park_name", "uname")
def _transforms(epsg: int) -> tuple[Any, Any]:
"""(경위도 → m, m → 경위도) shapely 변환 함수."""
from pyproj import Transformer
to_m = Transformer.from_crs("EPSG:4326", f"EPSG:{epsg}", always_xy=True).transform
to_ll = Transformer.from_crs(f"EPSG:{epsg}", "EPSG:4326", always_xy=True).transform
return to_m, to_ll
def _clip_metric(feature: dict[str, Any], frame: Any, to_m: Any) -> Any:
"""지물 면(경위도) → 작업 좌표계 m 로 바꿔 범위로 자름 · 못 읽거나 빈 것 = None."""
from shapely.geometry import shape
from shapely.ops import transform
try:
metric = transform(to_m, shape(feature["geometry"]))
except (KeyError, TypeError, ValueError, AttributeError):
return None
if not metric.is_valid:
metric = metric.buffer(0)
part = metric.intersection(frame)
return None if part.is_empty or part.area <= 0 else part
def clip_zones(
found: dict[str, list[dict[str, Any]]], box: Box, epsg: int
) -> tuple[list[dict[str, Any]], dict[str, dict[str, Any]]]:
"""종류별 지물 → 받은 범위로 자른 지물(종류 · 이름) + 종류별 {이름, 개수, 넓이 m²}."""
from shapely.geometry import box as rect
from shapely.geometry import mapping
from shapely.ops import transform
to_m, to_ll = _transforms(epsg)
frame = rect(*box)
features: list[dict[str, Any]] = []
summary: dict[str, dict[str, Any]] = {}
for kind, label, _ in ZONE_KINDS:
count, area = 0, 0.0
for feature in dedup_parcels(found.get(kind) or []):
part = _clip_metric(feature, frame, to_m)
if part is None:
continue
props = feature.get("properties") or {}
name = next((str(props[k]) for k in ZONE_NAME_PROPS if props.get(k)), label)
count += 1
area += part.area
features.append(
{
"type": "Feature",
"geometry": mapping(transform(to_ll, part)),
"properties": {"kind": kind, "name": name, "area_m2": round(part.area, 1)},
}
)
summary[kind] = {"label": label, "count": count, "area_m2": round(area, 1)}
return features, summary
def ensure_zones(processed: Path, box: Box, epsg: int) -> dict[str, Any]:
"""보호 · 규제 구역 — 열쇠 상자가 덮고 실패 종류 없으면 그대로 · 아니면 통째로 다시 받음."""
path, key_path = processed / ZONE_FILE, processed / ZONE_KEY
key = _read_json(key_path) or {}
if (
path.is_file()
and key.get("version") == ZONE_VERSION
and not key.get("failed")
and _covers(key.get("box"), box)
):
return {"path": path, "fetched": False, "failed": []}
jobs = [(kind, name) for kind, _, names in ZONE_KINDS for name in names]
# ponytail: 겹마다 범위 통째 한 번 · 1,000 꽉 차면 넷 — 구역은 필지보다 훨씬 적음
with ThreadPoolExecutor(max_workers=WFS_WORKERS) as pool:
results = list(pool.map(lambda job: fetch_tile(box, epsg, 0, job[1]), jobs))
failed = sorted({kind for (kind, _), (_, bad) in zip(jobs, results) if bad})
if len(failed) == len(ZONE_KINDS):
raise RuntimeError("보호 · 규제 구역을 받지 못했습니다 — 브이월드 WFS 실패")
found: dict[str, list[dict[str, Any]]] = {}
for (kind, _), (got, _) in zip(jobs, results):
found.setdefault(kind, []).extend(got)
features, summary = clip_zones(found, box, epsg)
for kind in failed:
summary[kind]["failed"] = True
processed.mkdir(parents=True, exist_ok=True)
atomic_write_json(path, {"type": "FeatureCollection", "features": features, "summary": summary})
atomic_write_json(
key_path,
{
"version": ZONE_VERSION,
"box": [float(v) for v in box],
"epsg": int(epsg),
"failed": failed,
},
)
logger.info("D 보호구역 받기: 지물 %d · 실패 종류 %s", len(features), failed)
return {"path": path, "fetched": True, "failed": failed}
# ── 산림입지도 ──────────────────────────────────────────
SOIL_TYPENAME = "lt_c_fsdifrsts"
SOIL_FILE = "산림입지도_bounds.geojson"
SOIL_KEY = "산림입지도_bounds.key.json"
SOIL_VERSION = 1
def clip_soils(
found: list[dict[str, Any]], box: Box, epsg: int
) -> tuple[list[dict[str, Any]], list[dict[str, Any]]]:
"""산림입지도 지물 → 받은 범위로 자른 지물(토양형 코드 · 이름 · 넓이) + 토양형별 {코드, 이름,
개수, 넓이 m²}(코드 차례). 이름 없으면 토양형 글자 그대로."""
from shapely.geometry import box as rect
from shapely.geometry import mapping
from shapely.ops import transform
to_m, to_ll = _transforms(epsg)
frame = rect(*box)
features: list[dict[str, Any]] = []
summary: dict[str, dict[str, Any]] = {}
for feature in dedup_parcels(found):
part = _clip_metric(feature, frame, to_m)
if part is None:
continue
props = feature.get("properties") or {}
code = str(props.get("toyanghyun") or props.get("sltp_cd") or "?")
name = str(props.get("name") or code)
row = summary.setdefault(
code,
{"code": code, "sort": str(props.get("sltp_cd") or code), "name": name, "count": 0},
)
row["count"] += 1
row["area_m2"] = row.get("area_m2", 0.0) + part.area
features.append(
{
"type": "Feature",
"geometry": mapping(transform(to_ll, part)),
"properties": {"code": code, "name": name, "area_m2": round(part.area, 1)},
}
)
rows = sorted(summary.values(), key=lambda row: (row["sort"], row["code"]))
for row in rows:
del row["sort"]
row["area_m2"] = round(row["area_m2"], 1)
return features, rows
def ensure_soils(processed: Path, box: Box, epsg: int) -> dict[str, Any]:
"""산림입지도 — 열쇠 상자가 덮으면 그대로 · 아니면 범위 통째로 다시 받음 · 실패 = 오류."""
path, key_path = processed / SOIL_FILE, processed / SOIL_KEY
key = _read_json(key_path) or {}
if path.is_file() and key.get("version") == SOIL_VERSION and _covers(key.get("box"), box):
return {"path": path, "fetched": False}
got, failed = fetch_tile(box, epsg, 0, SOIL_TYPENAME)
if failed:
raise RuntimeError("산림입지도를 받지 못했습니다 — 브이월드 WFS 실패")
features, summary = clip_soils(got, box, epsg)
processed.mkdir(parents=True, exist_ok=True)
atomic_write_json(path, {"type": "FeatureCollection", "features": features, "summary": summary})
atomic_write_json(
key_path, {"version": SOIL_VERSION, "box": [float(v) for v in box], "epsg": int(epsg)}
)
logger.info("D 산림입지도 받기: 지물 %d · 토양형 %d", len(features), len(summary))
return {"path": path, "fetched": True}
# ── 선 겹(단층 · 지질구조선) ─────────────────────────────
# 열쇠 → (이름, 브이월드 WFS 겹 이름) — 새 선 겹은 여기 한 줄(PLAN 584-1)
LINE_KINDS: dict[str, tuple[str, str]] = {
"faults": ("단층", "lt_l_gimsfault"),
"lineaments": ("지질구조선", "lt_l_gimslinea"),
}
LINE_VERSION = 1
def clip_lines(
found: list[dict[str, Any]], box: Box, epsg: int
) -> tuple[list[dict[str, Any]], dict[str, Any]]:
"""선 지물 → 받은 범위로 자른 선(길이 m · 원 속성 `legend`) + {가닥, 길이 m}."""
from shapely.geometry import box as rect
from shapely.geometry import mapping, shape
from shapely.ops import transform
to_m, to_ll = _transforms(epsg)
frame = rect(*box)
features: list[dict[str, Any]] = []
total = 0.0
for feature in dedup_parcels(found):
try:
part = transform(to_m, shape(feature["geometry"])).intersection(frame)
except (KeyError, TypeError, ValueError, AttributeError):
continue
if part.is_empty or part.length <= 0:
continue
total += part.length
props = feature.get("properties") or {}
features.append(
{
"type": "Feature",
"geometry": mapping(transform(to_ll, part)),
"properties": {"legend": props.get("legend"), "length_m": round(part.length, 1)},
}
)
return features, {"count": len(features), "length_m": round(total, 1)}
def ensure_lines(processed: Path, box: Box, epsg: int, kind: str) -> dict[str, Any]:
"""선 겹 — 열쇠 상자가 덮으면 그대로 · 아니면 범위 통째로 다시 받음 · 실패 = 오류."""
label, typename = LINE_KINDS[kind]
path, key_path = processed / f"{label}_bounds.geojson", processed / f"{label}_bounds.key.json"
key = _read_json(key_path) or {}
if path.is_file() and key.get("version") == LINE_VERSION and _covers(key.get("box"), box):
return {"path": path, "fetched": False}
got, failed = fetch_tile(box, epsg, 0, typename)
if failed:
raise RuntimeError(f"{label}을 받지 못했습니다 — 브이월드 WFS 실패")
features, summary = clip_lines(got, box, epsg)
processed.mkdir(parents=True, exist_ok=True)
atomic_write_json(path, {"type": "FeatureCollection", "features": features, "summary": summary})
atomic_write_json(
key_path, {"version": LINE_VERSION, "box": [float(v) for v in box], "epsg": int(epsg)}
)
logger.info("D %s 받기: 선 %d · %.0f m", label, summary["count"], summary["length_m"])
return {"path": path, "fetched": True}
# ── 계류 · 도로 띠 ───────────────────────────────────────
BAND_KIND = "plan_band"
BAND_VERSION = 1
# 타당성평가 「산지계류 100 m 이내」 · 「도로법 도로 300 m 이내」
STREAM_BAND_M = 100.0
ROAD_BAND_M = 300.0
# ponytail: 원문 = 1/25,000 주요 수계 · 도로법 도로 — 1/5,000 도엽에서 세류 · 소로급을 뺀
# 어림(협의 값)
SMALL_STREAMS = frozenset({"세류"})
LOCAL_ROADS = frozenset({"미분류", "소로", "면리간도로"})
BAND_SEGMENTS = 8 # 둥근 끝 사분원 꼭짓점
def band_union(lines: list[Any], dist_m: float, box: Box) -> Any:
"""선(작업 좌표계 m) → 거리 `dist_m` 띠 하나(겹친 곳 합침 · 범위로 자름) · 선 없으면 None."""
from shapely import union_all
from shapely.geometry import box as rect
if not lines:
return None
band = union_all(lines).buffer(dist_m, quad_segs=BAND_SEGMENTS).intersection(rect(*box))
return None if band.is_empty else band
def _band(
path: Path, epsg: int, box: Box, dist_m: float, prop: str, skip: frozenset[str]
) -> dict[str, Any]:
from shapely.geometry import mapping, shape
from shapely.ops import transform
to_m, to_ll = _transforms(epsg)
features = (_read_json(path) or {}).get("features") or []
lines = [
transform(to_m, shape(f["geometry"]))
for f in features
if f.get("geometry") and (f.get("properties") or {}).get(prop) not in skip
]
band = band_union(lines, dist_m, box)
frame = (box[2] - box[0]) * (box[3] - box[1])
area = band.area if band is not None else 0.0
shown = [] if band is None else [band.simplify(0.5)]
return {
"dist_m": dist_m,
"lines": len(lines),
"area_m2": round(area, 1),
"ratio": round(area / frame, 4) if frame > 0 else 0.0,
"collection": {
"type": "FeatureCollection",
"features": [
{"type": "Feature", "geometry": mapping(transform(to_ll, g)), "properties": {}}
for g in shown
],
},
}
def band_layer(project_root: Path, epsg: int, box: Box) -> dict[str, Any]:
"""계류 100 m · 도로 300 m 띠 — 두 도엽 겹 지문 · 범위 같으면 저장본 그대로."""
stream = ensure_sheet_layer(project_root, "하천중심선")
road = ensure_sheet_layer(project_root, "도로중심선")
key = derived_key(
BAND_KIND,
BAND_VERSION,
file_stamp(stream),
file_stamp(road),
int(epsg),
[round(v, 1) for v in box],
STREAM_BAND_M,
ROAD_BAND_M,
sorted(SMALL_STREAMS),
sorted(LOCAL_ROADS),
)
saved = read_derived(project_root, BAND_KIND, key)
if saved is not None:
return saved["value"]
value = {
"stream": _band(stream, epsg, box, STREAM_BAND_M, "구분", SMALL_STREAMS),
"road": _band(road, epsg, box, ROAD_BAND_M, "도로구분", LOCAL_ROADS),
}
write_derived(project_root, BAND_KIND, key, value)
return value