feat(B05): 세류망 흐름 새김 + 32방위 화살표 + 확장 루프 분리

1. 세류선 위 셀이 파랑으로 나오던 문제
   원인: 등고선 TIN 보간면은 실제 물골(thalweg)을 재현하지 못해, 세류선 위
   셀인데도 D8이 옆 사면으로 흘려보내 도로에 닿지 못했다.
   조치 (2가지 함께):
   - burn_stream_flow: 확정된 상류 세류망을 따라 격자 흐름 방향을 강제로 새긴다.
     세류망은 이미 도로를 건너 하류로 빠지는 물길로 확정된 자료다.
   - 사슬 추적의 종결 조건에 세류망 셀을 추가. 물이 세류에 합류한 시점에
     도로 도달이 결정된다(사슬 끝이 도로 셀에 정확히 닿지 않아도 된다).
   - StreamSplit.upstream 을 물 흐름 방향(상류->하류)으로 정렬. 방향은
     _spread_network 확산 시 진입 끝점을 기록해 한 번에 정한다.

2. 화살표 8방위 -> 32방위
   descent_azimuth: 지표면 기울기에서 연속 최급강하 방위를 뽑아 32단계로 양자화.
   D8은 연결(도로 도달 판정)에만 쓰고 표시는 실제 지형 방위를 따른다.
   세류망 새김 셀은 확정된 물길 방향을 그대로 쓴다.
   응답 바이트: 하위 6비트=32방위(32=제자리, 33=표고없음), 0x80=도로 도달.

3. 확장 루프 분리
   Basin._solve_grid 안에 있던 루프를 Flow.expand_until_closed 로 분리.
   단계 검증 미리보기는 이 경로를 타지 않는다.

4. 도로 미도달 셀 화살표를 백색으로 변경(파랑 채움 + 백색 화살표).

