From 10f175d6b400fbbaad5ffbb02bfe4413605a89f3 Mon Sep 17 00:00:00 2001
From: umsangdon <umsangdon@gmail.com>
Date: Sun, 30 Aug 2026 14:59:50 +0900
Subject: [PATCH] =?UTF-8?q?feat(B04):=20=EB=A7=88=EC=A7=80=EB=A7=89=20?=
 =?UTF-8?q?=EB=93=B1=EA=B3=A0=EC=84=A0=20=EC=95=88=EC=AA=BD=20=ED=8F=89?=
 =?UTF-8?q?=ED=83=84=20=EB=A7=88=EB=A3=A8=EB=A5=BC=20=EA=B3=A8=EA=B2=A9=20?=
 =?UTF-8?q?=EC=BA=A1=EC=9C=BC=EB=A1=9C=20=EB=8F=8B=EC=9A=B4=EB=8B=A4?=
MIME-Version: 1.0
Content-Type: text/plain; charset=UTF-8
Content-Transfer-Encoding: 8bit

다음 등고선이 없는 봉우리·능선 마루는 TIN이 마지막 등고 표고로 평탄하게
채워 'ㅜ'가 갑자기 'ㅡ'로 변하는 단이 생긴다(2026-08-30 사용자 지적).
등표고 평탄 컴포넌트 중 바깥 유효 이웃이 전부 낮은 것(=마루)만 골라,
경계 거리 비례로 중심 골격이 +간격/2(주곡선 5m 기준 +2.5m)가 되게 올린다.
실제 마루는 마지막 등고와 +간격 사이가 보장되므로 오차 상한이 데이터
불확실성 이내다. 골짜기 바닥은 배수 방향을 모르므로 불변.
값별 사전 필터와 컴포넌트 bbox 창 연산으로 전체 생성 11.8s 유지.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
---
 .../B04_PreProcess_Engine_SheetSurface.py     | 72 +++++++++++++++++++
 1 file changed, 72 insertions(+)

diff --git a/B04_PreProcess/B04_PreProcess_Engine_SheetSurface.py b/B04_PreProcess/B04_PreProcess_Engine_SheetSurface.py
index e475fbc3..14bdbd95 100644
--- a/B04_PreProcess/B04_PreProcess_Engine_SheetSurface.py
+++ b/B04_PreProcess/B04_PreProcess_Engine_SheetSurface.py
@@ -164,6 +164,73 @@ def _densify_between_contours(spec: Any, features: list[dict[str, Any]], cloud:
     )
 
 
