From 4fbaf63a602a2d3d90391dd81c31733da3e72415 Mon Sep 17 00:00:00 2001 From: umsangdon Date: Wed, 29 Jul 2026 20:00:50 +0900 Subject: [PATCH] auto: 2026-07-29 20:00 (EOMSANGDON-HOME) --- ...B05_wf2_Route_Engine_Drainage_Watershed.py | 129 ++++++++++++--- .../B05_wf2_Route_Engine_Watershed_Trace.py | 155 ++++++++++++++++++ 2 files changed, 258 insertions(+), 26 deletions(-) diff --git a/B05_wf2_Route/B05_wf2_Route_Engine_Drainage_Watershed.py b/B05_wf2_Route/B05_wf2_Route_Engine_Drainage_Watershed.py index c7697a8e..3a77ccdd 100644 --- a/B05_wf2_Route/B05_wf2_Route_Engine_Drainage_Watershed.py +++ b/B05_wf2_Route/B05_wf2_Route_Engine_Drainage_Watershed.py @@ -20,6 +20,7 @@ from typing import Any from shapely.geometry import LineString, Point, Polygon from shapely.ops import nearest_points, substring, unary_union +from shapely.strtree import STRtree from B05_wf2_Route.B05_wf2_Route_Engine_Drainage import ( StructureCandidate, @@ -33,8 +34,10 @@ from B05_wf2_Route.B05_wf2_Route_Engine_Watershed_Trace import ( ContourIndex, DividerStep, _explode_lines, + rim_walk, side_sign, trace_divider, + trace_ridge_march, trace_upstream_network, valley_region_polygon, ) @@ -373,6 +376,65 @@ def _clip_to_uphill( return merged +def _assemble_march_polygon( + vertices: list[Any], + start_m: float, + end_m: float, + contour_index: ContourIndex, + road_line: LineString, + uphill_sign: int, + network_union: Any, + stream_lines: list[Any], + valley_top_z: float, +) -> Polygon | None: + """개선 2안 — 도로 양끝 물갈림점에서 능선 행진으로 경계를 직접 조립한다. + + 좌·우 능선 행진 체인(세류 제약이 분수령을 강제) + 발원부 위 공통 등고선 아크로 + 폐합(기존 `_assemble_polygon` 재사용). 상류망 커버리지가 목표 미달이면 None을 + 돌려 1안(등거리+체인) 폴백을 태운다. + """ + others = [line for line in stream_lines if line.distance(network_union) > STREAM_JOIN_TOL_M] + other_tree = STRtree(others) if others else None + sx, sy, _ = _interpolate_vertex(vertices, start_m) + ex, ey, _ = _interpolate_vertex(vertices, end_m) + left = trace_ridge_march( + Point(sx, sy), contour_index, network_union, others, other_tree, road_line, uphill_sign + ) + right = trace_ridge_march( + Point(ex, ey), contour_index, network_union, others, other_tree, road_line, uphill_sign + ) + if len(left) < 5 or len(right) < 5: + return None + # 상측 폐합: 우측 정상 → 좌측 정상을 능선마루 보행으로 잇는다. + rim = rim_walk( + right[-1].point, + left[-1].point, + contour_index, + network_union, + others, + other_tree, + road_line, + uphill_sign, + ) + if rim is None: + return None + ring = _road_segment_coords(vertices, start_m, end_m) + ring.extend((step.point.x, step.point.y) for step in right[1:]) + ring.extend((point.x, point.y) for point in rim) + ring.extend((step.point.x, step.point.y) for step in reversed(left[1:])) + if len(ring) < 4: + return None + polygon = Polygon(ring).buffer(0) + if polygon.geom_type == "MultiPolygon": + polygon = max(polygon.geoms, key=lambda part: part.area) + if polygon.is_empty or polygon.geom_type != "Polygon": + return None + coverage = network_union.intersection(polygon).length / max(network_union.length, 1.0) + if coverage < CLOSING_COVERAGE_GOAL: + return None + return polygon + + def build_watershed_basins( vertices: list[Any], candidates: list[StructureCandidate], @@ -449,32 +511,47 @@ def build_watershed_basins( if network_union is not None else None ) - # 도로변 스트립(분할선 폐합) + 세류 계곡 영역(분수령=인접 세류 등거리)을 합친다. - base = _assemble_polygon( - vertices, - divides[position], - divides[position + 1], - dividers[position], - dividers[position + 1], - contour_index, - road_line, - network_union, - valley_top_z, - ) - valley = ( - valley_region_polygon(network_union, stream_features or [], outlet, contour_index) - if network_union is not None - else None - ) - if base is None and valley is None: - continue - if base is not None and valley is not None: - merged = base.union(valley).buffer(0) - if merged.geom_type == "MultiPolygon": - merged = max(merged.geoms, key=lambda part: part.area) - polygon = merged if merged.geom_type == "Polygon" and not merged.is_empty else base - else: - polygon = base if base is not None else valley + # 개선 2안: 도로 양끝 물갈림점에서 능선 행진으로 경계를 직접 그린다. + polygon = None + if network_union is not None and valley_top_z is not None: + polygon = _assemble_march_polygon( + vertices, + divides[position], + divides[position + 1], + contour_index, + road_line, + signs[position], + network_union, + stream_lines, + valley_top_z, + ) + if polygon is None: + # 개선 1안(폴백): 도로변 스트립 + 세류 계곡 영역(등거리+등고선 체인) 합집합. + base = _assemble_polygon( + vertices, + divides[position], + divides[position + 1], + dividers[position], + dividers[position + 1], + contour_index, + road_line, + network_union, + valley_top_z, + ) + valley = ( + valley_region_polygon(network_union, stream_features or [], outlet, contour_index) + if network_union is not None + else None + ) + if base is None and valley is None: + continue + if base is not None and valley is not None: + merged = base.union(valley).buffer(0) + if merged.geom_type == "MultiPolygon": + merged = max(merged.geoms, key=lambda part: part.area) + polygon = merged if merged.geom_type == "Polygon" and not merged.is_empty else base + else: + polygon = base if base is not None else valley if network_union is not None: coverage = network_union.intersection(polygon).length / max(network_union.length, 1.0) if coverage < CLOSING_COVERAGE_GOAL: diff --git a/B05_wf2_Route/B05_wf2_Route_Engine_Watershed_Trace.py b/B05_wf2_Route/B05_wf2_Route_Engine_Watershed_Trace.py index 8604b8ca..a6de8c24 100644 --- a/B05_wf2_Route/B05_wf2_Route_Engine_Watershed_Trace.py +++ b/B05_wf2_Route/B05_wf2_Route_Engine_Watershed_Trace.py @@ -33,6 +33,15 @@ VALLEY_CAP_M = 350.0 # 능선 스냅: 경계 정점에서 등고선 탐색 반경(m)과 등고선 위 능선 꼭짓점 탐색 폭(m). RIDGE_SNAP_M = 40.0 RIDGE_WALK_M = 80.0 +# 능선 행진(개선 2안): 다음 상위 등고선 탐색 반경(m)·등고선 위 꼭짓점 탐색 폭(m)·최대 단계. +# 탐색 폭을 좁게 유지해야 체인이 자기 능선을 국소 추종한다(넓으면 이웃 능선으로 가로 이탈). +MARCH_RADIUS_M = 120.0 +MARCH_WALK_M = 40.0 +MARCH_MAX_STEPS = 150 +# 첫 스텝(시드)만 넓게 탐색 — 물갈림점이 능선 위가 아닐 수 있어 국소 분수령을 먼저 찾는다. +MARCH_SEED_WALK_M = 250.0 +# 이탈 제약 완화 비율 — 등거리선(dn=do) 부근에서 체인이 멈추지 않게 소폭 허용. +MARCH_OTHER_RATIO = 0.7 @dataclass @@ -418,6 +427,152 @@ def _orthogonalize_chain( return entries +def trace_ridge_march( + start: Point, + contour_index: Any, + network_union: Any, + others: list[LineString], + other_tree: STRtree | None, + road_line: LineString, + uphill_sign: int, +) -> list[DividerStep]: + """능선 행진(개선 2안) — 상위 등고선마다 능선 꼭짓점을 한 칸씩 밟아 오른다. + + 수작업 유역도 작도법의 자동화: 물갈림점에서 출발해 매 단계 "반경 안 현재보다 높은 + 등고선 중 가장 낮은 것" 위에서 **|우리 상류망까지 거리 − 인접 세류까지 거리|가 + 최소인 지점**(분수령 = 두 세류망 등거리점)으로 이동한다. 두 세류에서 가장 먼 점을 + 고르면 우리 지류들 사이 내부 지능선으로 새므로, 등거리 조건이 바깥 분수령을 강제 + 한다. 꼭짓점 연결선은 등고선과 자연히 직교한다. + 제약: 도로 산측 유지, 우리 세류 침범(15m)·인접 세류 과이탈 금지. 더 높은 등고선이 + 없으면(능선 정상) 자연 종료. 반환은 DividerStep 목록 — 기존 폐합 로직과 호환. + """ + + def _score(point: Point) -> tuple[float, float]: + to_network = network_union.distance(point) + to_other = ( + others[int(other_tree.nearest(point))].distance(point) + if other_tree is not None + else float("inf") + ) + return to_network, to_other + + z = contour_index.nearest_elevation(start, MARCH_RADIUS_M) + if z is None: + return [] + steps_out = [DividerStep(point=start, z=z, geom_index=-1)] + current = start + for step_no in range(MARCH_MAX_STEPS): + walk_m = MARCH_SEED_WALK_M if step_no == 0 else MARCH_WALK_M + reach_m = max(MARCH_RADIUS_M, walk_m) + best: tuple[float, float, Point, int] | None = None # (레벨, |dn-do|, 지점, geom idx) + for index in contour_index.query(current.buffer(reach_m)): + level = contour_index.zs[index] + if level <= z + 0.01: + continue + if best is not None and level > best[0]: + continue + geom = contour_index.geoms[index] + if geom.distance(current) > reach_m: + continue + t0 = geom.project(current) + walk = int(walk_m / 10.0) + for offset in [0.0] + [sign * k * 10.0 for k in range(1, walk + 1) for sign in (1, -1)]: + t = min(max(t0 + offset, 0.0), geom.length) + candidate = geom.interpolate(t) + if candidate.distance(current) > reach_m: + continue + # 도로 하류측 이탈 금지 — 단, 노선 끝 너머(투영이 끝점에 걸림)는 좌우 + # 부호가 무의미하므로 세류 제약에만 맡긴다(끝을 감아 도는 분수령 허용). + projection = road_line.project(candidate) + if ( + 5.0 < projection < road_line.length - 5.0 + and side_sign(road_line, candidate) == -uphill_sign + ): + continue + to_network, to_other = _score(candidate) + if to_network < STREAM_JOIN_TOL_M: + continue # 우리 세류 침범 금지 + if to_other < to_network * MARCH_OTHER_RATIO: + continue # 인접 세류 쪽 과이탈 금지(등거리선 부근 소폭 허용) + balance = abs(to_network - to_other) + if best is None or level < best[0] or (level == best[0] and balance < best[1]): + best = (level, balance, candidate, index) + if best is None: + break + z = best[0] + current = best[2] + steps_out.append(DividerStep(point=current, z=z, geom_index=best[3])) + return steps_out + + +def rim_walk( + start: Point, + target: Point, + contour_index: Any, + network_union: Any, + others: list[LineString], + other_tree: STRtree | None, + road_line: LineString, + uphill_sign: int, +) -> list[Point] | None: + """능선마루를 따라 두 행진 정상을 잇는다(개선 2안 상측 폐합). + + 좌·우 능선 정상 높이가 달라 단일 등고선 아크로 못 닫는 경우, 매 단계 target에 + 가까워지는 등고선 위 지점 중 |우리 세류 거리 − 인접 세류 거리|가 최소인 곳 + (분수령)으로 이동한다. 레벨 제한 없음(마루는 오르내린다). 막히면 None. + """ + + def _score(point: Point) -> tuple[float, float]: + to_network = network_union.distance(point) + to_other = ( + others[int(other_tree.nearest(point))].distance(point) + if other_tree is not None + else float("inf") + ) + return to_network, to_other + + points: list[Point] = [] + current = start + remaining = current.distance(target) + for _ in range(MARCH_MAX_STEPS): + if remaining <= MARCH_RADIUS_M: + return points + best: tuple[float, Point] | None = None # (|dn-do|, 지점) + for index in contour_index.query(current.buffer(MARCH_RADIUS_M)): + geom = contour_index.geoms[index] + if geom.distance(current) > MARCH_RADIUS_M: + continue + t0 = geom.project(current) + walk = int(MARCH_WALK_M / 10.0) + for offset in [0.0] + [sign * k * 10.0 for k in range(1, walk + 1) for sign in (1, -1)]: + t = min(max(t0 + offset, 0.0), geom.length) + candidate = geom.interpolate(t) + if candidate.distance(current) > MARCH_RADIUS_M: + continue + if candidate.distance(target) > remaining - 5.0: + continue # target에 실질적으로 가까워지는 이동만 허용 + projection = road_line.project(candidate) + if ( + 5.0 < projection < road_line.length - 5.0 + and side_sign(road_line, candidate) == -uphill_sign + ): + continue + to_network, to_other = _score(candidate) + if to_network < STREAM_JOIN_TOL_M: + continue + if to_other < to_network * MARCH_OTHER_RATIO: + continue + balance = abs(to_network - to_other) + if best is None or balance < best[0]: + best = (balance, candidate) + if best is None: + return None + current = best[1] + remaining = current.distance(target) + points.append(current) + return None + + def trace_upstream_network( crossing: Point, stream_features: list[dict[str, Any]],