auto: 2026-07-29 19:31 (EOMSANGDON-HOME)

B05_wf2_Route/B05_wf2_Route_Engine_Watershed_Trace.py

@@ -10,6 +10,7 @@
 
 from __future__ import annotations
 
+import bisect
 from dataclasses import dataclass
 from typing import Any
 
@@ -256,23 +257,42 @@ def valley_region_polygon(
     if region.is_empty or region.geom_type != "Polygon":
         return None
     if contour_index is not None:
-        region = _snap_boundary_to_ridge(region, contour_index, network_union, others, other_tree)
+        region = _contour_chain_boundary(region, contour_index, network_union, others, other_tree)
     return region
 
 
-def _snap_boundary_to_ridge(
+def _collect_intersection_points(geometry: Any) -> list[Point]:
+    """교차 결과에서 대표 점들을 뽑는다."""
+    if geometry.is_empty:
+        return []
+    if geometry.geom_type == "Point":
+        return [geometry]
+    if geometry.geom_type in {"MultiPoint", "GeometryCollection"}:
+        points: list[Point] = []
+        for part in geometry.geoms:
+            points.extend(_collect_intersection_points(part))
+        return points
+    if geometry.geom_type in {"LineString", "MultiLineString"}:
+        return [geometry.interpolate(0.5, normalized=True)]
+    return []
+
+
+def _contour_chain_boundary(
     region: Any,
     contour_index: Any,
     network_union: Any,
     others: list[LineString],
     other_tree: STRtree | None,
 ) -> Any:
-    """등거리 경계 정점을 등고선 위 능선 꼭짓점으로 스냅한다.
+    """유역 경계를 **등고선마다 능선 포인트 1개씩 찍어 연결**한 체인으로 재구성한다.
 
-    능선 꼭짓점 = 근처 등고선을 따라 걸었을 때 두 세류(우리 상류망·인접 세류)
-    모두로부터의 최소거리가 최대가 되는 지점. 우리 세류 침범(근접)과 인접 세류
-    이탈(인접이 더 가까워짐)은 금지한다. 실패 시 원본 경계를 유지한다.
+    (2026-07-29 사용자 지시: 정점 스냅은 점이 듬성듬성해 등고선을 건너뛴다.)
+    등거리 1차 경계 링은 순서 뼈대로만 쓴다: 링을 가로지르는 모든 등고선 교차점마다
+    그 등고선 위에서 두 세류(우리 상류망·인접 세류) 모두로부터 가장 먼 지점(능선
+    꼭짓점)을 정제해 포인트를 얻고, 링 위 위치 순으로 연결한다. 등고선이 없는 구간은
+    원래 링 정점으로 메운다. 제약: 우리 세류 침범·인접 세류 이탈 금지.
     """
+    ring = LineString(region.exterior.coords)
 
     def _score(point: Point) -> tuple[float, float]:
         to_network = network_union.distance(point)
@@ -283,22 +303,26 @@ def _snap_boundary_to_ridge(
         )
         return to_network, to_other
 
-    snapped: list[tuple[float, float]] = []
-    for x, y in list(region.exterior.coords)[:-1]:
-        vertex = Point(x, y)
-        best = vertex
-        best_net, best_other = _score(vertex)
-        best_value = min(best_net, best_other)
-        for index in contour_index.query(vertex.buffer(RIDGE_SNAP_M)):
-            geom = contour_index.geoms[index]
-            if geom.distance(vertex) > RIDGE_SNAP_M:
-                continue
-            t0 = geom.project(vertex)
+    # 링을 가로지르는 등고선 교차점마다 능선 꼭짓점 1개.
+    chained: list[tuple[float, float, float]] = []  # (링 위치 s, x, y)
+    for index in contour_index.query(ring.buffer(1.0)):
+        geom = contour_index.geoms[index]
+        try:
+            crossings = _collect_intersection_points(geom.intersection(ring))
+        except Exception:  # noqa: BLE001
+            continue
+        for crossing in crossings:
+            s = ring.project(crossing)
+            t0 = geom.project(crossing)
+            best = None
+            best_value = -1.0
             steps = int(RIDGE_WALK_M / 10.0)
-            for offset in [0.0] + [s * 10.0 for k in range(1, steps + 1) for s in (k, -k)]:
+            for offset in [0.0] + [
+                sign * k * 10.0 for k in range(1, steps + 1) for sign in (1, -1)
+            ]:
                 t = min(max(t0 + offset, 0.0), geom.length)
                 candidate = geom.interpolate(t)
-                if candidate.distance(vertex) > RIDGE_SNAP_M + RIDGE_WALK_M:
+                if candidate.distance(crossing) > RIDGE_WALK_M + RIDGE_SNAP_M:
                     continue
                 to_network, to_other = _score(candidate)
                 if to_network < STREAM_JOIN_TOL_M:
@@ -308,10 +332,23 @@ def _snap_boundary_to_ridge(
                 if min(to_network, to_other) > best_value:
                     best_value = min(to_network, to_other)
                     best = candidate
-        snapped.append((best.x, best.y))
-    if len(snapped) < 4:
+            if best is not None:
+                chained.append((s, best.x, best.y))
+    if len(chained) < 4:
         return region
-    polygon = Polygon(snapped).buffer(0)
+    chained.sort()
+    # 등고선 공백 구간(교차점 사이가 먼 곳)은 원래 링 정점으로 메운다.
+    positions = [s for s, _, _ in chained]
+    filled: list[tuple[float, float, float]] = list(chained)
+    for x, y in list(region.exterior.coords)[:-1]:
+        s = ring.project(Point(x, y))
+        slot = bisect.bisect_left(positions, s)
+        before = positions[slot - 1] if slot > 0 else positions[-1] - ring.length
+        after = positions[slot] if slot < len(positions) else positions[0] + ring.length
+        if min(s - before, after - s) > VALLEY_CELL_M * 2.5:
+            filled.append((s, x, y))
+    filled.sort()
+    polygon = Polygon([(x, y) for _, x, y in filled]).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":