+def _cap_flat_summits(surface: np.ndarray, cell_m: float, interval_m: float) -> int:
+    """마지막 등고선 안쪽 평탄 마루를 완만한 능선 돔으로 올린다 (in-place).
+
+    다음 등고선이 없는 봉우리·능선 마루는 TIN이 마지막 등고 표고로 평탄하게 채운다
+    ('ㅜ'가 갑자기 'ㅡ'로 변하는 단 — 2026-08-30 사용자 지적). 실제 마루는 그 표고와
+    +간격 사이이므로, 경계에서 가장 먼 골격선이 +간격/2가 되도록 거리 비례로 올린다
+    (기대값 캡, 오차 상한 ±간격/2 = 원천 데이터 불확실성 이내).
+    적용 조건: 등표고 평탄 컴포넌트이고 바깥 유효 이웃이 전부 더 낮을 때(=마루).
+    골짜기 바닥은 배수 방향을 모르므로 건드리지 않는다. 돋운 컴포넌트 수를 반환.
+    """
+    from scipy.ndimage import binary_dilation, distance_transform_edt, label
+
+    finite = np.isfinite(surface)
+    quantized = np.round(surface * 100.0)
+    capped = 0
+    # 평탄 후보: 상하좌우 이웃 중 같은 표고가 있는 셀만 모아 컴포넌트를 짠다.
+    same_right = np.zeros_like(finite)
+    same_right[:, :-1] = finite[:, :-1] & finite[:, 1:] & (quantized[:, :-1] == quantized[:, 1:])
+    same_down = np.zeros_like(finite)
+    same_down[:-1, :] = finite[:-1, :] & finite[1:, :] & (quantized[:-1, :] == quantized[1:, :])
+    flat = same_right | same_down
+    flat[:, 1:] |= same_right[:, :-1]
+    flat[1:, :] |= same_down[:-1, :]
+
+    # 서로 다른 표고의 평탄면이 맞닿아 있을 수 있으므로 표고값별로 컴포넌트를 짠다.
+    # 값별 셀 수 25 미만은 컴포넌트도 25 미만이므로 라벨링 전에 거른다(속도).
+    flat_values, flat_counts = np.unique(quantized[flat], return_counts=True)
+    for value in flat_values[flat_counts >= 25]:
+        value_labels, value_count = label(flat & (quantized == value))
+        sizes = np.bincount(value_labels.reshape(-1))
+        for component_id in range(1, value_count + 1):
+            if sizes[component_id] < 25:  # 25㎡ 미만은 시각적으로 단이 아니다
+                continue
+            rows, cols = np.nonzero(value_labels == component_id)
+            if (
+                rows.min() == 0
+                or cols.min() == 0
+                or rows.max() == surface.shape[0] - 1
+                or cols.max() == surface.shape[1] - 1
+            ):
+                continue  # 격자 가장자리에 닿으면 바깥을 모른다
+            # 이후 연산은 컴포넌트 bbox(+1) 창에서만 — 전체 격자 반복을 피한다.
+            window = (
+                slice(rows.min() - 1, rows.max() + 2),
+                slice(cols.min() - 1, cols.max() + 2),
+            )
+            component = value_labels[window] == component_id
+            ring = binary_dilation(component) & ~component & finite[window]
+            if not ring.any():
+                continue
+            level = float(surface[rows[0], cols[0]])
+            if float(surface[window][ring].max()) >= level - 1e-3:
+                continue  # 더 높은 이웃이 있으면 마루가 아니다(사면 벤치·골짜기)
+            inner = distance_transform_edt(component, sampling=cell_m)
+            peak = float(inner.max())
+            if peak <= 0.0:
+                continue
+            patch = surface[window]
+            patch[component] += (interval_m / 2.0) * (inner[component] / peak).astype(surface.dtype)
+            capped += 1
+    if capped:
+        logger.info(
+            "도엽 서피스: 평탄 마루 %d곳을 +%.1fm 골격 캡으로 돋움", capped, interval_m / 2.0
+        )
+    return capped
+
+
 def build_sheet_surface_model(
     project_root: Path,
     processed_dir: Path,
@@ -193,9 +260,14 @@ def build_sheet_surface_model(
         return None
 
     spec = grid_spec_from_bounds(x_min, y_min, x_max, y_max, SHEET_SURFACE_GRID_M)
+    # 등고 간격(m) — 마루 캡 크기의 근거. 레벨이 하나뿐이면 5m(1:5,000 주곡선) 폴백.
+    contour_levels = np.unique(cloud.z)
+    interval_m = float(np.diff(contour_levels).min()) if len(contour_levels) > 1 else 5.0
     # 계단(평탄 삼각형) 해소 — 인접 등고선 사이 중간 보간선을 정점으로 추가.
     cloud = _densify_between_contours(spec, features, cloud)
     surface = interpolate_elevation(spec, cloud)  # (R, C), 북→남 행 순서, 외부 NaN
+    # 마지막 등고선 안쪽 평탄 마루를 완만한 돔으로 (2026-08-30 사용자 확정).
+    _cap_flat_summits(surface, spec.cell_m, interval_m)
 
     # DtmGridSampler 규약에 맞춰 y 오름차순으로 뒤집어 저장한다.
     x_coords = spec.cell_centers_x()