Files
Aislo/B05_Profile/B05_Profile_Engine_RouteInitial_Metrics.py
T
eomsangdonandClaude Opus 5.5 a718c2d67c feat(B05): 추천 갈래 다시 — 영선 · 안전 우선 · 계류 횡단 최소 · 공사비 어림 최소(PLAN 23-5a)
- 갈래 short · balance 뺌(요청 422) · gradeline · safe · crossings · cost 더함 — 각 목표로 제약 탐색한 후보를 묶음에 넣고 제 비교값으로 고름
- 탐색 목표 셋 더함: 목표 기울기 편차(영선 · 등경사선) · 위험(옆경사 35° 이상 · 계류 버퍼 안 통과 길이) · 계류 횡단 링크 · 공사비는 상대 단가 그대로
- 비교값 더함: side_slope_avg_pct · side_slope_max_pct · hazard_m · stream_crossings · stream_offset_avg_m · grade_dev_pct · cost_index
- 공사비 상대 단가(절성토 1 ㎥ = 1 · 연장 1 m = 2 · 계류 횡단 1 곳 = 400)는 임시 상대값 · 브레인 임의 판정
- 요청 criteria.target_grade_pct(영선 목표 기울기 · 없으면 상한) · 응답 criteria.target_grade_pct
- 시험 test_route_initial 31 · 왕복 시험은 화면이 아는 갈래만(새 넷은 23-5b 뒤)

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_019ACSGaHdLgnkEkoA4LDtMU
2026-09-30 09:12:05 +09:00

306 lines
15 KiB
Python