합성 검증: 세류 새김 389셀 전부 도로 도달, 32방위 전부 등장,
지류가 본류 중간에 합류하는 갈래를 미연결로 오탐하던 경고 제거.

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
This commit is contained in:
2026-07-31 18:50:21 +09:00
co-authored by Claude Fable 5
parent 06f2328b55
commit 638af3ecc6
7 changed files with 388 additions and 113 deletions
+9 -1
View File
@@ -330,14 +330,22 @@ export interface DrainagePrimaryRegion {
/** 셀별 흐름 방향과 도로 도달 여부. 등고선이 없어 판정을 못하면 null. */
flow: {
encoding: "base64-uint8";
/** 방위 분해능(32). 코드 0 = 화면 오른쪽, 시계방향 증가. */
azimuth_steps: number;
/** 제자리(더 낮은 이웃 없음)를 뜻하는 코드. */
sink_code: number;
/** 표고가 없어 판정 못한 셀 코드. */
invalid_code: number;
cells: number;
reaches_road: number;
no_road: number;
/** 격자에는 있으나 등고선 TIN 밖이라 표고가 없어 판정 못한 셀. */
unanalyzed: number;
/** 확정된 상류 세류망을 따라 흐름을 강제로 새긴 셀 수. */
burned: number;
outer_seeds: number;
interior_seeds: number;
/** 셀당 1바이트. 하위 4비트=3×3 방향코드(4=제자리, 15=무효), 0x80=도로 도달.
/** 셀당 1바이트. 하위 6비트=32방위 코드(32=제자리, 33=무효), 0x80=도로 도달.
* 순서는 grid.row_spans를 행 → 구간 → 열 오름차순으로 훑은 순서와 같다. */
data: string;
} | null;
@@ -38,20 +38,19 @@ from B05_wf2_Route.B05_wf2_Route_Engine_Drainage import (
from B05_wf2_Route.B05_wf2_Route_Engine_Watershed_Flow import (
FlowClassification,
RoadRaster,
border_contact,
burn_stream_flow,
classify_flow,
expand_until_closed,
largest_ring,
outer_boundary,
polygonize_labels,
rasterize_road,
trace_flow,
)
from B05_wf2_Route.B05_wf2_Route_Engine_Watershed_Grid import (
GridSpec,
TerrainGrid,
build_contour_cloud,
build_terrain_grid,
expand_grid_spec,
route_elevation_floor,
)
from B05_wf2_Route.B05_wf2_Route_Engine_Watershed_Stream import (
@@ -254,7 +253,10 @@ def preview_stages(
started = time.perf_counter()
terrain = build_terrain_grid(spec, cloud, region.cell_mask)
road = rasterize_road(spec, route_line, DRAINAGE_ROAD_WIDTH_M)
flow = classify_flow(terrain, road)
# 확정된 상류 세류망을 따라 흐름을 새긴다 — 세류선 위 셀과 그리로 흘러드는 셀은
# 반드시 도로에 도달해야 한다(TIN 보간면은 실제 물골을 재현하지 못한다).
terrain, burned = burn_stream_flow(terrain, road, region.split.upstream)
flow = classify_flow(terrain, road, burned)
logger.info("배수유역: 흐름 판정 %.1fs (셀 %d개)", time.perf_counter() - started, spec.size)
return StagePreview(region=region, terrain=terrain, road=road, flow=flow)
@@ -291,40 +293,18 @@ def _solve_grid(
logger.warning("배수유역: 1차 영역 안에 등고선이 없어 격자 해석을 건너뜁니다.")
return None
terrain = road = flow = None
for round_index in range(DRAINAGE_MAX_EXPAND_ROUNDS + 1):
started = time.perf_counter()
terrain = build_terrain_grid(spec, cloud)
road = rasterize_road(spec, route_line, DRAINAGE_ROAD_WIDTH_M)
flow = trace_flow(terrain, road)
# 격자 크기(config DRAINAGE_GRID_SIZE_M)를 조정할 근거가 되도록 회차별 소요를 남긴다.
logger.info(
"배수유역: %d회차 해석 %.1fs (셀 %d개)",
round_index + 1,
time.perf_counter() - started,
spec.size,
)
contact = border_contact(flow.active)
if not any(contact.values()):
break
if round_index == DRAINAGE_MAX_EXPAND_ROUNDS:
logger.warning(
"배수유역: 확장 상한(%d회)에 도달 — 경계 %s가 아직 활성입니다.",
DRAINAGE_MAX_EXPAND_ROUNDS,
[side for side, touched in contact.items() if touched],
)
break
widened = expand_grid_spec(spec, contact, DRAINAGE_EXPAND_STEP_M)
if widened == spec:
break
logger.info(
"배수유역: 경계 %s 활성 — %.0fm 확장",
[side for side, touched in contact.items() if touched],
DRAINAGE_EXPAND_STEP_M,
)
spec = widened
assert terrain is not None and road is not None and flow is not None
# 확장 루프는 `Watershed_Flow.expand_until_closed()`로 분리했다(2026-07-31 사용자 지시).
# 단계 검증 미리보기(`preview_stages`)는 이 경로를 타지 않는다 — 확장 자체가 아직 검증 대상.
started = time.perf_counter()
expansion = expand_until_closed(spec, cloud, route_line)
spec, terrain, road, flow = expansion.spec, expansion.terrain, expansion.road, expansion.flow
logger.info(
"배수유역: 격자 해석 %.1fs (확장 %d회, %s, 셀 %d개)",
time.perf_counter() - started,
expansion.rounds,
"닫힘" if expansion.closed else "미닫힘",
spec.size,
)
solution = _GridSolution(
spec=spec,
elevation=terrain.elevation.reshape(-1),
@@ -16,6 +16,7 @@ from __future__ import annotations
import logging
from dataclasses import dataclass
from typing import Any
import numpy as np
from rasterio.features import rasterize, shapes
@@ -24,12 +25,18 @@ from shapely.geometry import LineString, MultiPolygon, Polygon, shape
from shapely.ops import unary_union
from B05_wf2_Route.B05_wf2_Route_Engine_Watershed_Grid import (
ContourCloud,
GridSpec,
TerrainGrid,
direction_codes,
build_cell_mask,
build_terrain_grid,
descent_azimuth,
expand_grid_spec,
grid_transform,
)
from config.config_system import (
DRAINAGE_EXPAND_STEP_M,
DRAINAGE_MAX_EXPAND_ROUNDS,
DRAINAGE_MIN_BASIN_AREA_M2,
DRAINAGE_POLYGON_SIMPLIFY_M,
DRAINAGE_ROAD_WIDTH_M,
@@ -171,15 +178,103 @@ def trace_flow(terrain: TerrainGrid, road: RoadRaster) -> FlowResult:
)
def burn_stream_flow(
terrain: TerrainGrid, road: RoadRaster, streams: list[LineString]
) -> tuple[TerrainGrid, np.ndarray]:
"""확정된 상류 세류망을 따라 격자 흐름 방향을 강제로 새긴다.
등고선 TIN 보간면은 실제 물골(thalweg)을 그대로 재현하지 못한다. 그래서 세류선 위
셀인데도 D8이 옆 사면으로 흘려보내 도로에 닿지 못하는 일이 생긴다. 세류선은 이미
"도로를 건너 하류로 빠지는 물길"로 확정된 자료이므로, 그 위 셀의 흐름 방향은 추정할
것이 아니라 **그대로 따라야 한다**(2026-07-31 사용자 지시).
세류선 위 셀은 물길 방향의 다음 셀을 수신 셀로 삼는다. 그러면 세류선 셀은 물론,
세류선으로 흘러드는 사면 셀까지 전부 도로에 도달한다. 도로 셀은 흡수점이므로 건드리지 않는다.
돌려주는 값: (흐름이 새겨진 지형, 새긴 셀 마스크).
"""
spec = terrain.spec
receiver = terrain.receiver.copy()
step_length = terrain.step_length.copy()
valid = terrain.valid.reshape(-1)
road_cells = road.mask.reshape(-1)
burned = np.zeros(spec.size, dtype=bool)
tails: list[int] = []
for line in streams:
chain = _line_cell_chain(spec, line)
if len(chain) < 2:
continue
for current, following in zip(chain, chain[1:]):
if not valid[current] or not valid[following] or road_cells[current]:
continue
if not _is_neighbour(spec, current, following):
continue
receiver[current] = following
step_length[current] = _cell_distance(spec, current, following)
burned[current] = True
tails.append(chain[-1])
# 하류 끝이 도로에도, 다른 세류 갈래에도 닿지 않은 갈래만 진짜 문제다.
# 지류가 본류 중간에 합류하는 경우는 끝점이 본류 셀이므로 정상이다.
detached = sum(1 for tail in tails if not road_cells[tail] and not burned[tail])
if detached:
logger.warning(
"배수유역: 세류망 %d갈래의 하류 끝이 도로·다른 세류 어디에도 닿지 않습니다.",
detached,
)
logger.info("배수유역: 세류망 흐름 새김 %d셀 (세류 %d갈래)", int(burned.sum()), len(streams))
return (
TerrainGrid(
spec=spec,
elevation=terrain.elevation,
valid=terrain.valid,
receiver=receiver,
step_length=step_length,
),
burned,
)
def _line_cell_chain(spec: GridSpec, line: LineString) -> list[int]:
"""선을 따라 지나가는 셀을 순서대로 뽑는다(연속 중복 제거)."""
step = max(spec.cell_m / 2.0, 0.1)
positions = np.arange(0.0, line.length + step, step)
chain: list[int] = []
for position in positions:
point = line.interpolate(float(position))
col = int((point.x - spec.x_min) // spec.cell_m)
row = int((spec.y_max - point.y) // spec.cell_m)
if not (0 <= row < spec.n_rows and 0 <= col < spec.n_cols):
continue
index = row * spec.n_cols + col
if not chain or chain[-1] != index:
chain.append(index)
return chain
def _is_neighbour(spec: GridSpec, first: int, second: int) -> bool:
row_delta = abs(first // spec.n_cols - second // spec.n_cols)
col_delta = abs(first % spec.n_cols - second % spec.n_cols)
return max(row_delta, col_delta) == 1
def _cell_distance(spec: GridSpec, first: int, second: int) -> float:
row_delta = abs(first // spec.n_cols - second // spec.n_cols)
col_delta = abs(first % spec.n_cols - second % spec.n_cols)
return spec.cell_m * float(np.hypot(row_delta, col_delta))
@dataclass
class FlowClassification:
"""셀별 흐름 방향과 도로 도달 여부 — 확장 없이 현재 격자만 본 결과."""
direction: np.ndarray # (R*C,) int8 — 3×3 방향 코드(4=제자리), 무효 셀은 1
direction: np.ndarray # (R*C,) int16 — 32방위 코드(32=제자리, 33=무효)
reaches_road: np.ndarray # (R*C,) bool — 물길을 따라가면 도로에 닿는가
analyzed: np.ndarray # (R*C,) bool — 실제로 판정한 셀
outer_seeds: int # 최외곽에서 출발해 판정한 셀 수
interior_seeds: int # 최외곽 추적에 안 걸려 따로 출발시킨 내부 셀 수
burned: np.ndarray | None = None # (R*C,) bool — 세류망을 따라 흐름을 새긴 셀
def outermost_cells(domain: np.ndarray) -> np.ndarray:
@@ -196,33 +291,54 @@ def outermost_cells(domain: np.ndarray) -> np.ndarray:
return exposed & domain
def classify_flow(terrain: TerrainGrid, road: RoadRaster) -> FlowClassification:
def classify_flow(
terrain: TerrainGrid, road: RoadRaster, burned: np.ndarray | None = None
) -> FlowClassification:
"""최외곽 셀부터 물길을 따라가며 도로 도달 여부를 판정한다.
최외곽 셀에서 출발해 D8 수신 셀을 한 칸씩 따라가고, 경로가 끝나면 그 결과를 경로 전체에
되돌려 적는다. **한 번 판정한 셀은 다시 분석하지 않는다** — 다른 경로가 그 셀을 만나면
거기서 즉시 결론을 가져온다. 최외곽 추적으로 안 닿은 내부 셀은 그다음에 따로 출발시킨다
물이 흐르는 순서 그대로 따라간다 — 셀이 각자 도로를 바라보는 게 아니라, 한 셀에서
출발해 화살표가 가리키는 다음 셀, 그 셀이 가리키는 그다음 셀로 사슬처럼 이어 간다
(2026-07-31 사용자 지시).
경로가 도로 셀에 닿으면 경로 전체가 도로 도달(적색), 싱크에서 멈추거나 해석 영역 밖으로
나가면 미도달(파랑)이다. 격자 확장은 하지 않는다 — 다음 단계다.
사슬이 **도로 셀 또는 확정된 세류망 셀**에 닿으면 그 사슬 전체가 적색이다. 세류망은
이미 "도로를 건너 하류로 빠지는 물길"로 확정된 자료이므로, 물이 세류에 합류한 시점에
도로 도달이 결정된다. 사슬이 싱크에서 멈추거나 해석 영역 밖으로 나가면 전체가 미도달
(파랑 채움 + 백색 화살표)이다.
**한 번 판정한 셀은 다시 분석하지 않는다** — 사슬이 이미 색이 정해진 셀을 만나면 그
셀의 색을 그대로 물려받고 끝낸다. 최외곽 추적으로 안 닿은 내부 셀은 그다음에 따로
출발시킨다.
화살표 방위는 D8(8방위)이 아니라 지형 최급강하 32방위를 쓴다. `burned`(세류망을 따라
흐름을 새긴 셀)는 확정된 물길 방향을 그대로 쓴다.
"""
spec = terrain.spec
valid = terrain.valid.reshape(-1)
receiver = terrain.receiver
is_road = road.mask.reshape(-1) & valid
# 도로 셀과 세류망 셀 둘 다 "여기 닿으면 적색"인 종결점이다.
absorbing = road.mask.reshape(-1) & valid
if burned is not None:
absorbing = absorbing | (burned & valid)
direction = np.where(valid, direction_codes(spec, receiver), -1).astype(np.int8)
direction = descent_azimuth(spec, terrain.elevation, terrain.valid, receiver, burned)
reaches = np.zeros(spec.size, dtype=bool)
# 0=미방문, 1=경로에 올라 있음, 2=판정 완료
state = np.zeros(spec.size, dtype=np.int8)
outer = np.flatnonzero(outermost_cells(terrain.valid))
remaining = np.flatnonzero(valid)
outer_done = _walk_from(outer, receiver, valid, is_road, state, reaches)
interior_done = _walk_from(remaining, receiver, valid, is_road, state, reaches)
outer_done = _walk_from(outer, receiver, valid, absorbing, state, reaches)
interior_done = _walk_from(remaining, receiver, valid, absorbing, state, reaches)
analyzed = state == 2
if burned is not None:
stranded = int((burned & analyzed & ~reaches).sum())
if stranded:
logger.warning(
"배수유역: 세류망 새김 셀 %d개가 도로에 닿지 않습니다 — 하류 끝이 도로 셀과 "
"이어지지 않았는지 확인이 필요합니다.",
stranded,
)
logger.info(
"배수유역: 흐름 판정 %d셀 (최외곽 출발 %d / 내부 보충 %d) — 도로 도달 %d, 미도달 %d",
int(analyzed.sum()),
@@ -237,6 +353,7 @@ def classify_flow(terrain: TerrainGrid, road: RoadRaster) -> FlowClassification:
analyzed=analyzed,
outer_seeds=outer_done,
interior_seeds=interior_done,
burned=burned,
)
@@ -244,11 +361,15 @@ def _walk_from(
starts: np.ndarray,
receiver: np.ndarray,
valid: np.ndarray,
is_road: np.ndarray,
absorbing: np.ndarray,
state: np.ndarray,
reaches: np.ndarray,
) -> int:
"""출발 셀 목록에서 물길을 따라가며 판정한다. 새로 판정한 셀 수를 돌려준다."""
"""출발 셀에서 물길 사슬을 따라가며 판정하고, 사슬 전체에 같은 색을 적는다.
`absorbing`은 도로 셀과 확정된 세류망 셀 — 여기 닿으면 사슬 전체가 도로 도달이다.
새로 판정한 셀 수를 돌려준다.
"""
resolved = 0
path: list[int] = []
for start in starts.tolist():
@@ -258,15 +379,15 @@ def _walk_from(
node = start
while True:
if state[node] == 2:
verdict = bool(reaches[node])
verdict = bool(reaches[node]) # 이미 색이 정해진 셀 — 그 색을 물려받는다
break
if state[node] == 1: # 방어: 채움·평탄해소 후에는 순환이 없어야 한다
verdict = False
break
state[node] = 1
path.append(node)
if is_road[node]:
verdict = True
if absorbing[node]:
verdict = True # 도로 또는 세류망에 합류 — 여기서 하류로 빠진다
break
following = int(receiver[node])
if following == node or not valid[following]:
@@ -280,6 +401,77 @@ def _walk_from(
return resolved
# ── 격자 확장 (별도 단계) ───────────────────────────────────────────────────
@dataclass
class ExpansionResult:
"""확장 루프 결과. 확장을 쓰지 않는 경로에서는 이 모듈을 부르지 않는다."""
spec: GridSpec
terrain: TerrainGrid
road: RoadRaster
flow: FlowResult
rounds: int # 실제로 넓힌 횟수 (0 = 처음부터 닫혀 있었음)
closed: bool # 경계 링이 전부 비활성이 되어 스스로 멈췄는가
def expand_until_closed(
spec: GridSpec,
cloud: ContourCloud,
route_line: LineString,
region_area: Any = None,
road_width_m: float = DRAINAGE_ROAD_WIDTH_M,
max_rounds: int = DRAINAGE_MAX_EXPAND_ROUNDS,
step_m: float = DRAINAGE_EXPAND_STEP_M,
) -> ExpansionResult:
"""활성 셀이 격자 최외곽에 닿은 방향으로만 넓히며 유역이 닫힐 때까지 반복한다.
**본 계산 경로에서만 쓰는 별도 단계다**(2026-07-31 사용자 지시로 분리). 단계 검증
미리보기는 확장 없이 현재 격자만 본다 — 확장 로직 자체가 아직 검증 대상이기 때문이다.
종료 조건은 반경 상한이 아니라 **경계 링 전체가 비활성**이 되는 것이다. 비활성 셀이
최외곽에 띠로 완성되면 그 바깥은 볼 필요가 없다. `max_rounds`는 무한 반복 방지용이다.
`region_area`를 주면 확장된 격자에서도 그 영역에 걸치는 셀만 해석 대상으로 삼는다.
"""
terrain = road = flow = None
rounds = 0
closed = False
for attempt in range(max_rounds + 1):
domain = build_cell_mask(spec, region_area) if region_area is not None else None
terrain = build_terrain_grid(spec, cloud, domain)
road = rasterize_road(spec, route_line, road_width_m)
flow = trace_flow(terrain, road)
contact = border_contact(flow.active)
if not any(contact.values()):
closed = True
break
if attempt == max_rounds:
logger.warning(
"배수유역: 확장 상한(%d회) 도달 — 경계 %s가 아직 활성입니다.",
max_rounds,
[side for side, touched in contact.items() if touched],
)
break
widened = expand_grid_spec(spec, contact, step_m)
if widened == spec:
break
logger.info(
"배수유역: 경계 %s 활성 — %.0fm 확장 (%d회차)",
[side for side, touched in contact.items() if touched],
step_m,
attempt + 1,
)
spec = widened
rounds += 1
assert terrain is not None and road is not None and flow is not None
return ExpansionResult(
spec=spec, terrain=terrain, road=road, flow=flow, rounds=rounds, closed=closed
)
def border_contact(active: np.ndarray) -> dict[str, bool]:
"""활성 셀이 격자 최외곽에 닿은 방향. 전부 False면 유역이 능선 안에서 닫힌 것이다."""
return {
@@ -495,16 +495,58 @@ def build_terrain_grid(
)
def direction_codes(spec: GridSpec, receiver: np.ndarray) -> np.ndarray:
"""수신 셀 인덱스를 3×3 방향 코드로 바꾼다.
# 화살표 방위 분해능. 0 = 화면상 오른쪽(+열), 시계방향으로 증가(행이 아래로 증가하므로).
AZIMUTH_STEPS = 32
# 방위 코드 특수값.
AZIMUTH_SINK = AZIMUTH_STEPS # 32 = 제자리(더 낮은 이웃 없음)
AZIMUTH_INVALID = AZIMUTH_STEPS + 1 # 33 = 표고 없음(해석 불가)
코드 = (행증분 + 1) * 3 + (열증분 + 1) → 0~8. 4는 제자리(싱크)를 뜻한다.
화면이 셀마다 화살표를 그릴 수 있도록 방향만 뽑아낸 표현이다.
def descent_azimuth(
spec: GridSpec,
surface: np.ndarray,
valid: np.ndarray,
receiver: np.ndarray,
forced: np.ndarray | None = None,
) -> np.ndarray:
"""셀별 물 흐름 방위를 32방위 코드로 낸다.
D8은 연결(도로 도달 판정)에는 충분하지만 화면에 8방위밖에 못 그린다. 실제 지표수는
지형 최급강하 방향으로 흐르고 그 방향은 연속값이므로, **표시는 지표면 기울기에서 뽑은
연속 방위를 32단계로 양자화**해 보여 준다(2026-07-31 사용자 지시).
`forced`(세류망을 따라 흐름을 새긴 셀)는 기울기 대신 실제 수신 셀 방향을 쓴다 — 그
셀들은 지형 추정이 아니라 확정된 물길을 따르기 때문이다. 기울기가 0에 가까운 셀도
수신 셀 방향으로 대체한다.
"""
rows, cols = spec.n_rows, spec.n_cols
filled = np.where(valid, surface, np.nan)
# np.gradient는 NaN이 번지므로 무효 셀을 주변 유효값으로 임시 대체한 뒤 기울기를 잡는다.
working = np.where(np.isfinite(filled), filled, np.nanmean(filled) if valid.any() else 0.0)
grad_row, grad_col = np.gradient(working.astype(np.float64), spec.cell_m)
# 내리막 방향 = 기울기 반대. 행은 아래로 증가하므로 화면 좌표와 부호가 같다.
move_row = -grad_row
move_col = -grad_col
magnitude = np.hypot(move_row, move_col)
index = np.arange(receiver.size, dtype=np.int64)
row_delta = receiver // spec.n_cols - index // spec.n_cols
col_delta = receiver % spec.n_cols - index % spec.n_cols
return ((row_delta + 1) * 3 + (col_delta + 1)).astype(np.int8)
receiver_row = (receiver // cols - index // cols).reshape(rows, cols).astype(np.float64)
receiver_col = (receiver % cols - index % cols).reshape(rows, cols).astype(np.float64)
use_receiver = magnitude < 1e-9
if forced is not None:
use_receiver |= forced.reshape(rows, cols)
move_row = np.where(use_receiver, receiver_row, move_row)
move_col = np.where(use_receiver, receiver_col, move_col)
angle = np.arctan2(move_row, move_col)
code = np.rint(angle / (2.0 * math.pi / AZIMUTH_STEPS)).astype(np.int16) % AZIMUTH_STEPS
# 수신 셀이 자기 자신이거나 이동량이 없는 셀은 방향이 없다.
is_sink = (receiver.reshape(rows, cols) == index.reshape(rows, cols)) | (
(np.abs(move_row) < 1e-12) & (np.abs(move_col) < 1e-12)
)
code = np.where(is_sink, AZIMUTH_SINK, code)
code = np.where(valid, code, AZIMUTH_INVALID)
return code.reshape(-1).astype(np.int16)
def route_elevation_floor(route_z_values: list[float]) -> float | None:
@@ -37,7 +37,11 @@ logger = logging.getLogger(__name__)
@dataclass
class StreamSplit:
"""세류망을 도로 교차점에서 상·하류로 가른 결과. 검증 화면이 그대로 그린다."""
"""세류망을 도로 교차점에서 상·하류로 가른 결과. 검증 화면이 그대로 그린다.
`upstream`은 **물이 흐르는 방향(상류 → 하류)으로 정렬**돼 있다. 마지막 좌표가 도로에
가까운 끝이다. 격자 흐름에 세류 방향을 새겨 넣을 때 이 순서를 그대로 쓴다.
"""
upstream: list[LineString] = field(default_factory=list) # 채택 — 1차 영역의 기준
downstream: list[LineString] = field(default_factory=list) # 도로 아래로 이어진 망
@@ -87,20 +91,24 @@ def split_streams_at_road(
node_edges.setdefault(tail, []).append(index)
sampler = ElevationSampler(cloud)
upper_seeds: set[int] = set()
lower_seeds: set[int] = set()
# 씨앗 조각 → 그 조각의 하류쪽 끝점(= 도로 교차 노드). 이 값이 물 흐름 방향의 기준이 된다.
upper_seeds: dict[int, tuple[float, float]] = {}
lower_seeds: dict[int, tuple[float, float]] = {}
for index, piece in enumerate(pieces):
touching = [node for node in ends[index] if node in crossing_nodes]
if not touching:
continue
crossing_z = float(np.min(sampler.at(np.array(touching, dtype=np.float64))))
if _mean_elevation(piece, sampler) > crossing_z:
upper_seeds.add(index)
heights = sampler.at(np.array(touching, dtype=np.float64))
crossing_node = touching[int(np.argmin(heights))]
if _mean_elevation(piece, sampler) > float(np.min(heights)):
upper_seeds[index] = crossing_node
else:
lower_seeds.add(index)
lower_seeds[index] = crossing_node
upstream = _spread_network(upper_seeds, ends, node_edges, crossing_nodes)
downstream = _spread_network(lower_seeds, ends, node_edges, crossing_nodes) - upstream
upstream_flow = _spread_network(upper_seeds, ends, node_edges, crossing_nodes)
downstream_flow = _spread_network(lower_seeds, ends, node_edges, crossing_nodes)
upstream = set(upstream_flow)
downstream = set(downstream_flow) - upstream
logger.info(
"배수유역: 세류 조각 %d개 → 상류망 %d개 채택 / 하류망 %d개 · 미연결 %d개 제외",
len(pieces),
@@ -109,12 +117,26 @@ def split_streams_at_road(
len(pieces) - len(upstream) - len(downstream),
)
return StreamSplit(
upstream=[pieces[index] for index in sorted(upstream)],
# 상류망은 물 흐름 방향(상류 → 하류)으로 뒤집어 둔다 — 격자 흐름 새김에 그대로 쓴다.
upstream=[
_oriented(pieces[index], upstream_flow[index], ends[index])
for index in sorted(upstream)
],
downstream=[pieces[index] for index in sorted(downstream)],
no_contact=len(pieces) - len(upstream) - len(downstream),
)
def _oriented(
piece: LineString,
downstream_node: tuple[float, float],
piece_ends: tuple[tuple[float, float], tuple[float, float]],
) -> LineString:
"""조각을 하류쪽 끝이 마지막 좌표가 되도록 정렬한다."""
head, _tail = piece_ends
return LineString(list(piece.coords)[::-1]) if head == downstream_node else piece
def _cut_network_at_road(
lines: list[LineString], route_line: LineString
) -> tuple[list[LineString], set[tuple[float, float]]]:
@@ -151,13 +173,18 @@ def _cut_network_at_road(
def _spread_network(
seeds: set[int],
seeds: dict[int, tuple[float, float]],
ends: list[tuple[tuple[float, float], tuple[float, float]]],
node_edges: dict[tuple[float, float], list[int]],
blocked: set[tuple[float, float]],
) -> set[int]:
"""씨앗 조각에서 끝점을 타고 퍼진다. 도로 교차 노드는 통과하지 않는다."""
reached = set(seeds)
) -> dict[int, tuple[float, float]]:
"""씨앗 조각에서 끝점을 타고 퍼진다. 도로 교차 노드는 통과하지 않는다.
조각마다 **어느 끝점을 통해 도달했는지**를 함께 기록한다. 그 끝점이 도로에 더 가까운
쪽이므로 곧 그 조각의 하류 방향이다 — 세류망 전체의 물 흐름 방향이 이 한 번의 확산으로
같이 정해진다.
"""
downstream = dict(seeds)
queue = list(seeds)
while queue:
index = queue.pop()
@@ -165,10 +192,11 @@ def _spread_network(
if node in blocked:
continue
for neighbour in node_edges.get(node, ()):
if neighbour not in reached:
reached.add(neighbour)
queue.append(neighbour)
return reached
if neighbour in downstream:
continue
downstream[neighbour] = node
queue.append(neighbour)
return downstream
def _node_key(x: float, y: float) -> tuple[float, float]:
+28 -11
View File
@@ -30,7 +30,12 @@ from B05_wf2_Route.B05_wf2_Route_Engine_Watershed_Basin import (
preview_stages,
)
from B05_wf2_Route.B05_wf2_Route_Engine_Watershed_Export import write_grid_arrays, write_stage
from B05_wf2_Route.B05_wf2_Route_Engine_Watershed_Grid import mask_row_spans
from B05_wf2_Route.B05_wf2_Route_Engine_Watershed_Grid import (
AZIMUTH_INVALID,
AZIMUTH_SINK,
AZIMUTH_STEPS,
mask_row_spans,
)
from B05_wf2_Route.B05_wf2_Route_Repository import (
get_latest_route,
get_route_points,
@@ -282,19 +287,26 @@ def _write_stage_arrays(stored_path: str, preview: Any, region: Any, spec: Any)
flow = preview.flow
if flow is None:
return
arrays = {
"direction": flow.direction.reshape(spec.n_rows, spec.n_cols),
"reaches_road": flow.reaches_road.reshape(spec.n_rows, spec.n_cols),
"analyzed": flow.analyzed.reshape(spec.n_rows, spec.n_cols),
}
if flow.burned is not None:
arrays["burned"] = flow.burned.reshape(spec.n_rows, spec.n_cols)
write_grid_arrays(
stored_path,
"flow_direction",
spec,
arrays,
{
"direction": flow.direction.reshape(spec.n_rows, spec.n_cols),
"reaches_road": flow.reaches_road.reshape(spec.n_rows, spec.n_cols),
"analyzed": flow.analyzed.reshape(spec.n_rows, spec.n_cols),
},
{
"azimuth_steps": AZIMUTH_STEPS,
"sink_code": AZIMUTH_SINK,
"invalid_code": AZIMUTH_INVALID,
"analyzed": int(flow.analyzed.sum()),
"reaches_road": int((flow.reaches_road & flow.analyzed).sum()),
"no_road": int((~flow.reaches_road & flow.analyzed).sum()),
"burned": 0 if flow.burned is None else int(flow.burned.sum()),
"outer_seeds": flow.outer_seeds,
"interior_seeds": flow.interior_seeds,
},
@@ -305,26 +317,31 @@ def _flow_payload(preview: Any, region: Any) -> dict[str, Any] | None:
"""셀별 흐름 방향·도로 도달 여부를 바이트 배열로 압축한다.
셀이 수십만 개라 JSON 객체로는 못 보낸다. 셀 하나당 1바이트로 줄이고 base64로 싣는다:
하위 4비트 = 3×3 방향 코드(0~8, 4=제자리), 15 = 무효(표고 없음)
최상위 비트(0x80) = 도로 도달(적색). 꺼져 있으면 미도달(파랑).
하위 6비트(0x3F) = 32방위 코드(0~31, 0=화면 오른쪽·시계방향), 32=제자리, 33=표고 없음
최상위 비트(0x80) = 도로 도달(적색). 꺼져 있으면 미도달(백색 화살표).
바이트 순서는 `grid.row_spans`를 행 → 구간 → 열 오름차순으로 훑은 순서와 같다.
"""
flow = preview.flow
if flow is None or region.cell_mask is None:
return None
order = np.flatnonzero(region.cell_mask.reshape(-1))
codes = flow.direction[order]
analyzed = flow.analyzed[order]
reaches = flow.reaches_road[order]
packed = np.where(codes < 0, 15, codes).astype(np.uint8)
packed = np.clip(flow.direction[order], 0, AZIMUTH_INVALID).astype(np.uint8)
packed |= np.where(reaches, 0x80, 0).astype(np.uint8)
burned = flow.burned
return {
"encoding": "base64-uint8",
"azimuth_steps": AZIMUTH_STEPS,
"sink_code": AZIMUTH_SINK,
"invalid_code": AZIMUTH_INVALID,
"cells": int(order.size),
"reaches_road": int((reaches & analyzed).sum()),
"no_road": int((~reaches & analyzed).sum()),
# 격자에는 있으나 등고선 TIN 밖이라 표고가 없어 판정 못한 셀(방향 코드 15).
# 격자에는 있으나 등고선 TIN 밖이라 표고가 없어 판정 못한 셀.
"unanalyzed": int((~analyzed).sum()),
# 확정된 상류 세류망을 따라 흐름을 강제로 새긴 셀 수.
"burned": 0 if burned is None else int(burned[order].sum()),
"outer_seeds": flow.outer_seeds,
"interior_seeds": flow.interior_seeds,
"data": base64.b64encode(packed.tobytes()).decode("ascii"),
@@ -53,11 +53,11 @@ const COLLAPSED_KEY = "b05-route-drainage-collapsed";
// 해석 격자 셀 선 — 등고선·세류 위에 얹으므로 흰색으로 둔다(2026-07-31 사용자 지시).
const GRID_LINE_COLOR = "rgba(255, 255, 255, 0.55)";
// 흐름 판정 색 — 도로로 물이 오는 셀은 적색, 오지 않는 셀은 파랑.
// 흐름 판정 색 — 도로로 물이 오는 셀은 적색, 오지 않는 셀은 파랑 채움 + 백색 화살표.
const FLOW_TO_ROAD_FILL = "rgba(220, 38, 38, 0.28)";
const FLOW_TO_ROAD_LINE = "rgba(153, 27, 27, 0.95)";
const FLOW_AWAY_FILL = "rgba(37, 99, 235, 0.22)";
const FLOW_AWAY_LINE = "rgba(30, 64, 175, 0.9)";
const FLOW_AWAY_LINE = "rgba(255, 255, 255, 0.95)";
/** 등고선 TIN 밖이라 표고가 없어 판정하지 못한 셀 — 미도달(파랑)과 구분한다. */
const FLOW_UNKNOWN_FILL = "rgba(120, 113, 108, 0.18)";
/** 셀이 이보다 작으면 화살표가 뭉개져 읽히지 않으므로 채움색만 남긴다(px). */
@@ -357,6 +357,9 @@ export function createDrainagePanel(): DrainagePanel {
}
return;
}
const sink = region.flow?.sink_code ?? 32;
const invalid = region.flow?.invalid_code ?? 33;
const steps = region.flow?.azimuth_steps ?? 32;
for (let offset = 0; offset < count; offset += 1) {
drawFlowCell(
context,
@@ -366,13 +369,14 @@ export function createDrainagePanel(): DrainagePanel {
cellW,
cellH,
cellPx,
{ sink, invalid, steps },
);
}
});
context.restore();
}
/** 셀 하나 — 도달 여부로 칠하고, 충분히 크면 흐름 방향 화살표를 얹는다. */
/** 셀 하나 — 도달 여부로 칠하고, 충분히 크면 32방위 흐름 화살표를 얹는다. */
function drawFlowCell(
context: CanvasRenderingContext2D,
code: number,
@@ -381,11 +385,12 @@ export function createDrainagePanel(): DrainagePanel {
cellW: number,
cellH: number,
cellPx: number,
codes: { sink: number; invalid: number; steps: number },
): void {
const direction = code & 0x0f;
const azimuth = code & 0x3f;
const reaches = (code & 0x80) !== 0;
// 코드 15 = 등고선 TIN 밖이라 표고가 없어 판정 못한 셀. 파랑(미도달)과 구분한다.
const unanalyzed = direction === 15;
// 표고가 없어 판정 못한 셀 — 미도달(파랑 채움)과 구분해야 오독이 없다.
const unanalyzed = azimuth === codes.invalid;
context.fillStyle = unanalyzed
? FLOW_UNKNOWN_FILL
: reaches
@@ -397,11 +402,11 @@ export function createDrainagePanel(): DrainagePanel {
context.lineWidth = 0.5;
context.strokeRect(x, y, cellW, cellH);
}
if (cellPx < ARROW_MIN_PX || direction === 15) return;
if (cellPx < ARROW_MIN_PX || unanalyzed) return;
const stroke = reaches ? FLOW_TO_ROAD_LINE : FLOW_AWAY_LINE;
const midX = x + cellW / 2;
const midY = y + cellH / 2;
if (direction === 4) {
if (azimuth === codes.sink) {
// 제자리(싱크) — 방향이 없으므로 점으로 표시한다.
context.fillStyle = stroke;
context.beginPath();
@@ -409,27 +414,26 @@ export function createDrainagePanel(): DrainagePanel {
context.fill();
return;
}
const colDelta = (direction % 3) - 1;
const rowDelta = Math.floor(direction / 3) - 1;
const length = Math.hypot(colDelta, rowDelta) || 1;
const reach = (cellPx * 0.38) / length;
const tipX = midX + colDelta * reach;
const tipY = midY + rowDelta * reach;
// 코드 0 = 화면 오른쪽(+x), 시계방향(캔버스 y는 아래가 +).
const angle = (azimuth * 2 * Math.PI) / codes.steps;
const unitX = Math.cos(angle);
const unitY = Math.sin(angle);
const reach = cellPx * 0.38;
const tipX = midX + unitX * reach;
const tipY = midY + unitY * reach;
context.strokeStyle = stroke;
context.lineWidth = Math.max(0.6, cellPx * 0.09);
context.beginPath();
context.moveTo(midX - colDelta * reach, midY - rowDelta * reach);
context.moveTo(midX - unitX * reach, midY - unitY * reach);
context.lineTo(tipX, tipY);
context.stroke();
// 촉 — 진행 방향 기준 좌우로 짧게 접는다.
const head = cellPx * 0.16;
const unitX = (colDelta / length) * head;
const unitY = (rowDelta / length) * head;
const head = cellPx * 0.18;
context.beginPath();
context.moveTo(tipX, tipY);
context.lineTo(tipX - unitX - unitY * 0.7, tipY - unitY + unitX * 0.7);
context.lineTo(tipX - (unitX + unitY * 0.7) * head, tipY - (unitY - unitX * 0.7) * head);
context.moveTo(tipX, tipY);
context.lineTo(tipX - unitX + unitY * 0.7, tipY - unitY - unitX * 0.7);
context.lineTo(tipX - (unitX - unitY * 0.7) * head, tipY - (unitY + unitX * 0.7) * head);
context.stroke();
}
@@ -624,11 +628,15 @@ export function createDrainagePanel(): DrainagePanel {
region.flow && region.flow.unanalyzed > 0
? ` / 표고없음 ${region.flow.unanalyzed.toLocaleString()}(회)`
: "";
const burned =
region.flow && region.flow.burned > 0
? ` · 세류망 새김 ${region.flow.burned.toLocaleString()}`
: "";
const flow = region.flow
? ` · 흐름 도로도달 ${region.flow.reaches_road.toLocaleString()}(적) / ` +
`미도달 ${region.flow.no_road.toLocaleString()}(청)${unknown}, ` +
`최외곽 출발 ${region.flow.outer_seeds.toLocaleString()} + ` +
`내부 보충 ${region.flow.interior_seeds.toLocaleString()}`
`내부 보충 ${region.flow.interior_seeds.toLocaleString()}${burned}`
: " · 흐름 판정 없음";
return (
`1차 영역(반경 ${region.radius_m}m): 상류망 ${region.upstream_lines.length}조각 채택 / ` +