"""초기 계획노선 갈래의 비교값 · 제약 못 지킨 구간 목록(PLAN 20장).
`B05_Profile_Engine_RouteInitial` 이 만든 폴리라인 하나를 받아 연장 · 종단기울기 · 절성토량 ·
예상노선 이격을 재고, 공통 제약(종단기울기 상한 · 최소 곡선반지름 · 계류 이격)을 못 지킨
구간을 측점 · 값 · 사유로 모은다. `over_limit` 수는 이 목록 수와 같다.
"""
from __future__ import annotations
import math
from typing import TYPE_CHECKING, Any
import numpy as np
import shapely
from scipy.interpolate import RegularGridInterpolator
from shapely.geometry import LineString, MultiLineString, Point
if TYPE_CHECKING:
from B05_Profile.B05_Profile_Engine_RouteInitial import Terrain
from common_util.common_util_route_polyline import RouteCurve
SAMPLE_M = 2.0 # 측정 간격(m)
# 종단기울기를 잴 지반 = DEM 을 이 폭(가우스 σ, m)으로 고른 것 — 도상(1/5,000) 측정 가늠. 지반의
# 잔굴곡(수 m 턱 · 절개면)은 노면 계획고가 흡수하므로 노선 기울기에서 뺀다(사용자 확정 대기).
GRADE_SMOOTH_M = 5.0
STATION_M = 20.0 # 종단기울기 측정 간격(m) — 측점 체계 · 타당성평가 20 m 간격
GRADE_EPS = 1e-6 # 상한과 같은 기울기는 지킨 것으로 봄
# 안전 우선 — 옆경사가 이보다 급한 곳을 위험으로 봄(35° = 70 % · 지식DB 「35° 이상 전폭 절토」
# `산림과임업기술(임도)/2. 임도/4. 노선측량.md`)
SAFE_SIDE_SLOPE = math.tan(math.radians(35.0))
# 공사비 어림 상대 단가(절성토 1 ㎥ = 1) — ⚠ 임시 상대값 · 브레인 임의 판정(PLAN 23-5 · 지식DB
# 근거가 붙으면 바꿈). 연장 1 m = 벌개 · 노반 · 측구 등 절성토 밖 공종을 흙 2 ㎥ 로 어림(이 현장
# 절성토가 1 m 당 3~4 ㎥ 라 두 몫이 같은 자릿수) · 계류 횡단 1 곳 = 횡단 배수 구조물(관 ·
# 기슭막이) 400 ㎥.
COST_EARTH_M3 = 1.0
COST_LENGTH_M = 2.0
COST_CROSSING = 400.0
CROSSING_MERGE_M = 10.0 # 이 거리 안의 계류 교점은 한 횡단(갈라진 도엽 선 이음매)
REASON_TERRAIN = "지형상 불가능 — 복도 안 어떤 길도 {bound:.1f}% 이상 기울기를 지나야 이어짐"
REASON_GRADE_FIT = "곡선 맞춤 뒤 지반 차 — 측점 사이 지반이 상한을 넘음"
REASON_RADIUS = "곡선 둘 자리 모자람 — 앞뒤 직선이 짧아 최소 반지름 원호가 안 들어감"
REASON_STREAM = "계류 버퍼 안 — 건너기 · 기점 · 종점 봐줌을 넘음"
REASON_FOLLOW_GRADE = "예상노선 추적 — 예상노선이 지나는 지반이 상한을 넘음(제약 탐색 안 함)"
REASON_FOLLOW_STREAM = "예상노선 추적 — 예상노선이 계류 버퍼 안을 지남(제약 탐색 안 함)"
REASON_STREAM_TERRAIN = "지형상 불가능 — 기울기 · 반지름을 지키며 계류 버퍼를 비켜 잇는 길 없음"
REASON_WEIGHTED_GRADE = "비중 계산 — 예상노선 몫만큼 제약을 막지 않음(측점 사이 지반이 상한을 넘음)"
REASON_WEIGHTED_STREAM = "비중 계산 — 예상노선 몫만큼 제약을 막지 않음(계류 버퍼 안을 지남)"
def smoothed_ground(terrain: Terrain, grade_smooth_m: float = GRADE_SMOOTH_M) -> np.ndarray:
"""종단기울기를 잴 지반(`grade_smooth_m` 로 고른 DEM) — 탐색과 비교값이 같은 것을 본다.
폭은 화면에서 바꿀 수 있다(PLAN 20장 확정 · 기본값 `GRADE_SMOOTH_M`)."""
from scipy.ndimage import gaussian_filter
spacing = 0.5 * (_spacing(terrain.x) + _spacing(terrain.y))
return gaussian_filter(np.asarray(terrain.z, dtype=np.float64), grade_smooth_m / spacing)
def station_label(chainage_m: float) -> str:
"""측점 표기 — No.(20 m 번호)+나머지(m)."""
number = int(chainage_m // STATION_M)
return f"No.{number}+{chainage_m - number * STATION_M:.1f}"
def _resample(vertices: list[tuple[float, float]], step: float) -> tuple[np.ndarray, np.ndarray]:
"""폴리라인을 `step` 간격 점으로 — (점 [n,2], 누가거리 [n])."""
xy = np.asarray(vertices, dtype=np.float64)
seg = np.hypot(*np.diff(xy, axis=0).T)
chain = np.concatenate([[0.0], np.cumsum(seg)])
total = float(chain[-1])
at = np.append(np.arange(0.0, total, step), total) if total > 0 else np.array([0.0])
return np.column_stack([np.interp(at, chain, xy[:, 0]), np.interp(at, chain, xy[:, 1])]), at
def _runs(flags: np.ndarray) -> list[tuple[int, int]]:
"""참이 이어진 구간들의 (첫, 끝) 번호."""
runs: list[tuple[int, int]] = []
start = None
for i, flag in enumerate(flags):
if flag and start is None:
start = i
elif not flag and start is not None:
runs.append((start, i - 1))
start = None
if start is not None:
runs.append((start, len(flags) - 1))
return runs
def grade_line(ground: np.ndarray, step: float, max_grade: float) -> np.ndarray:
"""상한 안의 계획고 가늠 — 앞으로 · 뒤로 기울기를 누른 두 선의 평균(둘 다 상한 안 → 평균도)."""
rise = max_grade * step
forward, backward = ground.copy(), ground.copy()
for i in range(1, len(ground)):
forward[i] = min(max(ground[i], forward[i - 1] - rise), forward[i - 1] + rise)
for i in range(len(ground) - 2, -1, -1):
backward[i] = min(max(ground[i], backward[i + 1] - rise), backward[i + 1] + rise)
return 0.5 * (forward + backward)
def _item(kind: str, start: float, end: float, value: float, limit: float, reason: str) -> dict:
return {
"kind": kind,
"from_m": round(start, 1),
"to_m": round(end, 1),
"station": station_label(start),
"value": round(value, 2),
"limit": round(limit, 2),
"reason": reason,
}
def route_metrics(
vertices: list[tuple[float, float]],
expected: list[tuple[float, float]],
terrain: Terrain,
streams: list[list[tuple[float, float]]],
curves: list[RouteCurve],
*,
max_grade: float,
min_radius_m: float,
stream_offset_m: float,
road_width_m: float,
bound_pct: float | None,
grade_smooth_m: float = GRADE_SMOOTH_M,
follow: bool = False,
weighted: bool = False,
target_grade: float | None = None,
) -> dict[str, Any]:
"""비교값과 제약 못 지킨 구간 목록. 종단기울기 = 측점(20 m) 사이 고른 지반 기울기.
`bound_pct` 는 제약을 다 지키는 길이 없을 때만 — 복도 안 어떤 길도 넘어야 하는 기울기(%).
`follow` 는 예상노선 추적(제약 탐색 없이 예상노선 그대로) · `weighted` 는 비중 계산(제약을
벌점으로만 따짐) — 사유가 그 뜻으로 뜬다. `target_grade` = 영선 노선 목표 기울기(비율 · 없으면
상한) — 기울기 편차(`grade_dev_pct`)의 기준.
"""
grid = (terrain.y, terrain.x)
sampler = RegularGridInterpolator(grid, terrain.z, bounds_error=False, fill_value=None)
smoothed = smoothed_ground(terrain, grade_smooth_m)
smooth = RegularGridInterpolator(grid, smoothed, bounds_error=False, fill_value=None)
points, chain = _resample(vertices, SAMPLE_M)
length = float(chain[-1])
lookup = points[:, ::-1]
ground = sampler(lookup)
relaxed = bound_pct is not None
terrain_reason = REASON_TERRAIN.format(bound=bound_pct or 0.0)
violations: list[dict[str, Any]] = []
# 종단기울기 — 측점 사이 고른 지반 기울기 · 넘은 측점 구간을 이어 한 항목
stations = np.append(np.arange(0.0, length, STATION_M), length)
station_z = np.interp(stations, chain, smooth(lookup))
raw_z = np.interp(stations, chain, ground)
spans = np.diff(stations)
keep = spans > 0.5
grades = np.abs(np.diff(station_z))[keep] / spans[keep]
raw_grades = np.abs(np.diff(raw_z))[keep] / spans[keep]
starts, ends = stations[:-1][keep], stations[1:][keep]
for first, last in _runs(grades > max_grade + GRADE_EPS):
violations.append(
_item(
"grade",
float(starts[first]),
float(ends[last]),
float(grades[first : last + 1].max()) * 100.0,
max_grade * 100.0,
REASON_FOLLOW_GRADE
if follow
else REASON_WEIGHTED_GRADE
if weighted
else terrain_reason
if relaxed
else REASON_GRADE_FIT,
)
)
# 최소 곡선반지름 — 원호가 하한 아래로 줄어든 곡선
line = LineString(vertices)
for curve in curves:
if curve.violations:
violations.append(
_item(
"radius",
line.project(Point(curve.start)),
line.project(Point(curve.end)),
curve.radius_m,
min_radius_m,
REASON_RADIUS,
)
)
# 계류 이격 — 버퍼 안 구간 중 봐줄 길이를 넘는 것. 봐줌 = 건너기(계류에 닿음) 한 번의
# 버퍼 안 길이(45° 로 가로지름) + 기점 · 종점이 버퍼 안이면 45° 로 빠져나가는 길이
route_points = shapely.points(points)
near = np.full(len(points), np.inf)
crossings = 0
if streams:
network = MultiLineString(streams)
near = shapely.distance(route_points, network)
crossings = _crossings(LineString(vertices), network)
for first, last in _runs(near < stream_offset_m):
ends_in = int(first == 0) + int(last == len(near) - 1)
crosses = float(near[first : last + 1].min()) <= SAMPLE_M
allowance = (ends_in + 2 * crosses) * stream_offset_m * math.sqrt(2.0)
if chain[last] - chain[first] > allowance:
violations.append(
_item(
"stream",
float(chain[first]),
float(chain[last]),
float(near[first : last + 1].min()),
stream_offset_m,
REASON_FOLLOW_STREAM
if follow
else REASON_WEIGHTED_STREAM
if weighted
else REASON_STREAM_TERRAIN
if relaxed
else REASON_STREAM,
)
)
# 절 · 성토(가늠) — 상한 안 계획고와 옆경사 단면을 폭 W 로 적분(비탈면 제외)
grad_y, grad_x = np.gradient(terrain.z, _spacing(terrain.y), _spacing(terrain.x))
gx = RegularGridInterpolator(
(terrain.y, terrain.x), grad_x, bounds_error=False, fill_value=None
)
gy = RegularGridInterpolator(
(terrain.y, terrain.x), grad_y, bounds_error=False, fill_value=None
)
tangent = np.gradient(points, axis=0)
tangent /= np.maximum(np.hypot(tangent[:, 0], tangent[:, 1]), 1e-9)[:, None]
cross = np.abs(gx(lookup) * -tangent[:, 1] + gy(lookup) * tangent[:, 0])
height = grade_line(ground, SAMPLE_M, max_grade) - ground
across = np.linspace(-road_width_m / 2.0, road_width_m / 2.0, 9)
depth = cross[:, None] * across[None, :] - height[:, None] # +: 지반이 노면 위 = 절토
weight = np.gradient(chain) if len(chain) > 1 else np.zeros_like(chain)
cut = float((np.maximum(depth, 0.0).mean(axis=1) * road_width_m * weight).sum())
fill = float((np.maximum(-depth, 0.0).mean(axis=1) * road_width_m * weight).sum())
# 옆경사(고른 지반 · 노선에 수직) · 위험 통과 길이(옆경사 35° 이상 또는 계류 버퍼 안)
sy, sx = np.gradient(smoothed, _spacing(terrain.y), _spacing(terrain.x))
side = np.abs(
RegularGridInterpolator(grid, sx, bounds_error=False, fill_value=None)(lookup)
* -tangent[:, 1]
+ RegularGridInterpolator(grid, sy, bounds_error=False, fill_value=None)(lookup)
* tangent[:, 0]
)
hazard = float(weight[(side >= SAFE_SIDE_SLOPE) | (near < stream_offset_m)].sum())
target = max_grade if target_grade is None else target_grade
deviation = np.abs(grades - target)
offsets = shapely.distance(route_points, LineString(expected))
counts = {
kind: sum(1 for v in violations if v["kind"] == kind)
for kind in ("grade", "radius", "stream")
}
violations.sort(key=lambda item: item["from_m"])
return {
"length_m": round(length, 2),
"max_grade_pct": round(float(grades.max()) * 100.0 if grades.size else 0.0, 2),
"avg_grade_pct": round(
float((grades * spans[keep]).sum() / spans[keep].sum()) * 100.0 if grades.size else 0.0,
2,
),
# 참고 — 고르지 않은 DEM 측점 사이 최대 기울기(잔굴곡 포함 · 판정에 안 씀)
"raw_max_grade_pct": round(float(raw_grades.max()) * 100.0 if raw_grades.size else 0.0, 2),
"cut_m3": round(cut, 1),
"fill_m3": round(fill, 1),
# 흙 남음(절토 − 성토 · + = 남아 버릴 흙 · − = 모자라 들여올 흙) — 흙 균형 갈래의 목표
"surplus_m3": round(cut - fill, 1),
"over_limit": {**counts, "total": len(violations)},
"offset_avg_m": round(float(offsets.mean()), 2),
"offset_max_m": round(float(offsets.max()), 2),
# PLAN 23-5 — 옆경사 · 위험 통과 길이(안전 우선) · 계류 횡단(계류 횡단 최소) · 계류 평균
# 거리(계류 없으면 None) · 목표 기울기 편차(영선 노선 · 측점 기울기와 목표의 차 평균) ·
# 공사비 어림
"side_slope_avg_pct": round(float((side * weight).sum() / max(length, 1e-9)) * 100.0, 2),
"side_slope_max_pct": round(float(side.max()) * 100.0, 2),
"hazard_m": round(hazard, 1),
"stream_crossings": crossings,
"stream_offset_avg_m": round(float(near.mean()), 2) if streams else None,
"grade_dev_pct": round(
float((deviation * spans[keep]).sum() / spans[keep].sum()) * 100.0
if grades.size
else 0.0,
2,
),
"cost_index": round(
COST_EARTH_M3 * (cut + fill) + COST_LENGTH_M * length + COST_CROSSING * crossings
),
"violations": violations,
}
def _crossings(route: LineString, network: MultiLineString) -> int:
"""노선이 계류를 건너는 횟수 — 교점을 노선 따라 세고 `CROSSING_MERGE_M` 안은 한 번."""
hit = route.intersection(network)
points = [geom for geom in getattr(hit, "geoms", [hit]) if not geom.is_empty]
at = sorted(route.project(point.centroid) for point in points)
return sum(1 for i, value in enumerate(at) if i == 0 or value - at[i - 1] > CROSSING_MERGE_M)
def _spacing(coords: np.ndarray) -> float:
return float((coords[-1] - coords[0]) / (len(coords) - 1)) if len(coords) > 1 else 1.0