Merge remote-tracking branches 'origin/sub_laptop_3' and 'origin/sub_laptop_8' into main_laptop_1

This commit is contained in:
2026-11-03 16:36:06 +09:00
69 changed files with 267 additions and 2675 deletions
+1 -1
View File
@@ -66,7 +66,7 @@ async def trigger_wf1_analysis_and_email(
project_info: dict[str, Any] | None = None
try:
# 자동 전처리는 기본 조합만 만든다 — 필터는 엔진이 입력 LAS를 보고 정하고
# (분류 있으면 classification, 없으면 csf), 표현은 SURFACE_AUTO_METHODS.
# (분류 있으면 업체 분류, 없으면 점진 TIN · 둘 다 다듬기 판), 표현은 SURFACE_AUTO_METHODS.
# 다른 조합은 관리자가 B04 드롭다운에서 요청할 때 만든다 (2026-09-01 사용자 확정).
# 확정 표현(implicit)도 함께 지음 — 엔진이 기본 필터에만 지음(요구 558).
methods = list(dict.fromkeys([*SURFACE_AUTO_METHODS, SURFACE_CONFIRM_DEFAULT_METHOD]))
@@ -97,12 +97,6 @@ export interface SurfacePointCloudSampleResponse {
rgb?: [number, number, number][];
}
export interface SurfaceGroundStatsResponse {
status: string;
project_id: string;
filters: Record<string, Record<string, unknown>>;
}
export interface SurfaceStatusResponse {
project_id: string;
status: "pending" | "in_progress" | "completed" | "failed";
@@ -298,14 +292,6 @@ export async function fetchConfirmedSurface(projectId: string): Promise<SurfaceC
});
}
export async function fetchSurfaceGroundStats(
projectId: string,
): Promise<SurfaceGroundStatsResponse> {
return requestJson<SurfaceGroundStatsResponse>(`/projects/${projectId}/surface/ground-stats`, {
method: "GET",
});
}
export async function fetchSurfaceStatus(projectId: string): Promise<SurfaceStatusResponse> {
return requestJson<SurfaceStatusResponse>(`/projects/${projectId}/surface/status`, {
method: "GET",
+7 -3
View File
@@ -14,12 +14,14 @@ from typing import Any
import numpy as np
from B04_PreProcess.B04_PreProcess_Engine_Ground import (
available_filters,
base_filter,
build_ground_masks,
resolve_auto_source_filters,
run_ground_filter,
summarize_masks,
)
from B04_PreProcess.B04_PreProcess_Engine_ModelBuild import BUILDERS
from B04_PreProcess.B04_PreProcess_Engine_Pipeline import build_all_terrain_models
from B04_PreProcess.B04_PreProcess_Engine_Structurize import structurize_las
from common_util.common_util_atomic import atomic_write_npz
@@ -608,8 +610,10 @@ def _collect_analysis_result(
# manifest에서 저장된 모델별 정보 추출 (필터별 대표 모델)
models: list[dict[str, Any]] = []
for filter_key, filter_entry in manifest.get("source_filters", {}).items():
if filter_key not in available_filters():
continue # 지운 갈래(608 · csf 등) 옛 줄은 DB 로 되살리지 않음
for method, meta in filter_entry.get("methods", {}).items():
if meta.get("status") != "completed":
if meta.get("status") != "completed" or method not in BUILDERS:
continue
model_file = meta.get("model_file")
model_path = (models_dir / model_file) if model_file else None
@@ -618,11 +622,11 @@ def _collect_analysis_result(
layers.append(
{
"layer_name": f"{method}_{filter_key}_preview",
"geometry_type": "MESH" if method != "meshfree" else "POINTCLOUD",
"geometry_type": "MESH",
"file_path": _relative_to_project(
project_root, models_dir / meta["preview_file"]
),
"file_format": "glb" if method != "meshfree" else "ply",
"file_format": "glb",
}
)
models.append(
+14 -327
View File
@@ -1,16 +1,13 @@
"""B04 등고선 추출 엔진.
5종 표현(regular_grid/triangular_mesh/bspline_surface/local_rbf_height_field/
meshfree_surfels)의 npz 모델에서 표고 격자를 환원하고, marching squares로
지정 간격 등고선 라인을 추출한다. DTM valid_mask를 footprint로 사용해
경계 누출을 차단한다.
표현(regular_grid · local_rbf_height_field 굳힌 격자)의 표고 격자에서 marching
squares로 지정 간격 등고선 라인을 추출한다. 격자 valid_mask가 경계 누출을 막는다.
"""
from pathlib import Path
from typing import Any
import numpy as np
from scipy.interpolate import RBFInterpolator
from skimage import measure
# 등고선 캐시 형식/추출 규칙이 바뀔 때 증가시킨다.
@@ -139,233 +136,6 @@ def contour_cache_is_stale(contour_path: Path, model_path: Path) -> bool:
return any(path.is_file() and path.stat().st_mtime_ns > built for path in sources)
def _load_footprint_mask(
model_npz_path: Path, x_coords: np.ndarray, y_coords: np.ndarray
) -> np.ndarray | None:
"""같은 source filter의 DTM valid_mask를 현재 격자에 최근접 리샘플한다."""
stem = Path(model_npz_path).stem
if stem.endswith("_smooth"):
stem = stem[:-7]
parts = stem.split("_", 1)
if len(parts) < 2:
return None
filter_key = parts[1]
dtm_path = Path(model_npz_path).parent / f"dtm_{filter_key}.npz"
if not dtm_path.exists():
return None
try:
with np.load(dtm_path) as d:
dtm_x = np.asarray(d["x"]).ravel()
dtm_y = np.asarray(d["y"]).ravel()
dtm_mask = np.asarray(d["valid_mask"], dtype=bool)
except Exception:
return None
if len(dtm_x) < 2 or len(dtm_y) < 2:
return None
def _nearest_idx(axis: np.ndarray, coords: np.ndarray) -> np.ndarray:
ascending = bool(axis[0] <= axis[-1])
a = axis if ascending else axis[::-1]
idx = np.clip(np.searchsorted(a, coords), 1, len(a) - 1)
idx = np.where(np.abs(a[idx - 1] - coords) <= np.abs(a[idx] - coords), idx - 1, idx)
return idx if ascending else (len(axis) - 1 - idx)
xi = _nearest_idx(dtm_x, np.asarray(x_coords, dtype=np.float64))
yi = _nearest_idx(dtm_y, np.asarray(y_coords, dtype=np.float64))
return dtm_mask[np.ix_(yi, xi)]
def _apply_footprint(
model_npz_path: Path, x_coords: np.ndarray, y_coords: np.ndarray, valid_mask: np.ndarray
) -> np.ndarray:
"""valid_mask에 DTM footprint를 교집합으로 적용한다 (형상 다르면 최근접 리샘플)."""
fp = _load_footprint_mask(model_npz_path, x_coords, y_coords)
if fp is not None:
if fp.shape == valid_mask.shape:
return valid_mask & fp
from scipy.ndimage import zoom
zoom_y = valid_mask.shape[0] / fp.shape[0]
zoom_x = valid_mask.shape[1] / fp.shape[1]
fp_resized = zoom(fp.astype(float), (zoom_y, zoom_x), order=0) > 0.5
if fp_resized.shape == valid_mask.shape:
return valid_mask & fp_resized
return valid_mask
def _tin_face_coverage_mask(
vertices: np.ndarray, faces: np.ndarray, xx: np.ndarray, yy: np.ndarray
) -> np.ndarray:
"""저장된 TIN 면이 실제로 덮는 XY 영역만 True로 반환한다."""
vertices = np.asarray(vertices)
faces = np.asarray(faces, dtype=np.int64)
if vertices.ndim != 2 or vertices.shape[1] < 2 or not len(faces):
return np.zeros(xx.shape, dtype=bool)
edges = np.vstack((faces[:, [0, 1]], faces[:, [1, 2]], faces[:, [2, 0]]))
edges = np.sort(edges, axis=1)
unique_edges, counts = np.unique(edges, axis=0, return_counts=True)
boundary_edges = unique_edges[counts == 1]
if not len(boundary_edges):
return np.zeros(xx.shape, dtype=bool)
import affine
import rasterio.features
from shapely import get_parts, linestrings, polygonize
boundary_lines = linestrings(vertices[boundary_edges, :2])
polygons = list(get_parts(polygonize(boundary_lines)))
if not polygons:
return np.zeros(xx.shape, dtype=bool)
x_coords = xx[0, :]
y_coords = yy[:, 0]
dx = float(x_coords[1] - x_coords[0]) if len(x_coords) > 1 else 1.0
dy = float(y_coords[1] - y_coords[0]) if len(y_coords) > 1 else 1.0
transform = affine.Affine(dx, 0.0, x_coords[0] - dx / 2.0, 0.0, dy, y_coords[0] - dy / 2.0)
mask = rasterio.features.rasterize(
polygons, out_shape=xx.shape, transform=transform, fill=0, default_value=1, dtype="uint8"
)
return mask.astype(bool)
def _grid_axes(x_min: float, x_max: float, y_min: float, y_max: float, target_grid_m: float):
cols = max(2, int(np.ceil((x_max - x_min) / target_grid_m)) + 1)
rows = max(2, int(np.ceil((y_max - y_min) / target_grid_m)) + 1)
x_coords = np.linspace(x_min, x_max, cols, dtype=np.float32)
y_coords = np.linspace(y_min, y_max, rows, dtype=np.float32)
return x_coords, y_coords
# 후보 셀 배열이 한 번에 커지지 않도록 삼각형을 나눠 처리하는 단위.
_RASTER_CHUNK = 200_000
# 격자점이 변 위에 놓였을 때 양쪽 삼각형 모두에서 탈락하지 않도록 두는 여유.
_BARYCENTRIC_EPS = 1e-9
def rasterize_triangle_mesh(
vertices: np.ndarray,
triangles: np.ndarray,
x_coords: np.ndarray,
y_coords: np.ndarray,
) -> np.ndarray:
"""삼각망을 규칙격자에 직접 래스터화해 표고 격자를 만든다.
삼각형마다 자기가 덮는 격자 셀만 무게중심 선형보간으로 채우므로, 질의점별
위치탐색(scipy `find_simplex`)이 사라진다. 삼각망 밖 셀은 NaN으로 남아
기존 `griddata` 결과와 같은 유효 판정을 유지한다.
"""
rows, cols = len(y_coords), len(x_coords)
z_grid = np.full((rows, cols), np.nan, dtype=np.float64)
vertices = np.asarray(vertices, dtype=np.float64)
triangles = np.asarray(triangles, dtype=np.int64)
if not len(triangles) or vertices.ndim != 2 or vertices.shape[1] < 3:
return z_grid
axis_x = np.asarray(x_coords, dtype=np.float64)
axis_y = np.asarray(y_coords, dtype=np.float64)
def _index_range(axis: np.ndarray, low: np.ndarray, high: np.ndarray, size: int):
"""축 배열을 직접 탐색해 [low, high]를 덮는 인덱스 구간을 낸다.
`(좌표 - 원점) / 간격` 식은 쓰지 않는다. `_grid_axes`가 float32 축을 만드는데
UTM 좌표(18만대)에서 float32 해상도가 0.015625m라 격자 간격이 0.984~1.0으로
흔들리고, 등간격으로 환산하면 최대 0.81셀까지 어긋난다. 삼각형이 평균 0.3m라
그만큼 어긋나면 격자점을 통째로 놓친다 (2026-08-17 실측).
"""
ascending = size < 2 or axis[-1] >= axis[0]
sorted_axis = axis if ascending else axis[::-1]
lo = np.searchsorted(sorted_axis, low, side="left")
hi = np.searchsorted(sorted_axis, high, side="right") - 1
if ascending:
return lo, hi
return size - 1 - hi, size - 1 - lo
for begin in range(0, len(triangles), _RASTER_CHUNK):
chunk = triangles[begin : begin + _RASTER_CHUNK]
corners = vertices[chunk]
ax, ay, az = corners[:, 0, 0], corners[:, 0, 1], corners[:, 0, 2]
bx, by, bz = corners[:, 1, 0], corners[:, 1, 1], corners[:, 1, 2]
cx, cy, cz = corners[:, 2, 0], corners[:, 2, 1], corners[:, 2, 2]
denom = (by - cy) * (ax - cx) + (cx - bx) * (ay - cy)
col_lo, col_hi = _index_range(
axis_x, np.minimum(np.minimum(ax, bx), cx), np.maximum(np.maximum(ax, bx), cx), cols
)
row_lo, row_hi = _index_range(
axis_y, np.minimum(np.minimum(ay, by), cy), np.maximum(np.maximum(ay, by), cy), rows
)
# 클립 전에 격자 밖 삼각형을 걸러야 경계에 눌러붙지 않는다.
selected = (
np.isfinite(denom)
& (np.abs(denom) > 1e-12)
& (col_hi >= col_lo)
& (row_hi >= row_lo)
& (col_hi >= 0)
& (row_hi >= 0)
& (col_lo <= cols - 1)
& (row_lo <= rows - 1)
)
if not selected.any():
continue
col_lo = np.clip(col_lo[selected], 0, cols - 1).astype(np.int64)
col_hi = np.clip(col_hi[selected], 0, cols - 1).astype(np.int64)
row_lo = np.clip(row_lo[selected], 0, rows - 1).astype(np.int64)
row_hi = np.clip(row_hi[selected], 0, rows - 1).astype(np.int64)
width = col_hi - col_lo + 1
counts = width * (row_hi - row_lo + 1)
total = int(counts.sum())
if total <= 0:
continue
# 삼각형별 후보 셀 수가 제각각이라 반복 인덱스를 펼쳐서 한 번에 계산한다.
starts = np.concatenate(([0], np.cumsum(counts)[:-1]))
local = np.arange(total, dtype=np.int64) - np.repeat(starts, counts)
width_rep = np.repeat(width, counts)
offset_row = local // width_rep
cell_row = np.repeat(row_lo, counts) + offset_row
cell_col = np.repeat(col_lo, counts) + (local - offset_row * width_rep)
point_x = axis_x[cell_col]
point_y = axis_y[cell_row]
ax_r, ay_r, az_r = (
np.repeat(ax[selected], counts),
np.repeat(ay[selected], counts),
np.repeat(az[selected], counts),
)
bx_r, by_r, bz_r = (
np.repeat(bx[selected], counts),
np.repeat(by[selected], counts),
np.repeat(bz[selected], counts),
)
cx_r, cy_r, cz_r = (
np.repeat(cx[selected], counts),
np.repeat(cy[selected], counts),
np.repeat(cz[selected], counts),
)
denom_r = np.repeat(denom[selected], counts)
weight_a = ((by_r - cy_r) * (point_x - cx_r) + (cx_r - bx_r) * (point_y - cy_r)) / denom_r
weight_b = ((cy_r - ay_r) * (point_x - cx_r) + (ax_r - cx_r) * (point_y - cy_r)) / denom_r
weight_c = 1.0 - weight_a - weight_b
hit = (
(weight_a >= -_BARYCENTRIC_EPS)
& (weight_b >= -_BARYCENTRIC_EPS)
& (weight_c >= -_BARYCENTRIC_EPS)
)
if not hit.any():
continue
z_grid[cell_row[hit], cell_col[hit]] = (
weight_a[hit] * az_r[hit] + weight_b[hit] * bz_r[hit] + weight_c[hit] * cz_r[hit]
)
return z_grid
def extract_contours(
model_npz_path: Path,
representation: str,
@@ -415,102 +185,19 @@ def _extract_loaded(
x_coords, y_coords, z_grid, valid_mask, interval, scene_center=scene_center
)
if representation == "triangular_mesh":
vertices, faces = data["vertices"], data["faces"]
x_min, x_max = float(np.min(vertices[:, 0])), float(np.max(vertices[:, 0]))
y_min, y_max = float(np.min(vertices[:, 1])), float(np.max(vertices[:, 1]))
x_coords, y_coords = _grid_axes(x_min, x_max, y_min, y_max, target_grid_m)
xx, yy = np.meshgrid(x_coords, y_coords)
# 저장된 faces를 그대로 쓴다 — griddata는 정점을 다시 삼각분할하느라 긴 변
# 제거·외곽 클리핑이 반영되지 않은 면까지 되살렸다 (PLAN 2026-08-17).
z_grid = rasterize_triangle_mesh(vertices, faces, x_coords, y_coords)
face_mask = _tin_face_coverage_mask(vertices, faces, xx, yy)
valid_mask = np.isfinite(z_grid) & face_mask
valid_mask = _apply_footprint(model_npz_path, x_coords, y_coords, valid_mask)
return extract_contours_from_grid(
x_coords, y_coords, z_grid, valid_mask, interval, scene_center=scene_center
)
if representation == "bspline_surface":
# 모델 빌더와 같은 곡면을 다시 만든다 — 평활 계수·발산 방어를 공유한다.
from B04_PreProcess.B04_PreProcess_Engine_ModelContext import (
evaluate_nurbs_spline,
fit_nurbs_spline,
)
control_x, control_y, control_z = data["control_x"], data["control_y"], data["control_z"]
degree = int(data["degree"][0])
spline, z_range = fit_nurbs_spline(control_x, control_y, control_z, degree)
x_coords, y_coords = _grid_axes(
float(control_x[0]),
float(control_x[-1]),
float(control_y[0]),
float(control_y[-1]),
target_grid_m,
)
z_grid = evaluate_nurbs_spline(
spline, y_coords, x_coords, z_range, Path(model_npz_path).stem
).astype(np.float32)
valid_mask = _apply_footprint(
model_npz_path, x_coords, y_coords, np.ones_like(z_grid, dtype=bool)
)
return extract_contours_from_grid(
x_coords, y_coords, z_grid, valid_mask, interval, scene_center=scene_center
)
if representation == "local_rbf_height_field":
# 굳힌 격자(641-1)가 있으면 그 면 그대로 — 등고선 = 서버 표본 · 종횡단과 같은 면(643)
# 굳힌 격자(641-1) 그대로 — 등고선 = 서버 표본 · 종횡단과 같은 면(643) · dtm 늘 먼저(754-1)
from common_util.common_util_surface_sampler import implicit_grid_npz
filter_key = model_npz_path.stem.removeprefix("implicit_")
if (model_npz_path.parent / f"dtm_{filter_key}.npz").is_file():
from common_util.common_util_surface_sampler import implicit_grid_npz
with np.load(implicit_grid_npz(model_npz_path.parent, filter_key)) as grid:
return extract_contours_from_grid(
grid["x"],
grid["y"],
grid["z"],
grid["valid_mask"],
interval,
scene_center=scene_center,
)
centers_xy, center_z = data["centers_xy"], data["center_z"]
smoothing = float(data["smoothing"][0])
interpolator = RBFInterpolator(
centers_xy.astype(np.float64),
center_z.astype(np.float64),
neighbors=min(64, len(centers_xy)),
smoothing=smoothing,
kernel="thin_plate_spline",
)
x_min, x_max = float(np.min(centers_xy[:, 0])), float(np.max(centers_xy[:, 0]))
y_min, y_max = float(np.min(centers_xy[:, 1])), float(np.max(centers_xy[:, 1]))
x_coords, y_coords = _grid_axes(x_min, x_max, y_min, y_max, target_grid_m)
xx, yy = np.meshgrid(x_coords, y_coords)
z_values = interpolator(np.column_stack([xx.ravel(), yy.ravel()])).astype(np.float32)
z_grid = z_values.reshape(len(y_coords), len(x_coords))
valid_mask = _apply_footprint(
model_npz_path, x_coords, y_coords, np.ones_like(z_grid, dtype=bool)
)
return extract_contours_from_grid(
x_coords, y_coords, z_grid, valid_mask, interval, scene_center=scene_center
)
if representation == "meshfree_surfels":
from scipy.spatial import Delaunay
points = data["points"]
x_min, x_max = float(np.min(points[:, 0])), float(np.max(points[:, 0]))
y_min, y_max = float(np.min(points[:, 1])), float(np.max(points[:, 1]))
x_coords, y_coords = _grid_axes(x_min, x_max, y_min, y_max, target_grid_m)
# Delaunay 1회로 삼각분할과 hull 판정을 함께 해결한다. 기존에는 griddata가
# 내부에서 한 번, hull 판정용으로 또 한 번 만들었다 (PLAN 2026-08-17).
triangulation = Delaunay(np.asarray(points[:, :2], dtype=np.float64))
z_grid = rasterize_triangle_mesh(points, triangulation.simplices, x_coords, y_coords)
# hull 밖은 어느 삼각형에도 안 닿아 NaN으로 남는다 = 기존 hull_inside와 같다.
valid_mask = np.isfinite(z_grid)
valid_mask = _apply_footprint(model_npz_path, x_coords, y_coords, valid_mask)
return extract_contours_from_grid(
x_coords, y_coords, z_grid, valid_mask, interval, scene_center=scene_center
)
with np.load(implicit_grid_npz(model_npz_path.parent, filter_key)) as grid:
return extract_contours_from_grid(
grid["x"],
grid["y"],
grid["z"],
grid["valid_mask"],
interval,
scene_center=scene_center,
)
raise ValueError(f"지원하지 않는 표현 방식입니다: {representation}")
@@ -1,122 +0,0 @@
"""B04 CSF(Cloth Simulation Filter) 지면 필터.
Pure NumPy 기반 CSF. 지형을 반전시킨 뒤 가상의 천을 중력으로 낙하시켜
원 지면(반전 최하단)에 밀착시키고, 천과의 오차가 임계값 이내인 포인트를
지면으로 분류한다.
"""
import math
from typing import Any
import numpy as np
from config.config_system import (
SURFACE_CSF_CLASS_THRESHOLD_M,
SURFACE_CSF_CLOTH_RESOLUTION_M,
SURFACE_CSF_ITERATIONS,
SURFACE_CSF_MAX_ITERATIONS,
SURFACE_CSF_RIGIDNESS,
SURFACE_CSF_SLOPE_SMOOTH,
SURFACE_CSF_SLOPE_SMOOTH_THRESHOLD_M,
SURFACE_CSF_TIME_STEP,
)
# rigidness 단계별 스프링 완화 계수 (1: 산악 밀착 ~ 3: 부드럽게 덮음)
_RIGIDNESS_SPRING_COEFF = {1: 0.25, 2: 0.45, 3: 0.65}
# 중력 하강 계수 (time_step에 곱해 반복당 낙하량 산출)
_GRAVITY_BASE = 9.8 * 0.05
def filter_csf(
structured_data: dict[str, Any] | np.lib.npyio.NpzFile,
cloth_resolution: float = SURFACE_CSF_CLOTH_RESOLUTION_M,
rigidness: int = SURFACE_CSF_RIGIDNESS,
time_step: float = SURFACE_CSF_TIME_STEP,
class_threshold: float = SURFACE_CSF_CLASS_THRESHOLD_M,
iterations: int = SURFACE_CSF_ITERATIONS,
slope_smooth: bool = SURFACE_CSF_SLOPE_SMOOTH,
slope_smooth_threshold: float = SURFACE_CSF_SLOPE_SMOOTH_THRESHOLD_M,
max_iterations: int = SURFACE_CSF_MAX_ITERATIONS,
) -> np.ndarray:
"""CSF로 지면 포인트의 불리언 마스크를 반환한다."""
if not math.isfinite(cloth_resolution) or cloth_resolution <= 0:
raise ValueError("천 해상도는 0보다 큰 유한한 값이어야 합니다.")
if rigidness not in _RIGIDNESS_SPRING_COEFF:
raise ValueError("rigidness는 1, 2, 3 중 하나여야 합니다.")
if not math.isfinite(time_step) or time_step <= 0:
raise ValueError("time_step은 0보다 큰 유한한 값이어야 합니다.")
if not math.isfinite(class_threshold) or class_threshold < 0:
raise ValueError("분류 임계값은 0 이상의 유한한 값이어야 합니다.")
if iterations <= 0:
raise ValueError("반복 횟수는 1 이상이어야 합니다.")
if max_iterations < iterations:
raise ValueError("반복 상한은 정착 여유 반복 수보다 작을 수 없습니다.")
xyz = np.asarray(structured_data["xyz"], dtype=np.float64)
if xyz.ndim != 2 or xyz.shape[1] != 3:
raise ValueError("xyz 배열은 (N, 3) 형태여야 합니다.")
if xyz.shape[0] == 0:
return np.zeros(0, dtype=bool)
xs, ys, zs = xyz[:, 0], xyz[:, 1], xyz[:, 2]
# 1. 지형 반전 — 지표면 추출을 위해 높이를 뒤집는다.
z_max = float(np.max(zs))
inverted_zs = (z_max - zs).astype(np.float32)
# 2. 2D 가상 천 격자 설정 (바운더리 밀착 매핑)
x_min, x_max = float(np.min(xs)), float(np.max(xs))
y_min, y_max = float(np.min(ys)), float(np.max(ys))
cols = int(np.ceil((x_max - x_min) / cloth_resolution)) + 1
rows = int(np.ceil((y_max - y_min) / cloth_resolution)) + 1
# 천 노드 초기 높이 — 반전 지형 최고점보다 약간 높은 곳에서 낙하 시작
start_height = float(np.max(inverted_zs)) + 1.0
cloth_z = np.full((rows, cols), start_height, dtype=np.float32)
# 3. 격자 충돌 타겟 구성 (Drape Target)
collision_grid = np.full((rows, cols), -np.inf, dtype=np.float32)
gx = np.clip(((xs - x_min) / cloth_resolution).astype(np.int32), 0, cols - 1)
gy = np.clip(((ys - y_min) / cloth_resolution).astype(np.int32), 0, rows - 1)
np.maximum.at(collision_grid, (gy, gx), inverted_zs)
collision_grid[collision_grid == -np.inf] = 0.0
# 4. 천 시뮬레이션 반복 루프 (물리 하강)
# 반복 수는 현장 기복에서 뽑는다 — 고정 횟수는 곧 고정 하강 예산이라
# (0.3185m × 150회 = 47.8m) 기복이 그보다 큰 산악지에서 천이 지면에 닿지
# 못하고 공중에 멈춘다. 낙하 구간 + 스프링 정착 여유(iterations)로 잡는다.
gravity = _GRAVITY_BASE * time_step
spring_coeff = _RIGIDNESS_SPRING_COEFF[rigidness]
required_drop = start_height - float(np.min(collision_grid))
total_iterations = int(math.ceil(max(required_drop, 0.0) / gravity)) + iterations
total_iterations = min(total_iterations, max_iterations)
for _ in range(total_iterations):
cloth_z -= gravity
cloth_z = np.maximum(cloth_z, collision_grid)
# 노드 간 스프링 제약 완화 (가로/세로 인접 교정)
diff_h = cloth_z[:, 1:] - cloth_z[:, :-1]
correction_h = diff_h * spring_coeff * 0.5
cloth_z[:, :-1] += correction_h
cloth_z[:, 1:] -= correction_h
diff_v = cloth_z[1:, :] - cloth_z[:-1, :]
correction_v = diff_v * spring_coeff * 0.5
cloth_z[:-1, :] += correction_v
cloth_z[1:, :] -= correction_v
cloth_z = np.maximum(cloth_z, collision_grid)
# 5. 시뮬레이션 천 높이와 원본 대조 → 오차 이내면 지면
simulated_inverted_z = cloth_z[gy, gx]
height_diff = np.abs(inverted_zs - simulated_inverted_z)
mask = height_diff <= class_threshold
# 6. 수목 노이즈 2차 필터 보정 — 셀 지면 후보(반전고 최대 = 원표고 최저)보다
# 얼마나 위에 떠 있는지를 본다. 반전 공간에서 높이는 (지면 후보 - 자기 반전고)다.
# 반대로 빼면 항상 0 이하가 나와 조건이 늘 참이 된다(= 필터가 무력화).
if slope_smooth:
local_ground_z = collision_grid[gy, gx]
mask = mask & ((local_ground_z - inverted_zs) < slope_smooth_threshold)
return np.asarray(mask, dtype=bool)
@@ -1,50 +0,0 @@
"""B04 grid minimum-Z 지면 필터."""
import math
from typing import Any
import numpy as np
from config.config_system import SURFACE_GRID_CELL_SIZE_M, SURFACE_GRID_HEIGHT_THRESHOLD_M
def filter_grid_min_z(
structured_data: dict[str, Any] | np.lib.npyio.NpzFile,
cell_size: float = SURFACE_GRID_CELL_SIZE_M,
height_threshold: float = SURFACE_GRID_HEIGHT_THRESHOLD_M,
) -> np.ndarray:
"""격자 최저 표면에서 높이 임계값 이내인 포인트의 불리언 마스크를 반환한다."""
if not math.isfinite(cell_size) or cell_size <= 0:
raise ValueError("격자 크기는 0보다 큰 유한한 값이어야 합니다.")
if not math.isfinite(height_threshold) or height_threshold < 0:
raise ValueError("높이 임계값은 0 이상의 유한한 값이어야 합니다.")
xyz = np.asarray(structured_data["xyz"], dtype=np.float64)
bounds = np.asarray(structured_data["bounds"], dtype=np.float64)
if xyz.ndim != 2 or xyz.shape[1] != 3:
raise ValueError("xyz 배열은 (N, 3) 형태여야 합니다.")
if bounds.shape != (3, 2):
raise ValueError("bounds 배열은 (3, 2) 형태여야 합니다.")
if xyz.shape[0] == 0:
return np.zeros(0, dtype=bool)
x_min, y_min = bounds[0, 0], bounds[1, 0]
x_max, y_max = bounds[0, 1], bounds[1, 1]
grid_width = int(np.ceil((x_max - x_min) / cell_size)) + 2
grid_height = int(np.ceil((y_max - y_min) / cell_size)) + 2
# float32로 담으면 셀 최저점(가장 확실한 지면점)이 반올림 때문에 자기 기준면보다
# 낮아져 height_above < 0 으로 탈락한다 — 원본 좌표와 같은 float64로 둔다.
minimum_z = np.full((grid_height, grid_width), np.inf, dtype=np.float64)
grid_x = np.clip(((xyz[:, 0] - x_min) / cell_size).astype(np.int64), 0, grid_width - 1)
grid_y = np.clip(((xyz[:, 1] - y_min) / cell_size).astype(np.int64), 0, grid_height - 1)
np.minimum.at(minimum_z, (grid_y, grid_x), xyz[:, 2])
# 기준면은 점이 실제로 들어 있는 셀에서만 읽는다 — 3x3 최소필터로 이웃 셀 값을
# 끌어오면 급경사에서 기준면이 경사만큼 파고들어(중앙 1.88m > 임계 1.5m) 지면점이
# 통째로 탈락한다. 빈 셀은 아래 인덱싱에서 참조되지 않으므로 채울 필요가 없다.
height_above = xyz[:, 2] - minimum_z[grid_y, grid_x]
return np.asarray(
(height_above >= 0.0) & (height_above <= height_threshold),
dtype=bool,
)
@@ -1,91 +0,0 @@
"""B04 PMF(Progressive Morphological Filter) 지면 필터.
XY 평면을 격자로 투영해 Z-min 지형 맵을 만든 뒤, 윈도우 폭을 단계적으로
키워가며 형태학적 열림(Opening) 연산으로 수목·구조물을 제거하고, 최종 지면
대비 높이차가 임계값 이내인 포인트를 지면으로 분류한다. Pure NumPy 구현.
"""
import math
from typing import Any
import numpy as np
from numpy.lib.stride_tricks import sliding_window_view
from config.config_system import (
SURFACE_PMF_CELL_SIZE_M,
SURFACE_PMF_INITIAL_WINDOW_SIZE,
SURFACE_PMF_MAX_DISTANCE_M,
SURFACE_PMF_MAX_WINDOW_SIZE,
SURFACE_PMF_SLOPE,
)
def _min_max_window_filter(grid: np.ndarray, w_size: int, mode: str = "min") -> np.ndarray:
"""순수 NumPy 2D 이동 윈도우 최솟값/최댓값 필터 (경계는 edge 패딩)."""
pad_val = w_size // 2
padded = np.pad(grid, pad_val, mode="edge")
windows = sliding_window_view(padded, (w_size, w_size))
if mode == "min":
return np.min(windows, axis=(2, 3))
return np.max(windows, axis=(2, 3))
def filter_pmf(
structured_data: dict[str, Any] | np.lib.npyio.NpzFile,
cell_size: float = SURFACE_PMF_CELL_SIZE_M,
max_window_size: int = SURFACE_PMF_MAX_WINDOW_SIZE,
slope: float = SURFACE_PMF_SLOPE,
initial_window_size: int = SURFACE_PMF_INITIAL_WINDOW_SIZE,
max_distance: float = SURFACE_PMF_MAX_DISTANCE_M,
) -> np.ndarray:
"""PMF로 지면 포인트의 불리언 마스크를 반환한다."""
if not math.isfinite(cell_size) or cell_size <= 0:
raise ValueError("격자 크기는 0보다 큰 유한한 값이어야 합니다.")
if initial_window_size < 1 or max_window_size < initial_window_size:
raise ValueError("윈도우 크기는 1 이상이고 최대가 초기값 이상이어야 합니다.")
if not math.isfinite(max_distance) or max_distance < 0:
raise ValueError("최대 거리 임계값은 0 이상의 유한한 값이어야 합니다.")
xyz = np.asarray(structured_data["xyz"], dtype=np.float64)
bounds = np.asarray(structured_data["bounds"], dtype=np.float64)
if xyz.ndim != 2 or xyz.shape[1] != 3:
raise ValueError("xyz 배열은 (N, 3) 형태여야 합니다.")
if bounds.shape != (3, 2):
raise ValueError("bounds 배열은 (3, 2) 형태여야 합니다.")
if xyz.shape[0] == 0:
return np.zeros(0, dtype=bool)
xs, ys, zs = xyz[:, 0], xyz[:, 1], xyz[:, 2]
x_min, y_min, z_min = bounds[0, 0], bounds[1, 0], bounds[2, 0]
x_max, y_max = bounds[0, 1], bounds[1, 1]
grid_w = int(np.ceil((x_max - x_min) / cell_size)) + 2
grid_h = int(np.ceil((y_max - y_min) / cell_size)) + 2
z_grid = np.full((grid_h, grid_w), np.inf, dtype=np.float32)
gx = np.clip(((xs - x_min) / cell_size).astype(np.int32), 0, grid_w - 1)
gy = np.clip(((ys - y_min) / cell_size).astype(np.int32), 0, grid_h - 1)
# 1. Z-min 지형 구성
np.minimum.at(z_grid, (gy, gx), zs.astype(np.float32))
z_grid[z_grid == np.inf] = z_min
# 2. 점진적 형태학 필터 (Opening = Dilation of Erosion)
current_grid = z_grid.copy()
w_sizes = []
w = initial_window_size
while w <= max_window_size:
w_sizes.append(w)
w = w * 2 + 1
for w_size in w_sizes:
eroded = _min_max_window_filter(current_grid, w_size, mode="min")
opened = _min_max_window_filter(eroded, w_size, mode="max")
t_dist = min(slope * w_size * cell_size * 0.15 + 0.5, max_distance)
mask_elev = (current_grid - opened) > t_dist
current_grid[mask_elev] = opened[mask_elev]
# 3. 마스크 매핑 및 원본 비교
simulated_z = current_grid[gy, gx]
mask = (zs >= simulated_z - 0.4) & (zs <= simulated_z + max_distance)
return np.asarray(mask, dtype=bool)
@@ -1,116 +0,0 @@
"""B04 RANSAC 지면 필터 (Local plane fitting).
공간을 로컬 격자로 분할하고, 각 격자에서 RANSAC 평면 피팅을 수행해
평면과의 거리가 임계값 이내인 인라이어(지면)를 취합한다. Pure NumPy 구현.
"""
import math
from collections.abc import Callable
from typing import Any
import numpy as np
from config.config_system import (
SURFACE_RANSAC_DISTANCE_THRESHOLD_M,
SURFACE_RANSAC_ITERATIONS,
SURFACE_RANSAC_LOCAL_GRID_SIZE_M,
SURFACE_RANSAC_N,
SURFACE_RANSAC_SEED,
)
def fit_plane_ransac(
points: np.ndarray,
distance_threshold: float = SURFACE_RANSAC_DISTANCE_THRESHOLD_M,
ransac_n: int = SURFACE_RANSAC_N,
num_iterations: int = SURFACE_RANSAC_ITERATIONS,
seed: int = SURFACE_RANSAC_SEED,
) -> np.ndarray:
"""평면 ax+by+cz+d=0을 RANSAC으로 피팅하고 인라이어 마스크를 반환한다."""
n_points = len(points)
if n_points < ransac_n:
return np.ones(n_points, dtype=bool)
best_inliers = np.zeros(n_points, dtype=bool)
max_inlier_count = -1
rng = np.random.default_rng(seed)
for _ in range(num_iterations):
idx = rng.choice(n_points, ransac_n, replace=False)
p0, p1, p2 = points[idx]
normal = np.cross(p1 - p0, p2 - p0)
norm = np.linalg.norm(normal)
if norm < 1e-6:
continue # 세 점이 일직선상 → 스킵
normal = normal / norm
d = -np.dot(normal, p0)
distances = np.abs(np.dot(points, normal) + d)
inliers = distances < distance_threshold
inlier_count = int(np.sum(inliers))
if inlier_count > max_inlier_count:
max_inlier_count = inlier_count
best_inliers = inliers
if max_inlier_count > 0:
return best_inliers
return np.ones(n_points, dtype=bool)
def filter_ransac(
structured_data: dict[str, Any] | np.lib.npyio.NpzFile,
distance_threshold: float = SURFACE_RANSAC_DISTANCE_THRESHOLD_M,
ransac_n: int = SURFACE_RANSAC_N,
num_iterations: int = SURFACE_RANSAC_ITERATIONS,
local_grid_size: float = SURFACE_RANSAC_LOCAL_GRID_SIZE_M,
seed: int = SURFACE_RANSAC_SEED,
progress_callback: Callable[[int], None] | None = None,
) -> np.ndarray:
"""로컬 격자별 RANSAC 평면 분할로 지면 포인트 마스크를 반환한다."""
if not math.isfinite(local_grid_size) or local_grid_size <= 0:
raise ValueError("로컬 격자 크기는 0보다 큰 유한한 값이어야 합니다.")
if ransac_n < 3:
raise ValueError("RANSAC 샘플 수는 3 이상이어야 합니다.")
xyz = np.asarray(structured_data["xyz"], dtype=np.float64)
bounds = np.asarray(structured_data["bounds"], dtype=np.float64)
if xyz.ndim != 2 or xyz.shape[1] != 3:
raise ValueError("xyz 배열은 (N, 3) 형태여야 합니다.")
if bounds.shape != (3, 2):
raise ValueError("bounds 배열은 (3, 2) 형태여야 합니다.")
n_points = len(xyz)
mask = np.zeros(n_points, dtype=bool)
if n_points == 0:
return mask
x_min, y_min = bounds[0, 0], bounds[1, 0]
x_max, y_max = bounds[0, 1], bounds[1, 1]
grid_w = int(np.ceil((x_max - x_min) / local_grid_size)) + 1
grid_h = int(np.ceil((y_max - y_min) / local_grid_size)) + 1
xs, ys = xyz[:, 0], xyz[:, 1]
gx = np.clip(((xs - x_min) / local_grid_size).astype(np.int32), 0, grid_w - 1)
gy = np.clip(((ys - y_min) / local_grid_size).astype(np.int32), 0, grid_h - 1)
grid_indices = gy * grid_w + gx
unique_grids = np.unique(grid_indices)
total_grids = len(unique_grids)
for i, grid_id in enumerate(unique_grids):
cell_points_idx = np.where(grid_indices == grid_id)[0]
if len(cell_points_idx) < ransac_n:
mask[cell_points_idx] = True
continue
cell_inliers = fit_plane_ransac(
xyz[cell_points_idx],
distance_threshold=distance_threshold,
ransac_n=ransac_n,
num_iterations=num_iterations,
seed=seed,
)
mask[cell_points_idx[cell_inliers]] = True
if progress_callback:
progress_callback(int(((i + 1) / total_grids) * 100))
if progress_callback:
progress_callback(100) # 모든 격자 처리 완료 보장
return mask
@@ -1,6 +1,6 @@
"""B04 지면 필터 오케스트레이션.
구조화된 포인트클라우드(structured.npz)에 대해 grid_min_z/csf/pmf/ransac/ptd/classification
구조화된 포인트클라우드(structured.npz)에 대해 ptd(점진 TIN)/classification(업체 분류)
(+ `_refine` 다듬기 판)
필터를 실행하여 지면 마스크 딕셔너리를 만든다. 필터 선택은 config의
source_filters를 따른다.
@@ -16,11 +16,7 @@ from B04_PreProcess.B04_PreProcess_Engine_Filter_Classification import (
has_classified_coverage,
has_classified_ground,
)
from B04_PreProcess.B04_PreProcess_Engine_Filter_CSF import filter_csf
from B04_PreProcess.B04_PreProcess_Engine_Filter_Grid import filter_grid_min_z
from B04_PreProcess.B04_PreProcess_Engine_Filter_PMF import filter_pmf
from B04_PreProcess.B04_PreProcess_Engine_Filter_PTD import filter_ptd
from B04_PreProcess.B04_PreProcess_Engine_Filter_RANSAC import filter_ransac
from B04_PreProcess.B04_PreProcess_Engine_Filter_Refine import refine_ground_mask
from config.config_system import SURFACE_AUTO_FALLBACK_FILTER
@@ -30,10 +26,6 @@ _REFINABLE = ("classification", "ptd")
# 필터 키 → 함수 매핑
_FILTERS = {
"grid_min_z": filter_grid_min_z,
"csf": filter_csf,
"pmf": filter_pmf,
"ransac": filter_ransac,
"ptd": filter_ptd,
"classification": filter_classification,
"classification_refine": lambda data: refine_ground_mask(
@@ -1,7 +1,7 @@
"""B04 지표면 5종 표현 빌더 (TIN/DTM/NURBS/implicit/meshfree).
"""B04 지표면 표현 빌더 (DTM/implicit).
TerrainContext에서 파생 격자·샘플을 받아 각 표현의 모델(npz)과 프리뷰
(GLB/PLY)를 저장하고 메타데이터를 반환한다.
(GLB)를 저장하고 메타데이터를 반환한다.
"""
from contextlib import suppress
@@ -10,7 +10,6 @@ from typing import Any
import numpy as np
from scipy.interpolate import RBFInterpolator
from scipy.spatial import Delaunay
from B04_PreProcess.B04_PreProcess_Engine_ModelContext import (
ProgressCallback,
@@ -18,63 +17,14 @@ from B04_PreProcess.B04_PreProcess_Engine_ModelContext import (
artifact_size,
atomic_npz,
clip_and_compact_mesh,
evaluate_nurbs_spline,
fit_nurbs_spline,
grid_faces,
grid_vertices,
with_footprint,
write_binary_ply,
write_glb,
)
from common_util.common_util_surface_sampler import implicit_grid_npz
def build_tin(
context: TerrainContext, output_dir: Path, stem: str, progress: ProgressCallback
) -> dict[str, Any]:
progress(5)
points = context.sample(int(context.config["tin_max_input_points"]))
unique_xy, unique_indices = np.unique(points[:, :2], axis=0, return_index=True)
points = points[unique_indices]
if len(points) < 3:
raise ValueError("TIN 생성에 필요한 포인트가 부족합니다.")
progress(25)
faces = np.asarray(Delaunay(unique_xy).simplices, dtype=np.uint32)
if len(faces):
triangle_xy = points[faces, :2]
edges = np.stack(
[
np.linalg.norm(triangle_xy[:, 0] - triangle_xy[:, 1], axis=1),
np.linalg.norm(triangle_xy[:, 1] - triangle_xy[:, 2], axis=1),
np.linalg.norm(triangle_xy[:, 2] - triangle_xy[:, 0], axis=1),
],
axis=1,
)
faces = faces[np.max(edges, axis=1) <= float(context.config["tile_size_meters"]) * 2]
valid_vertices = context.contains_xy(points[:, 0], points[:, 1])
points, faces = clip_and_compact_mesh(points, faces, valid_vertices)
if not len(faces):
raise ValueError("외곽 안쪽 기준 적용 후 TIN 면이 남지 않았습니다.")
progress(65)
model_path = output_dir / f"{stem}.npz"
preview_path = output_dir / f"{stem}_preview.glb"
atomic_npz(model_path, vertices=points, faces=faces)
write_glb(preview_path, points, faces, context.bounds)
progress(100)
return with_footprint(
context,
{
"representation": "triangular_mesh",
"model_file": model_path.name,
"preview_file": preview_path.name,
"preview_media_type": "model/gltf-binary",
"vertex_count": int(len(points)),
"face_count": int(len(faces)),
"artifact_bytes": artifact_size(model_path, preview_path),
},
)
def build_dtm(
context: TerrainContext, output_dir: Path, stem: str, progress: ProgressCallback
) -> dict[str, Any]:
@@ -121,56 +71,6 @@ def build_dtm(
)
def build_nurbs(
context: TerrainContext, output_dir: Path, stem: str, progress: ProgressCallback
) -> dict[str, Any]:
degree = max(1, min(5, int(context.config["nurbs_degree"])))
patch_size = float(context.config["nurbs_patch_size_meters"])
controls = max(degree + 1, int(context.config["nurbs_control_points_per_axis"]))
control_resolution = max(patch_size / max(controls - 1, 1), 0.25)
x_control, y_control, z_control = context.grid(control_resolution)
progress(30)
spline, z_range = fit_nurbs_spline(x_control, y_control, z_control, degree)
x_preview, y_preview, _ = context.preview_grid(
float(context.config["dtm_grid_resolution_meters"])
)
z_preview = evaluate_nurbs_spline(spline, y_preview, x_preview, z_range, stem).astype(
np.float32
)
progress(65)
vertices = grid_vertices(x_preview, y_preview, z_preview)
faces = grid_faces(len(y_preview), len(x_preview))
valid_preview = context.contains_xy(vertices[:, 0], vertices[:, 1])
vertices, faces = clip_and_compact_mesh(vertices, faces, valid_preview)
model_path = output_dir / f"{stem}.npz"
preview_path = output_dir / f"{stem}_preview.glb"
atomic_npz(
model_path,
control_x=x_control,
control_y=y_control,
control_z=z_control,
degree=np.array([degree], np.int16),
patch_size_meters=np.array([patch_size], np.float32),
)
write_glb(preview_path, vertices, faces, context.bounds)
progress(100)
return with_footprint(
context,
{
"representation": "bspline_surface",
"model_file": model_path.name,
"preview_file": preview_path.name,
"preview_media_type": "model/gltf-binary",
"degree": degree,
"control_rows": int(len(y_control)),
"control_columns": int(len(x_control)),
"vertex_count": int(len(vertices)),
"face_count": int(len(faces)),
"artifact_bytes": artifact_size(model_path, preview_path),
},
)
def implicit_point_count(ground_count: int, config: dict[str, Any]) -> int:
"""implicit 중심점 수 = clip(지면점 × 비율, 하한, 상한) · 지면점이 하한보다 적으면 전부."""
minimum = int(config["implicit_min_points"])
@@ -240,57 +140,7 @@ def build_implicit(
)
def build_meshfree(
context: TerrainContext, output_dir: Path, stem: str, progress: ProgressCallback
) -> dict[str, Any]:
points = context.sample(int(context.config["meshfree_max_model_points"]))
points = points[context.contains_xy(points[:, 0], points[:, 1])]
if not len(points):
raise ValueError("외곽 안쪽 기준 적용 후 Meshfree 포인트가 남지 않았습니다.")
resolution = float(context.config["dtm_grid_resolution_meters"])
x_grid, y_grid, z_grid = context.grid(resolution)
dz_dy, dz_dx = np.gradient(z_grid, resolution, resolution)
gx = np.clip(np.searchsorted(x_grid, points[:, 0]), 0, len(x_grid) - 1)
gy = np.clip(np.searchsorted(y_grid, points[:, 1]), 0, len(y_grid) - 1)
normals = np.column_stack([-dz_dx[gy, gx], -dz_dy[gy, gx], np.ones(len(points), np.float32)])
normals /= np.maximum(np.linalg.norm(normals, axis=1, keepdims=True), 1e-9)
progress(55)
preview_max = int(context.config["max_preview_vertices"])
if len(points) > preview_max:
selection = np.linspace(0, len(points) - 1, preview_max, dtype=np.int64)
preview_points, preview_normals = points[selection], normals[selection]
else:
preview_points, preview_normals = points, normals
model_path = output_dir / f"{stem}.npz"
preview_path = output_dir / f"{stem}_preview.ply"
radius = float(context.config["meshfree_point_radius_meters"])
atomic_npz(
model_path,
points=points,
normals=normals.astype(np.float32),
radius=np.array([radius], np.float32),
)
write_binary_ply(preview_path, preview_points, preview_normals, context.bounds)
progress(100)
return with_footprint(
context,
{
"representation": "meshfree_surfels",
"model_file": model_path.name,
"preview_file": preview_path.name,
"preview_media_type": "application/octet-stream",
"point_count": int(len(points)),
"preview_point_count": int(len(preview_points)),
"point_radius_meters": radius,
"artifact_bytes": artifact_size(model_path, preview_path),
},
)
BUILDERS = {
"tin": build_tin,
"dtm": build_dtm,
"nurbs": build_nurbs,
"implicit": build_implicit,
"meshfree": build_meshfree,
}
@@ -1,8 +1,8 @@
"""B04 지표면 모델 공통 컨텍스트 및 메시 유틸리티.
지면 마스크가 적용된 포인트에서 footprint(외곽), 격자, 프리뷰 격자를 만들고,
GLB/PLY 프리뷰 및 npz 모델을 원자적으로 저장하는 공통 기능을 제공한다.
5개 표현(TIN/DTM/NURBS/implicit/meshfree) 빌더가 이 컨텍스트를 공유한다.
GLB 프리뷰 및 npz 모델을 원자적으로 저장하는 공통 기능을 제공한다.
표현(DTM/implicit) 빌더가 이 컨텍스트를 공유한다.
"""
import hashlib
@@ -16,7 +16,6 @@ from typing import Any, Callable
import numpy as np
import trimesh
from scipy import ndimage
from scipy.interpolate import RectBivariateSpline
from skimage.restoration import inpaint_biharmonic
from common_util.common_util_atomic import atomic_write_bytes, atomic_write_npz
@@ -24,8 +23,6 @@ from common_util.common_util_atomic import atomic_write_bytes, atomic_write_npz
logger = logging.getLogger(__name__)
MODEL_VERSION = 1
MODEL_METHODS = ("tin", "dtm", "nurbs", "implicit", "meshfree")
SOURCE_FILTER_LABELS = {"grid_min_z": "Grid Min-Z", "csf": "CSF", "pmf": "PMF"}
ProgressCallback = Callable[[int], None]
# 대용량 포인트 배치 처리 크기
@@ -33,63 +30,6 @@ _BATCH_SIZE = 500_000
# implicit 격자 칸 뽑기 씨앗 — 고정이라 다시 지어도 같은 중심점
IMPLICIT_SAMPLE_SEED = 558
# NURBS 평활 계수 — 제어점당 허용 잔차제곱(m²).
# 0.01(RMS 0.1m)은 요구가 너무 빡빡해 FITPACK이 제어점 수(122x113)보다 많은
# knot(126x117)을 밀어넣고 계수가 발산했다 — csf 표고 1e105, grid_min_z 1e10
# (2026-08-17 실측). 1.0(RMS 1m)이면 knot이 17~35로 떨어지고 세 지면 필터 모두
# 데이터 표고 범위 안에 머문다. NURBS는 평활 곡면 표현이라 이 정도 완화가 맞다.
NURBS_RESIDUAL_PER_CONTROL_POINT = 1.0
def fit_nurbs_spline(
x_control: np.ndarray, y_control: np.ndarray, z_control: np.ndarray, degree: int
) -> tuple[RectBivariateSpline, tuple[float, float]]:
"""제어 격자에 B-spline 곡면을 맞추고, 발산 판정용 표고 허용 범위를 함께 준다.
모델 빌더와 등고선 엔진이 같은 곡면을 각각 만들므로 여기서 한 번에 정의한다.
"""
spline = RectBivariateSpline(
y_control,
x_control,
z_control,
kx=min(degree, len(y_control) - 1),
ky=min(degree, len(x_control) - 1),
s=float(len(x_control) * len(y_control)) * NURBS_RESIDUAL_PER_CONTROL_POINT,
)
z_min = float(np.min(z_control))
z_max = float(np.max(z_control))
margin = max((z_max - z_min) * 0.5, 1.0)
return spline, (z_min - margin, z_max + margin)
def evaluate_nurbs_spline(
spline: RectBivariateSpline,
y_coords: np.ndarray,
x_coords: np.ndarray,
z_range: tuple[float, float],
label: str = "",
) -> np.ndarray:
"""스플라인을 평가하고 표고 허용 범위를 벗어난 값을 잘라낸다.
평활 스플라인은 제어점이 성긴 구석에서 발산할 수 있다. 그대로 두면 float32
캐스팅에서 inf가 되고 프리뷰 색상(nan)과 등고선 레벨 산출(`np.arange`가
"Maximum allowed size exceeded"로 실패)까지 번진다.
"""
z_values = np.asarray(spline(y_coords, x_coords), dtype=np.float64)
low, high = z_range
outside = np.count_nonzero(~np.isfinite(z_values) | (z_values < low) | (z_values > high))
if outside:
logger.warning(
"NURBS 곡면이 표고 범위를 벗어나 잘라냈습니다: %s%d개 셀 (허용 %.1f~%.1fm)",
f"{label} " if label else "",
int(outside),
low,
high,
)
z_values = np.nan_to_num(z_values, nan=low, posinf=high, neginf=low)
return np.clip(z_values, low, high)
# 서명에서 빼는 키 — 지오메트리를 만드는 방법이 아니라 "무엇을 만들지"를 고르는 값들.
# source_filters·precompute를 넣으면 조합을 바꿀 때마다 manifest가 통째로 폐기돼
# 이전에 만들어 둔 조합이 사라진다. manifest는 이미 필터·방식별로 키가 갈린다.
@@ -155,44 +95,6 @@ def write_glb(path: Path, vertices: np.ndarray, faces: np.ndarray, bounds: np.nd
atomic_write_bytes(path, payload)
def write_binary_ply(
path: Path, vertices: np.ndarray, normals: np.ndarray, bounds: np.ndarray
) -> None:
verts = scene_vertices(vertices, bounds)
scene_normals = np.empty_like(normals, dtype=np.float32)
scene_normals[:, 0] = normals[:, 0]
scene_normals[:, 1] = normals[:, 2]
scene_normals[:, 2] = -normals[:, 1]
colors = height_colors(vertices)
dtype = np.dtype(
[
("x", "<f4"),
("y", "<f4"),
("z", "<f4"),
("nx", "<f4"),
("ny", "<f4"),
("nz", "<f4"),
("red", "u1"),
("green", "u1"),
("blue", "u1"),
("alpha", "u1"),
]
)
records = np.empty(len(vertices), dtype=dtype)
records["x"], records["y"], records["z"] = verts.T
records["nx"], records["ny"], records["nz"] = scene_normals.T
records["red"], records["green"], records["blue"], records["alpha"] = colors.T
header = (
"ply\nformat binary_little_endian 1.0\n"
f"element vertex {len(vertices)}\n"
"property float x\nproperty float y\nproperty float z\n"
"property float nx\nproperty float ny\nproperty float nz\n"
"property uchar red\nproperty uchar green\nproperty uchar blue\nproperty uchar alpha\n"
"end_header\n"
).encode("ascii")
atomic_write_bytes(path, header + records.tobytes())
def grid_faces(rows: int, cols: int) -> np.ndarray:
"""정규 격자의 삼각형 면 인덱스를 만든다."""
if rows < 2 or cols < 2:
@@ -1,6 +1,6 @@
"""B04 지표면 모델 파이프라인 오케스트레이터.
세 지면 필터(grid_min_z/csf/pmf)와 다섯 표현(TIN/DTM/NURBS/implicit/meshfree)의
지면 필터(업체 분류 · 점진 TIN · 둘의 다듬기 판)와 표현(DTM/implicit)의
캐시를 만들고 manifest.json을 관리한다. 캐시 유효성 검증, 스무딩/등고선 연동,
동일 출력 폴더의 중복 실행 취소를 포함한다.
"""
@@ -88,10 +88,9 @@ def _cache_contours(
smooth_model_path = output_dir / f"{stem}_smooth.npz"
smooth_meta = metadata.get("smooth", {})
if smooth_model_path.exists() and smooth_meta.get("status") == "completed":
smooth_rep = "regular_grid" if method == "dtm" else "triangular_mesh"
smooth_contours = extract_contours(
smooth_model_path,
representation=smooth_rep,
representation="regular_grid",
interval=interval,
target_grid_m=target_grid_m,
scene_center=None,
@@ -143,10 +142,7 @@ def _cache_contours_safely(
_REPRESENTATIONS = {
"meshfree": "meshfree_surfels",
"dtm": "regular_grid",
"tin": "triangular_mesh",
"nurbs": "bspline_surface",
"implicit": "local_rbf_height_field",
}
@@ -157,17 +153,16 @@ def _cache_is_valid(
"""디스크의 결과 파일과 스무딩 메타데이터가 유효한지 검사한다.
manifest 항목이 완료여야 함 — 설정 지문이 바뀌면 manifest 가 비므로 파일만 남은
옛 결과(평활 없는 implicit · nurbs · meshfree)를 다시 짓는다(요구 558).
옛 결과(평활 없는 implicit)를 다시 짓는다(요구 558).
"""
if entry.get("status") != "completed":
return False
ext = "ply" if method == "meshfree" else "glb"
files_exist = (output_dir / f"{stem}_preview.{ext}").exists() and (
files_exist = (output_dir / f"{stem}_preview.glb").exists() and (
output_dir / f"{stem}.npz"
).exists()
if not files_exist:
return False
if method in config.get("smoothing_methods", ("dtm", "tin")):
if method in config.get("smoothing_methods", ("dtm",)):
smooth_entry = entry.get("smooth", {})
smooth_exist = (output_dir / f"{stem}_smooth.npz").exists() and (
output_dir / f"{stem}_smooth_preview.glb"
@@ -190,14 +185,14 @@ def _build_all_terrain_models(
force: bool = False,
progress: ProgressReporter | None = None,
) -> dict[str, Any]:
"""세 지면 필터와 다섯 표현의 캐시를 만들고 manifest를 반환한다."""
"""지면 필터와 표현의 캐시를 만들고 manifest를 반환한다."""
output_dir.mkdir(parents=True, exist_ok=True)
manifest_path = output_dir / "manifest.json"
filters = tuple(key for key in config["source_filters"] if key in ground_masks)
# dtm 은 늘 맨 앞 — 같은 필터 외곽 · implicit 굳힌 격자 · 등고선 자르기의 바탕이라
# 어느 표현을 지어도 낡았으면(manifest 항목 없음 · 지문 다름) 같이 다시 지음(754-1).
methods = ("dtm", *(key for key in config["precompute"] if key in BUILDERS and key != "dtm"))
# 무거운 implicit 은 기본 필터(첫 필터)에만 — 분류 LAS 의 csf 까지 두 번 짓지 않음(요구 558).
# 무거운 implicit 은 기본 필터(첫 필터)에만 — 둘째 필터엔 안 지음(요구 558).
later_methods = tuple(method for method in methods if method != "implicit")
signature = config_signature(config)
bounds = np.asarray(structured_data["bounds"], dtype=np.float64)
@@ -248,7 +243,7 @@ def _build_all_terrain_models(
interval = float(config.get("contour_interval_meters", 1.0))
contour_path = output_dir / f"contour_{filter_key}_{method}_{interval}m.json"
needs_contours = contour_cache_is_stale(contour_path, output_dir / f"{stem}.npz")
if method in config.get("smoothing_methods", ("dtm", "tin")):
if method in config.get("smoothing_methods", ("dtm",)):
smooth_contour_path = (
output_dir / f"contour_{filter_key}_{method}_smooth_{interval}m.json"
)
@@ -302,7 +297,7 @@ def _build_all_terrain_models(
"error": None,
}
)
if method in config.get("smoothing_methods", ("dtm", "tin")):
if method in config.get("smoothing_methods", ("dtm",)):
original_model_path = output_dir / f"{stem}.npz"
if original_model_path.exists():
try:
@@ -1,57 +1,23 @@
"""도엽등고선 → 표고 격자 보간 방식 모음 (비교용).
"""도엽등고선 → 표고 격자 보간 — 라플라스(거리 비례 초기값 + 이웃 평균 다듬기).
문헌(Hutchinson 1988/89 ANUDEM; Chaplot 2006; Arun 2013 등)은 지형 복잡도·자료 밀도에
따라 우열이 갈리며 단일 최적해가 없다고 본다. 그래서 방식을 하나로 고르지 않고 여기
모아 두고 B04 화면에서 바꿔 가며 보게 한다(2026-08-30 사용자 지시).
비교용 다섯 방식(TIN 둘 · TPS · ANUDEM형 · 다중해상도)은 걷음(608) — 확정 · 자동 = 라플라스.
각 builder는 `(spec, burned, features, cell_m) -> (R, C) float32` 격자를 돌려준다.
builder는 `(spec, burned, features, cell_m) -> (R, C) float32` 격자를 돌려준다.
`burned`는 등고 라인이 구워진 격자(라인 셀 = 표고, 그 외 NaN)다. 폐합 링 안쪽 처리와
프리뷰·저장은 호출측(`_SheetSurface`)이 방식과 무관하게 똑같이 해 준다.
"""
import logging
import warnings
from typing import Any, Callable
import numpy as np
logger = logging.getLogger(__name__)
# TIN(도엽선)이 물어 오는 격자 밖 여유 — 테두리가 삼각망 밖으로 나가지 않을 만큼만.
SHEET_TIN_CLIP_MARGIN_M = 100.0
# 화면 버튼에 쓰는 이름 — 키는 surface_models.generation_params.source_filter 접미사다.
SHEET_METHOD_LABELS: dict[str, str] = {
"tin_sheet": "TIN(도엽선)",
"tin": "TIN 격자",
"biharmonic": "TPS(박판)",
"anudem": "ANUDEM형",
"multires": "다중해상도",
"laplace": "라플라스",
}
def _laplacian(values: np.ndarray) -> np.ndarray:
padded = np.pad(values, 1, mode="edge")
return (
padded[:-2, 1:-1] + padded[2:, 1:-1] + padded[1:-1, :-2] + padded[1:-1, 2:] - 4.0 * values
)
def _contour_vertices(burned: np.ndarray, spec: Any) -> tuple[np.ndarray, np.ndarray]:
"""등고 라인 셀을 (N,2) 세계좌표와 표고로 바꾼다."""
rows, cols = np.nonzero(np.isfinite(burned))
xs = spec.cell_centers_x()[cols]
ys = spec.cell_centers_y()[rows]
return np.column_stack([xs, ys]), burned[rows, cols].astype(np.float64)
def _grid_points(spec: Any) -> tuple[np.ndarray, np.ndarray]:
grid_x, grid_y = np.meshgrid(spec.cell_centers_x(), spec.cell_centers_y())
return grid_x, grid_y
# ── ① 거리 비례 (버튼에서는 뺐지만 다른 방식의 초기추정으로 계속 쓴다) ───────
# ── 거리 비례 (라플라스 초기값) ─────────────────────────────────────────────
def build_distance(spec: Any, burned: np.ndarray, features: Any, cell_m: float) -> np.ndarray:
"""가장 가까운 서로 다른 표고 두 라인 사이를 거리 비례로 나눈다.
@@ -120,394 +86,17 @@ def relax_laplace(surface: np.ndarray, fixed: np.ndarray, iterations: int) -> No
surface[colour] = neighbours[colour]
# ── ② 라플라스(조화) ─────────────────────────────────────────────────────────
# ── 라플라스(조화) ─────────────────────────────────────────────────────────────
def build_laplace(spec: Any, burned: np.ndarray, features: Any, cell_m: float) -> np.ndarray:
"""등고선을 경계값으로 두고 Δz=0을 푼다.
면이 매끈해지지만 z=r(원뿔·능선)은 조화함수가 아니라 마루가 눌린다. 비교 기준으로
남겨 둔다 — ANUDEM이 라플라스 대신 박판 스플라인을 쓰는 이유를 눈으로 보기 위함.
면이 매끈해지지만 z=r(원뿔·능선)은 조화함수가 아니라 마루가 조금 눌린다.
"""
surface = build_distance(spec, burned, features, cell_m)
relax_laplace(surface, np.isfinite(burned), 400)
return surface
# ── ③ 박판 스플라인(중조화) ──────────────────────────────────────────────────
def solve_min_curvature(constrained: np.ndarray, guess: np.ndarray) -> np.ndarray:
"""제약 셀을 고정하고 Δ²z=0(박판 스플라인)을 최소곡률 최소제곱으로 푼다.
ANUDEM/Topo to Raster가 쓰는 박판 스플라인과 같은 연산자다. 라플라스와 달리 z=r을
그대로 통과시켜 능선·마루가 눌리지 않고, 경사가 등고선 너머로 자연스럽게 이어진다.
`constrained`: 값이 고정된 셀(등고 라인, 필요하면 구조선 앵커) — 그 외는 NaN.
`guess`: 시작값이자 감쇠 기준(보통 거리 보간 결과).
"""
from scipy.sparse.linalg import LinearOperator, lsmr
guess = guess.astype(np.float64)
burned = constrained
# 제약 셀 + **격자 테두리**를 고정한다. 최외곽 등고선 바깥이 통째로 자유면
# 1차함수가 Δ²의 영공간에 남아 해가 하나로 정해지지 않고 켤레기울기가 발산한다
# (2026-08-30 실측: |Δz| 2092m). 테두리는 거리 보간값으로 묶는다.
fixed = np.isfinite(burned)
fixed[0, :] = fixed[-1, :] = True
fixed[:, 0] = fixed[:, -1] = True
free = ~fixed & np.isfinite(guess)
if not free.any():
return guess.astype(np.float32)
index = np.flatnonzero(free.ravel())
values = np.where(np.isfinite(burned), burned, guess)
base = np.where(fixed, np.nan_to_num(values), 0.0)
# Δ²z=0을 정규방정식(CG)으로 풀면 조건수가 격자변 4제곱이라 발산한다(실측).
# 대신 **최소곡률** 최소제곱으로 세운다 — 자유 셀에 대해 ‖Δz‖를 최소화하며,
# 그 정상해가 곧 Δ²z=0이다. 조건수가 제곱으로 줄어 LSMR이 안정적으로 푼다
# (Briggs 1974의 최소곡률 격자화와 같은 목적함수).
# 최소곡률만으로는 제약(등고선)에서 먼 영역이 정해지지 않아 해가 폭주한다.
# ANUDEM의 거칠기 벌점과 같은 취지로 감쇠항을 붙여 거리 보간값에 묶어 둔다:
# minimize ‖Δz‖² + λ‖z − 거리보간‖²
# λ가 작을수록 더 매끈하고 클수록 거리 보간에 가깝다.
total_cells = base.size
free_count = len(index)
damping = np.float64(np.sqrt(0.02))
anchor = guess.ravel()[index]
def forward(vector: np.ndarray) -> np.ndarray:
# 반드시 **선형**이어야 한다 — 고정 셀 기여(base)를 여기서 더하면 아핀이 되어
# LSMR의 전제가 깨지고 해가 폭주한다. base 몫은 우변으로만 넘긴다.
scattered = np.zeros_like(base)
scattered.ravel()[index] = vector
return np.concatenate([_laplacian(scattered).ravel(), damping * vector])
def adjoint(vector: np.ndarray) -> np.ndarray:
curvature = _laplacian(vector[:total_cells].reshape(base.shape)).ravel()[index]
return curvature + damping * vector[total_cells:]
rhs = np.concatenate([-_laplacian(base).ravel(), damping * anchor])
linear = LinearOperator(
(total_cells + free_count, free_count),
matvec=forward,
rmatvec=adjoint,
dtype=np.float64,
)
result = lsmr(linear, rhs, x0=anchor, maxiter=400, atol=1e-8, btol=1e-8)
solution, info = result[0], result[1]
surface = base.copy()
surface.ravel()[index] = solution
surface[fixed] = values[fixed]
# 안전장치 — 발산하면 조용히 틀린 지형을 넘기지 말고 거리 보간으로 되돌린다.
drift = float(np.nanmax(np.abs(surface - guess)))
span = float(np.nanmax(guess) - np.nanmin(guess))
if not np.isfinite(drift) or drift > max(span, 1.0):
logger.warning(
"도엽 서피스(TPS): 해가 발산해(최대 %.1fm) 거리 보간으로 되돌립니다 (info=%s).",
drift,
info,
)
return guess.astype(np.float32)
logger.info("도엽 서피스(TPS): 최소제곱 info=%s, 최대 변화 %.2fm", info, drift)
return surface.astype(np.float32)
def build_biharmonic(spec: Any, burned: np.ndarray, features: Any, cell_m: float) -> np.ndarray:
"""등고선만 제약으로 둔 박판 스플라인."""
return solve_min_curvature(burned, build_distance(spec, burned, features, cell_m))
# ── ④ TIN ───────────────────────────────────────────────────────────────────
def _triangulate(
spec: Any, shape_: tuple[int, int], points: np.ndarray, values: np.ndarray
) -> np.ndarray:
"""정점 구름을 Delaunay 삼각망 선형 보간해 격자로 편다. 삼각망 밖은 NaN."""
from scipy.interpolate import LinearNDInterpolator
if len(points) < 3:
return np.full(shape_, np.nan, dtype=np.float32)
interpolator = LinearNDInterpolator(points, values)
grid_x, grid_y = _grid_points(spec)
surface = np.empty(shape_, dtype=np.float32)
chunk = max(1, int(4_000_000 // max(spec.n_cols, 1)))
for start in range(0, spec.n_rows, chunk):
stop = min(start + chunk, spec.n_rows)
surface[start:stop] = interpolator(grid_x[start:stop], grid_y[start:stop]).astype(
np.float32
)
return surface
def build_tin(spec: Any, burned: np.ndarray, features: Any, cell_m: float) -> np.ndarray:
"""격자에 구운 등고 라인 셀을 Delaunay 삼각망 선형 보간한다.
정점이 셀 중심에 맞춰져 있어 1m 계단이 삼각망에 그대로 실린다. 같은 표고 정점
3개로 이루어진 평탄 삼각형이 굴곡부·마루에 계단을 만든다. 비교 기준으로 남긴다.
"""
points, values = _contour_vertices(burned, spec)
if len(points) > 120_000: # 삼각망 비용은 정점 수에 비례한다
step = int(np.ceil(len(points) / 120_000))
points, values = points[::step], values[::step]
return _triangulate(spec, burned.shape, points, values)
def build_tin_sheet(spec: Any, burned: np.ndarray, features: Any, cell_m: float) -> np.ndarray:
"""**원본 5m 도엽 등고선 정점**을 그대로 이은 고전 TIN (2026-08-30 사용자 지시).
`tin`은 격자에 구운 라인 셀(=1m 계단으로 뭉개진 정점)을 쓰지만, 이쪽은 벡터
등고선의 정점을 좌표 그대로 쓴다 — 도면 등고선을 삼각망으로 잇는 측량 관행 그대로다.
격자 밖 등고선은 물지 않는다(삼각망 비용만 커지고 결과는 같다). 다만 격자 테두리가
삼각망 밖으로 나가지 않도록 여유를 두고 자른다.
길이 필터는 두지 않는다 — 짧은 봉우리 폐합 등고선을 버리면 마루가 통째로 평평해진다.
"""
from shapely.geometry import shape as to_shape
from B04_PreProcess.B04_PreProcess_Engine_Watershed_Grid import (
ELEVATION_KEYS,
iter_linestrings,
)
xs, ys = spec.cell_centers_x(), spec.cell_centers_y()
margin = max(SHEET_TIN_CLIP_MARGIN_M, cell_m * 2.0)
x_lo, x_hi = xs[0] - margin, xs[-1] + margin
y_lo, y_hi = ys[-1] - margin, ys[0] + margin
coords: list[np.ndarray] = []
levels: list[np.ndarray] = []
for feature in features or []:
properties = feature.get("properties") or {}
elevation = next(
(float(properties[key]) for key in ELEVATION_KEYS if properties.get(key) is not None),
None,
)
geometry = feature.get("geometry")
if elevation is None or not geometry:
continue
try:
parsed = to_shape(geometry)
except Exception: # noqa: BLE001 — 손상된 피처는 건너뛴다
continue
for line in iter_linestrings(parsed):
point = np.asarray(line.coords, dtype=np.float64)[:, :2]
inside = (
(point[:, 0] >= x_lo)
& (point[:, 0] <= x_hi)
& (point[:, 1] >= y_lo)
& (point[:, 1] <= y_hi)
)
if not inside.any():
continue
coords.append(point[inside])
levels.append(np.full(int(inside.sum()), elevation, dtype=np.float64))
if not coords:
logger.warning("도엽 서피스: TIN(도엽선)에 쓸 등고선 정점이 없습니다.")
return np.full(burned.shape, np.nan, dtype=np.float32)
points = np.vstack(coords)
values = np.concatenate(levels)
# 도엽 이음매에서 같은 정점이 겹쳐 들어온다 — Qhull 비용만 늘어 미리 접는다.
_, unique = np.unique(np.round(points, 3), axis=0, return_index=True)
points, values = points[unique], values[unique]
logger.info("도엽 서피스: TIN(도엽선) 정점 %d개", len(points))
return _triangulate(spec, burned.shape, points, values)
# ── ⑦ ANUDEM형 (구조선 + 배수 강제) ─────────────────────────────────────────
def _contour_corner_anchors(
spec: Any, burned: np.ndarray, guess: np.ndarray, cell_m: float
) -> tuple[np.ndarray, np.ndarray] | None:
"""등고선의 국소 최대 곡률점(코너)에서 능선·계곡 구조선 앵커를 만든다.
ANUDEM은 등고선 자체의 곡률에서 능선·계곡망을 먼저 뽑아 흐름 구조를 세운다
(Hutchinson 1988/89). 여기서도 같은 순서를 따른다.
① 라인마다 정점 곡률을 재 국소 최대점(V자 꼭짓점)을 고른다
② 굽은 안쪽이 더 높으면 **계곡**(등고선 V가 상류를 가리킴), 낮으면 **능선**
③ 같은 종류의 코너를 이웃 표고끼리 이어 그 사이를 선형 보간해 앵커로 심는다
앵커는 박판 해의 제약으로 들어가 계곡 바닥이 이어져 내려가고 능선 마루가 선다.
"""
from scipy.ndimage import label
from scipy.spatial import cKDTree
levels = np.unique(burned[np.isfinite(burned)])
if len(levels) < 2:
return None
interval = float(np.diff(levels).min())
xs = spec.cell_centers_x()
ys = spec.cell_centers_y()
corners: list[tuple[float, float, float, int]] = [] # x, y, level, +1 계곡 / -1 능선
span = 6 # 곡률을 재는 정점 간격(px) — 짧으면 노이즈, 길면 꼭짓점을 놓친다
for level in levels:
labelled, count = label(burned == level)
for component_id in range(1, count + 1):
line_rows, line_cols = np.nonzero(labelled == component_id)
if len(line_rows) < 3 * span:
continue
# 라인 셀을 한 줄로 세운다 — 좌표 정렬로 근사한다(정밀 추적은 과하다).
order = np.argsort(line_cols + line_rows * 1e-3)
path = np.column_stack([line_cols[order], line_rows[order]]).astype(np.float64)
before = np.roll(path, span, axis=0)
after = np.roll(path, -span, axis=0)
first = path - before
second = after - path
first_len = np.hypot(first[:, 0], first[:, 1])
second_len = np.hypot(second[:, 0], second[:, 1])
valid = (first_len > 1e-6) & (second_len > 1e-6)
cosine = np.ones(len(path))
cosine[valid] = (first[valid] * second[valid]).sum(axis=1) / (
first_len[valid] * second_len[valid]
)
sharp = np.flatnonzero(valid & (cosine < 0.3)) # 70도 이상 꺾인 자리
for i in sharp[:: max(1, span)]:
# 굽은 안쪽 방향 = 두 변 단위벡터 합의 반대
inward = -(first[i] / first_len[i] + second[i] / second_len[i])
norm = float(np.hypot(inward[0], inward[1]))
if norm < 1e-6:
continue
probe = path[i] + inward / norm * 6.0
probe_col = int(round(probe[0]))
probe_row = int(round(probe[1]))
if not (0 <= probe_row < guess.shape[0] and 0 <= probe_col < guess.shape[1]):
continue
inside = float(guess[probe_row, probe_col])
if not np.isfinite(inside) or abs(inside - level) < interval * 0.15:
continue
corners.append(
(
float(xs[int(path[i, 0])]),
float(ys[int(path[i, 1])]),
float(level),
1 if inside > level else -1,
)
)
if len(corners) < 4:
logger.info("도엽 서피스(ANUDEM형): 등고선 코너가 부족해 구조선을 건너뜁니다.")
return None
array = np.asarray(corners, dtype=np.float64)
anchor_xy: list[np.ndarray] = []
anchor_z: list[np.ndarray] = []
reach = interval * 20.0 # 이보다 먼 코너는 같은 구조선으로 보지 않는다
for level in levels[:-1]:
upper = level + interval
lower_set = array[np.abs(array[:, 2] - level) < 1e-6]
upper_set = array[np.abs(array[:, 2] - upper) < 1e-6]
if not len(lower_set) or not len(upper_set):
continue
tree = cKDTree(upper_set[:, :2])
distance, index = tree.query(lower_set[:, :2], k=1)
for i in range(len(lower_set)):
j = int(index[i])
if distance[i] > reach or lower_set[i, 3] != upper_set[j, 3]:
continue
start = lower_set[i, :2]
end = upper_set[j, :2]
steps = max(2, int(distance[i] / max(cell_m, 1e-6) / 4))
fraction = np.linspace(0.0, 1.0, steps + 1)[1:-1]
if not len(fraction):
continue
anchor_xy.append(start + (end - start) * fraction[:, None])
anchor_z.append(level + interval * fraction)
if not anchor_xy:
return None
return np.vstack(anchor_xy), np.concatenate(anchor_z)
def _enforce_drainage(surface: np.ndarray, epsilon: float = 0.01) -> int:
"""가짜 웅덩이를 메운다 — ANUDEM의 배수 강제와 같은 목적.
등고선만으로 만든 면에는 흐름이 끊기는 웅덩이가 남는다. 형태학적 재구성(erosion)
으로 채우되 완전 평탄해지지 않게 아주 작은 값을 얹는다. 채운 셀 수를 반환한다.
"""
from skimage.morphology import reconstruction
if not np.isfinite(surface).all():
return 0
seed = np.full(surface.shape, float(surface.max()), dtype=np.float64)
seed[0, :] = surface[0, :]
seed[-1, :] = surface[-1, :]
seed[:, 0] = surface[:, 0]
seed[:, -1] = surface[:, -1]
filled = reconstruction(seed, surface.astype(np.float64), method="erosion")
raised = filled > surface + 1e-6
if not raised.any():
return 0
# 그냥 채우면 웅덩이가 통째로 평탄해져 흐름 방향이 없어진다. 채운 영역 안쪽으로
# 갈수록 아주 조금 높아지게 해서 물이 가장자리(넘침점)로 빠져나가게 둔다
# (Garbrecht·Martz의 평탄면 해소를 간단히 옮긴 것 — 표고 변화는 cm 단위다).
from scipy.ndimage import distance_transform_edt
inner = distance_transform_edt(raised)
surface[raised] = (filled[raised] + epsilon * inner[raised]).astype(surface.dtype)
return int(raised.sum())
def build_anudem(spec: Any, burned: np.ndarray, features: Any, cell_m: float) -> np.ndarray:
"""ANUDEM형 — 등고선 곡률에서 능선·계곡 구조선을 뽑아 제약에 더하고, 박판으로 풀고,
가짜 웅덩이를 메운다. Topo to Raster가 밟는 세 단계를 그대로 옮긴 것이다."""
guess = build_distance(spec, burned, features, cell_m).astype(np.float64)
constrained = burned.astype(np.float64).copy()
anchors = _contour_corner_anchors(spec, burned, guess, cell_m)
if anchors is not None:
xy, z = anchors
row, col = spec.world_to_rc(xy[:, 0].copy(), xy[:, 1].copy())
inside = (row >= 0) & (col >= 0)
row, col, z = row[inside], col[inside], z[inside]
free = ~np.isfinite(constrained[row, col])
constrained[row[free], col[free]] = z[free]
logger.info("도엽 서피스(ANUDEM형): 구조선 앵커 %d셀", int(free.sum()))
surface = solve_min_curvature(constrained, guess).astype(np.float32)
logger.info("도엽 서피스(ANUDEM형): 가짜 웅덩이 %d셀 메움", _enforce_drainage(surface))
return surface
# ── ⑧ 다중해상도 (coarse → fine) ────────────────────────────────────────────
def build_multires(spec: Any, burned: np.ndarray, features: Any, cell_m: float) -> np.ndarray:
"""성긴 격자에서 풀고 점차 세밀화한다 — ANUDEM의 다중해상도 전략.
전체 형상은 성긴 격자에서 싸게 잡고, 세밀한 격자에서는 등고선 근처만 다듬는다.
한 해상도에서만 풀 때보다 넓은 밴드가 고르게 퍼지고 값싸게 수렴한다.
"""
from scipy.ndimage import zoom
surface: np.ndarray | None = None
for factor in (8, 4, 2, 1):
if factor == 1:
coarse = burned
else:
# 성긴 격자의 제약 — 블록 안 등고 라인의 평균 표고를 대표로 쓴다.
rows = burned.shape[0] // factor * factor
cols = burned.shape[1] // factor * factor
blocks = burned[:rows, :cols].reshape(rows // factor, factor, cols // factor, factor)
# 라인이 하나도 없는 블록은 NaN이 정상이라 경고를 삼킨다.
with warnings.catch_warnings():
warnings.simplefilter("ignore", RuntimeWarning)
coarse = np.nanmean(blocks, axis=(1, 3)).astype(np.float32)
level = build_distance(spec, coarse, features, cell_m * factor)
if surface is not None:
# 앞 단계 해를 지금 해상도로 올려 절반씩 섞는다 — 성긴 단계의 넓은 추세를
# 이어받되 이번 해상도의 등고선 정보를 덮지 않는다.
scale = (level.shape[0] / surface.shape[0], level.shape[1] / surface.shape[1])
upscaled = zoom(surface, scale, order=1)
free = ~np.isfinite(coarse)
level[free] = (level[free] + upscaled[free]) * 0.5
relax_laplace(level, np.isfinite(coarse), 8)
surface = level
assert surface is not None
if surface.shape != burned.shape: # 블록 자르기로 남은 가장자리 보정
scale = (burned.shape[0] / surface.shape[0], burned.shape[1] / surface.shape[1])
surface = zoom(surface, scale, order=1)
line = np.isfinite(burned)
surface[line] = burned[line]
return surface.astype(np.float32)
SHEET_METHOD_BUILDERS: dict[str, Callable[[Any, np.ndarray, Any, float], np.ndarray]] = {
"tin_sheet": build_tin_sheet,
"tin": build_tin,
"biharmonic": build_biharmonic,
"anudem": build_anudem,
"multires": build_multires,
"laplace": build_laplace,
}
+2 -41
View File
@@ -1,7 +1,7 @@
"""B04 지표면 스무딩 엔진 (DTM 가우시안+B-Spline, TIN Taubin).
"""B04 지표면 스무딩 엔진 (DTM 가우시안+B-Spline).
원본 모델(npz)을 로드해 표고 Z만 평활화한 스무딩 모델과 프리뷰(GLB)를 만든다.
XY 위치와 삼각형 연결은 원본을 유지한다.
XY 격자는 원본을 유지한다.
"""
import hashlib
@@ -10,7 +10,6 @@ from pathlib import Path
from typing import Any
import numpy as np
import trimesh
from scipy import ndimage
from scipy.interpolate import RectBivariateSpline
@@ -104,48 +103,10 @@ def smooth_dtm(
}
def smooth_tin(
context: TerrainContext, output_dir: Path, stem: str, original_model_path: Path
) -> dict[str, Any]:
"""TIN 삼각망에 Taubin 저수축 스무딩을 적용한다 (Z만 평활화)."""
data = np.load(original_model_path)
vertices, faces = data["vertices"], data["faces"]
if len(vertices) < 3 or not len(faces):
raise ValueError("TIN 스무딩을 수행할 삼각망 메쉬 데이터가 올바르지 않습니다.")
iterations = int(context.config.get("smoothing_tin_taubin_iterations", 10))
lamb = float(context.config.get("smoothing_tin_taubin_lambda", 0.5))
mu = float(context.config.get("smoothing_tin_taubin_mu", -0.53))
mesh = trimesh.Trimesh(vertices=vertices, faces=faces, process=False)
if iterations > 0:
trimesh.smoothing.filter_taubin(mesh, lamb=lamb, nu=mu, iterations=iterations)
vertices_smooth = np.asarray(mesh.vertices, dtype=np.float32)
# 지표면 스무딩은 리토폴로지가 아니다: XY·연결은 원본, Z만 평활화.
vertices_smooth[:, :2] = np.asarray(vertices[:, :2], dtype=np.float32)
faces_smooth = np.asarray(mesh.faces, dtype=np.uint32)
model_path = output_dir / f"{stem}_smooth.npz"
preview_path = output_dir / f"{stem}_smooth_preview.glb"
atomic_npz(model_path, vertices=vertices_smooth, faces=faces_smooth)
write_glb(preview_path, vertices_smooth, faces_smooth, context.bounds)
return {
"model_file": model_path.name,
"preview_file": preview_path.name,
"preview_media_type": "model/gltf-binary",
"vertex_count": int(len(vertices_smooth)),
"face_count": int(len(faces_smooth)),
"artifact_bytes": artifact_size(model_path, preview_path),
"smoothing_signature": compute_smoothing_signature(context.config),
}
def run_smoothing(
method: str, context: TerrainContext, output_dir: Path, stem: str, original_model_path: Path
) -> dict[str, Any]:
"""방식에 따라 적절한 스무딩 엔진을 수행한다."""
if method == "dtm":
return smooth_dtm(context, output_dir, stem, original_model_path)
if method == "tin":
return smooth_tin(context, output_dir, stem, original_model_path)
raise ValueError(f"스무딩이 불가능한 지형 표현 방식입니다: {method}")
@@ -4,7 +4,7 @@
만들고, 뒤 단계(지면필터·모델·등고선·배수)는 받는 형식이 그대로라 손대지 않는다.
점이 임계를 넘으면 **칸(기본 0.5m)마다 최저점 하나만** 남긴다(씨닝). 설계가 쓰는 격자가
1m(지면필터 2m·CSF 천 1.5m)라 0.5m 는 설계보다 촘촘해 결과 표고가 사실상 같고, 30GB 두
1m(다듬기 칸 1m)라 0.5m 는 설계보다 촘촘해 결과 표고가 사실상 같고, 30GB 두
장이 메모리 29GB → 1.4GB 로 내려간다. 임계 아래면 원본 점을 그대로 쓴다 — 작은 자료의
결과는 바뀌지 않는다.
"""
@@ -157,8 +157,8 @@ class _ThinGrid:
버리지 않으므로 **지면 결과는 씨닝 전과 완전히 같다**(2026-09-06 용화 실측: 1m 지면
격자 141,969칸 전부 표고 차이 0). 지면점은 원본의 2~3%뿐이라 남겨도 가볍다.
나머지(수목·구조물·잡음)는 칸마다 최저점만 남긴다 — 분류가 없는 자료에서 CSF·PMF가
지면을 찾을 밑그림으로 충분하다(CSF 천 간격 1.5m > 칸 0.5m).
나머지(수목·구조물·잡음)는 칸마다 최저점만 남긴다 — 분류가 없는 자료에서 점진 TIN이
지면을 찾을 밑그림으로 충분하다(씨앗 간격 10m · 다듬기 칸 1m > 칸 0.5m).
한 청크 안에서 같은 칸이 여러 번 나오면 뒤에 쓴 값이 이겨 최저점이 아니게 된다.
그래서 청크를 (칸, 표고)로 정렬해 **칸마다 첫 점**만 골라 낸 뒤 격자와 견준다.
+1 -80
View File
@@ -1,7 +1,6 @@
"""B04 지표면 분석 FastAPI 라우터."""
import asyncio
import json
import logging
from pathlib import Path
from typing import Any
@@ -22,7 +21,7 @@ from B04_PreProcess.B04_PreProcess_Engine_Extent import (
project_epsg_from_prj,
)
from B04_PreProcess.B04_PreProcess_Engine_Ground import available_filters
from B04_PreProcess.B04_PreProcess_Engine_PointPreview import ensure_point_preview, sample_points
from B04_PreProcess.B04_PreProcess_Engine_PointPreview import ensure_point_preview
from B04_PreProcess.B04_PreProcess_Repository import (
get_input_file,
list_project_point_cloud_inputs,
@@ -54,12 +53,10 @@ from B04_PreProcess.B04_PreProcess_Schema import (
SurfaceConfirmedResponse,
SurfaceConfirmRequest,
SurfaceConfirmResponse,
SurfaceGroundStatsResponse,
SurfaceInputFileListResponse,
SurfaceInputFileSummary,
SurfaceModelListResponse,
SurfaceModelSummary,
SurfacePointCloudSampleResponse,
)
from B04_PreProcess.B04_PreProcess_Service import confirm_surface_selection
from common_util.common_util_auth import require_system_admin
@@ -356,51 +353,6 @@ async def get_surface_point_cloud_binary(
)
@router.get("/{project_id}/surface/point-cloud", response_model=SurfacePointCloudSampleResponse)
async def get_surface_point_cloud(
project_id: UUID,
filter: str | None = None,
) -> SurfacePointCloudSampleResponse | JSONResponse:
"""원본 또는 필터링된 B04 3D 미리보기 포인트를 반환한다(JSON · 화면은 이진 길을 씀)."""
try:
source_path = await _point_cloud_source(project_id, filter)
if isinstance(source_path, JSONResponse):
return source_path
with np.load(source_path) as structured:
xyz = np.asarray(structured["xyz"], dtype=np.float32)
bounds = np.asarray(structured["bounds"], dtype=np.float64)
point_count = (
int(structured["point_count"]) if "point_count" in structured else int(len(xyz))
)
rgb = np.asarray(structured["rgb"]) if "rgb" in structured else None
sample, rgb_sample = sample_points(xyz, rgb)
return SurfacePointCloudSampleResponse(
project_id=str(project_id),
point_count=point_count,
sampled_count=int(len(sample)),
bounds={
"x_min": float(bounds[0, 0]),
"x_max": float(bounds[0, 1]),
"y_min": float(bounds[1, 0]),
"y_max": float(bounds[1, 1]),
"z_min": float(bounds[2, 0]),
"z_max": float(bounds[2, 1]),
},
points=sample.astype(float).tolist(),
rgb=rgb_sample.astype(int).tolist() if rgb_sample is not None else None,
)
except LookupError as exc:
return JSONResponse(status_code=404, content={"status": "error", "message": str(exc)})
except Exception:
logger.exception("B04 포인트클라우드 샘플 조회 실패: project_id=%s", project_id)
return JSONResponse(
status_code=500,
content={"status": "error", "message": "포인트클라우드 조회 중 오류가 발생했습니다."},
)
@router.get("/{project_id}/surface/confirmed", response_model=SurfaceConfirmedResponse)
async def get_confirmed_surface(project_id: UUID) -> SurfaceConfirmedResponse | JSONResponse:
"""확정 지표면 구성과 지형 가장자리만 반환한다(포인트 배열 없음).
@@ -496,34 +448,3 @@ async def get_confirmed_surface(project_id: UUID) -> SurfaceConfirmedResponse |
status_code=500,
content={"status": "error", "message": "확정 지표면 정보를 불러오지 못했습니다."},
)
@router.get("/{project_id}/surface/ground-stats", response_model=SurfaceGroundStatsResponse)
async def get_surface_ground_stats(project_id: UUID) -> SurfaceGroundStatsResponse | JSONResponse:
"""manifest에서 필터별 지면 포인트 통계를 반환한다."""
pool = get_db_pool()
try:
async with pool.acquire() as connection:
stored_path = await get_project_storage_relative_path(connection, project_id)
project_root = Path(resolve_stored_project_path(stored_path))
manifest_path = project_root / "B04_PreProcess" / "models" / "manifest.json"
if not manifest_path.is_file():
return SurfaceGroundStatsResponse(project_id=str(project_id), filters={})
manifest = json.loads(manifest_path.read_text(encoding="utf-8"))
filters = {
key: {
"source_point_count": value.get("source_point_count"),
"methods": value.get("methods", {}),
}
for key, value in manifest.get("source_filters", {}).items()
if isinstance(value, dict)
}
return SurfaceGroundStatsResponse(project_id=str(project_id), filters=filters)
except LookupError as exc:
return JSONResponse(status_code=404, content={"status": "error", "message": str(exc)})
except Exception:
logger.exception("B04 지면 통계 조회 실패: project_id=%s", project_id)
return JSONResponse(
status_code=500,
content={"status": "error", "message": "지면 통계 조회 중 오류가 발생했습니다."},
)
@@ -29,10 +29,7 @@ logger = logging.getLogger(__name__)
router = APIRouter(prefix="/api/projects", tags=["B04 Surface Contour"])
MODEL_REPRESENTATIONS = {
"meshfree": "meshfree_surfels",
"dtm": "regular_grid",
"tin": "triangular_mesh",
"nurbs": "bspline_surface",
"implicit": "local_rbf_height_field",
}
@@ -79,7 +76,7 @@ async def get_surface_model_contour(
status_code=404,
content={"status": "error", "message": "해당 모델을 찾을 수 없습니다."},
)
model_type, model_file_path = row[0], row[1]
model_file_path = row[1]
project_root = Path(resolve_stored_project_path(stored_path))
if not model_file_path:
@@ -92,21 +89,16 @@ async def get_surface_model_contour(
models_dir = model_path.parent
structured_path = project_root / "B04_PreProcess" / "processed" / "structured.npz"
stem = model_path.stem
parts = stem.split("_")
if len(parts) >= 2:
method = parts[0]
filter_key = "_".join(parts[1:])
else:
method = model_type
filter_key = "csf"
if smooth and method in ("dtm", "tin"):
# 파일 이름 = `{표현}_{필터}` · 지운 표현(608 · tin 등)이면 404
method, _, filter_key = stem.partition("_")
if smooth and method == "dtm":
contour_filename = f"contour_{filter_key}_{method}_smooth_{interval}m.json"
contour_model_path = models_dir / f"{stem}_smooth.npz"
representation = "regular_grid" if method == "dtm" else "triangular_mesh"
representation = "regular_grid"
else:
contour_filename = f"contour_{filter_key}_{method}_{interval}m.json"
contour_model_path = model_path
representation = MODEL_REPRESENTATIONS.get(method)
representation = MODEL_REPRESENTATIONS.get(method) if filter_key else None
contour_path = models_dir / contour_filename
if recalculate or not _is_contour_cache_current(contour_path, contour_model_path):
if not math.isfinite(interval) or interval < 0.5:
@@ -188,11 +188,10 @@ async def get_surface_model_preview(
models_dir = model_path.parent
stem = model_path.stem
ext = "ply" if model_type == "meshfree" else "glb"
if smooth and model_type in ("dtm", "tin"):
if smooth and model_type == "dtm":
preview_filename = f"{stem}_smooth_preview.glb"
else:
preview_filename = f"{stem}_preview.{ext}"
preview_filename = f"{stem}_preview.glb"
preview_path = models_dir / preview_filename
if not preview_path.is_file():
@@ -204,11 +203,7 @@ async def get_surface_model_preview(
},
)
media_type = "application/octet-stream"
if ext == "glb":
media_type = "model/gltf-binary"
elif ext == "ply":
media_type = "application/ply"
media_type = "model/gltf-binary"
# 브라우저가 이미 같은 파일을 갖고 있으면 304만 돌려준다(새로고침이 빨라진다).
return cached_file_response(request, preview_path, media_type, preview_filename)
+3 -28
View File
@@ -4,19 +4,14 @@ from typing import Any
from pydantic import BaseModel, ConfigDict, Field
from B04_PreProcess.B04_PreProcess_Engine_Ground import available_filters
from config.config_system import (
SURFACE_MODEL_PRECOMPUTE,
SURFACE_MODEL_SOURCE_FILTERS,
)
# classification은 LAS가 지면분류를 싣고 있을 때 엔진이 자동으로 붙이지만,
# 재계산 요청에서 이름으로 지정하는 것도 허용한다.
_ALLOWED_FILTERS = set(SURFACE_MODEL_SOURCE_FILTERS) | {
"ransac",
"classification",
"classification_refine",
"ptd_refine",
}
# 허용 필터 = 지면 필터 등록부 하나(Ground `_FILTERS`)
_ALLOWED_FILTERS = set(available_filters())
_ALLOWED_METHODS = set(SURFACE_MODEL_PRECOMPUTE)
@@ -105,18 +100,6 @@ class SurfaceInputFileListResponse(BaseModel):
files: list[SurfaceInputFileSummary]
class SurfacePointCloudSampleResponse(BaseModel):
"""B04 3D 미리보기용 포인트클라우드 샘플."""
status: str = "success"
project_id: str
point_count: int
sampled_count: int
bounds: dict[str, float]
points: list[list[float]]
rgb: list[list[int]] | None = None
class SurfaceConfirmedResponse(BaseModel):
"""확정 지표면 요약 — 화면 진입 판정·준비화면·B05가 공통으로 쓰는 단일 출처.
@@ -138,14 +121,6 @@ class SurfaceConfirmedResponse(BaseModel):
route_bounds: dict[str, float] | None = None
class SurfaceGroundStatsResponse(BaseModel):
"""필터별 지면 포인트 통계 응답."""
status: str = "success"
project_id: str
filters: dict[str, dict[str, Any]]
class SurfaceAnalyzeResponse(BaseModel):
"""지표면 분석 실행 결과."""
@@ -18,9 +18,6 @@ export const SOURCE_FILTERS = [
"classification_refine",
"ptd",
"ptd_refine",
"grid_min_z",
"csf",
"pmf",
] as const;
const FILTER_LABELS: Record<string, string> = {
+2 -2
View File
@@ -50,12 +50,12 @@ import "./B04_PreProcess_UI_Style.css";
import "./B04_PreProcess_UI_Style_Map.css";
import { currentProjectId, guardProjectLeave } from "../A00_Common/b_project_switch";
const MODEL_METHODS = ["tin", "dtm", "nurbs", "implicit", "meshfree"] as const;
const MODEL_METHODS = ["dtm", "implicit"] as const;
const DEFAULT_FILTER = "ptd";
const DEFAULT_METHOD = "implicit"; // 확정 기본 표현(요구 558) — 확정본이 있으면 그것이 먼저
const ROUTE_STAGE = ROUTES.B05_PROFILE;
// 도엽 서피스 보간 방식 버튼 순서 — 백엔드 SHEET_SURFACE_METHODS와 같은 차례로 둔다.
const SHEET_METHOD_ORDER = ["tin_sheet", "tin", "biharmonic", "anudem", "multires", "laplace"];
const SHEET_METHOD_ORDER = ["laplace"];
function L(key: keyof typeof ui_locales): string {
return ui_locales[key][currentLanguageIndex];
@@ -1,7 +1,6 @@
import * as THREE from "three";
import { OrbitControls } from "three/examples/jsm/controls/OrbitControls.js";
import { GLTFLoader } from "three/examples/jsm/loaders/GLTFLoader.js";
import { PLYLoader } from "three/examples/jsm/loaders/PLYLoader.js";
import { API_BASE_URL } from "@config/config_frontend";
import { fetchCachedBytes } from "../A00_Common/b_asset_cache";
import { currentLanguageIndex, ui_locales } from "@ui/ui_template_locale";
@@ -75,7 +74,7 @@ export function createSurfaceTerrainViewer(): SurfaceTerrainViewer {
minValSpan,
progress,
} = buildTerrainViewerChrome();
let activeFilter = "csf";
let activeFilter = "";
let activeMethod = "dtm";
/** 진행률(0~1, 모르면 null)과 문구를 표시한다. label이 null이면 서클을 감춘다. */
@@ -348,7 +347,7 @@ export function createSurfaceTerrainViewer(): SurfaceTerrainViewer {
}
function syncSmoothingSupport(): void {
const supported = activeMethod === "tin" || activeMethod === "dtm";
const supported = activeMethod === "dtm";
smoothSelect.disabled = !supported;
smoothSelect.value = supported && smoothPreferred ? "on" : "off";
smoothLabel.classList.toggle("is-disabled", !supported);
@@ -372,8 +371,8 @@ export function createSurfaceTerrainViewer(): SurfaceTerrainViewer {
showProgress(null, "모델 조회 중…");
// 1. Find matching model in list
// model_type is TIN / DTM / NURBS / Implicit / Meshfree (we match activeMethod)
// model_file_path contains the activeFilter (e.g. csf, pmf, grid_min_z)
// model_type is DTM / Implicit (we match activeMethod)
// model_file_path contains the activeFilter (e.g. ptd_refine, classification)
const match = currentModelsList.find((m) => {
const typeMatches = m.model_type.toLowerCase() === activeMethod.toLowerCase();
const configuredFilter = m.generation_params?.source_filter;
@@ -391,7 +390,7 @@ export function createSurfaceTerrainViewer(): SurfaceTerrainViewer {
}
const modelId = match.id;
const isSmooth = (activeMethod === "tin" || activeMethod === "dtm") && smoothingOn();
const isSmooth = activeMethod === "dtm" && smoothingOn();
currentModelId = modelId;
currentModelSmooth = isSmooth;
const generation = ++loadGeneration;
@@ -410,52 +409,34 @@ export function createSurfaceTerrainViewer(): SurfaceTerrainViewer {
});
if (generation !== loadGeneration) return;
if (activeMethod === "meshfree") {
const geometry = new PLYLoader().parse(buffer);
geometry.computeBoundingSphere();
const material = new THREE.PointsMaterial({
size: 0.35,
vertexColors: geometry.hasAttribute("color"),
sizeAttenuation: true,
});
const points = new THREE.Points(geometry, material);
points.visible = surfCheck.checked;
terrainMesh = points;
scene.add(points);
fitCamera(points);
showProgress(1, "등고선을 그리는 중…");
await loadSelectedContours(modelId, isSmooth);
showProgress(null, null);
} else {
new GLTFLoader().parse(
buffer,
"",
async (gltf) => {
if (generation !== loadGeneration) {
disposeObject(gltf.scene);
return;
new GLTFLoader().parse(
buffer,
"",
async (gltf) => {
if (generation !== loadGeneration) {
disposeObject(gltf.scene);
return;
}
gltf.scene.traverse((child) => {
if (child instanceof THREE.Mesh) {
child.material.side = THREE.DoubleSide;
child.material.vertexColors = child.geometry.hasAttribute("color");
}
gltf.scene.traverse((child) => {
if (child instanceof THREE.Mesh) {
child.material.side = THREE.DoubleSide;
child.material.vertexColors = child.geometry.hasAttribute("color");
}
});
gltf.scene.visible = surfCheck.checked;
terrainMesh = gltf.scene;
scene.add(gltf.scene);
fitCamera(gltf.scene);
showProgress(1, "등고선을 그리는 중…");
await loadSelectedContours(modelId, isSmooth);
showProgress(null, null);
},
() => {
if (generation !== loadGeneration) return;
statusSpan.textContent = "3D 메쉬 파일이 없거나 로드할 수 없습니다.";
showProgress(null, null);
},
);
}
});
gltf.scene.visible = surfCheck.checked;
terrainMesh = gltf.scene;
scene.add(gltf.scene);
fitCamera(gltf.scene);
showProgress(1, "등고선을 그리는 중…");
await loadSelectedContours(modelId, isSmooth);
showProgress(null, null);
},
() => {
if (generation !== loadGeneration) return;
statusSpan.textContent = "3D 메쉬 파일이 없거나 로드할 수 없습니다.";
showProgress(null, null);
},
);
} catch (e) {
statusSpan.textContent = "3D 파일 로드에 실패했습니다.";
showProgress(null, null);
+5 -31
View File
@@ -28,7 +28,8 @@
가파른 지형(옆경사 중앙값 67 %)에서는 제약 안의 길이 폭 10~20 m 한 띠뿐이라 갈래끼리 가까울 수
있다 — 응답 `similar`(10 m 안인 다른 갈래) · `band_width_m`(띠 폭 어림)로 알린다(PLAN 23-4).
공통 제약(follow 밖 모든 갈래 · **넘지 못함**) — 최소 곡선반지름 · 계류 이격 · 곡선 MC = 일반측점:
공통 제약(follow 밖 모든 갈래 · **넘지 못함** · 측점 맞춤 · 합치기 없음 = 가이드 선 · 요구 611) —
최소 곡선반지름 · 계류 이격:
· 종단기울기 — 지반 기울기가 상한을 넘는 링크는 탐색에서 뺌(비탈 곧장 오르는 길 막이). 곡선
맞춤 뒤 측점(20 m) 사이 **원지반 기울기는 참고만**(계획선 높이로 푸는 값 · 요구 253-1 · 419) —
못 지킨 곳 수(`over_limit.total`)에 안 듦 · 목록엔 남음. 시종점 높이 차가 상한을 넘으면
@@ -64,11 +65,8 @@ import numpy as np
from B05_Profile.B05_Profile_Engine_RouteInitial_CrossEarth import cross_earthwork
from B05_Profile.B05_Profile_Engine_RouteInitial_Follow import (
FollowLand,
cross_delta,
cross_total,
merged_follow,
)
from B05_Profile.B05_Profile_Engine_RouteInitial_Merge import RouteEarth
from B05_Profile.B05_Profile_Engine_RouteInitial_Metrics import GRADE_SMOOTH_M, route_metrics
from B05_Profile.B05_Profile_Engine_RouteInitial_Pool import (
Pool,
@@ -94,7 +92,6 @@ from B05_Profile.B05_Profile_Engine_RouteInitial_Search import (
)
from B05_Profile.B05_Profile_Engine_RouteInitial_Station import (
off_station_count,
snap_curves_to_stations,
with_station_misses,
)
@@ -352,12 +349,6 @@ def generate_initial_route(
)
)
def fit(outline: Any, guard: Any = None) -> Any:
# R = 일반측점(383-1) — 갈래 후보는 직접 잼(`guard` · 406 · 거리 한도 없음) · 추적은 제약을
# 안 따짐(새 초기 노선 `snapped_initial` 과 같은 잣대 · 거리 한도 그대로).
moved = snap_curves_to_stations(outline, min_radius_m, guard=guard)
return outline if moved is None else moved
def follow_line(grid: Any = None) -> tuple[Any, dict[str, Any]]:
# 예상노선 추적 = 전처리 초기값(`planned_route_initial.csv`)과 같은 함수(곡선 최적 · 측점
# 맞춤 · 758-1 `..._Follow.merged_follow`) · 같은 입력 → 같은 선. 넘은 곳만 목록에 뜸.
@@ -392,13 +383,6 @@ def generate_initial_route(
design_curves=[curve.as_dict() for curve in line.curves],
)
def route_earth(line: Any) -> float:
# 합치기 멈춤 잣대 = 계획고 지반고 횡단(1차 종단 규칙 · 종단 끔과 같음)
return cross_total(line, terrain, standard)
def merge_delta(before: Any, after: Any, window: tuple[float, float]) -> float:
return cross_delta(before, after, window, terrain, standard)
def cross_score(route: Any, capped: bool) -> float:
# 후보 줄 세움(419) — 켬 = 상한 누른 계획선(`_reference`) 횡단 · 끔 = 계획고 지반고 횡단
held = pool.scores if pool is not None else {}
@@ -446,16 +430,7 @@ def generate_initial_route(
graph = RouteGraph(grid, expected, limit=max_grade, min_radius_m=min_radius_m)
while len(_POOLS) >= POOL_CACHE_SIZE:
_POOLS.pop(next(iter(_POOLS)))
# 합치기 멈춤 토공 = 노선 횡단 + 바뀐 측점 횡단 차이(420 `..._CrossEarth` · 419)
_POOLS[key] = build_pool(
graph,
follow_line(grid),
max_grade,
min_radius_m,
measure,
fit,
cross=RouteEarth(route_earth, merge_delta),
)
_POOLS[key] = build_pool(graph, follow_line(grid), max_grade, min_radius_m, measure)
pool = _POOLS[key]
bound = pool.bound
if target == "weighted":
@@ -464,8 +439,7 @@ def generate_initial_route(
# 토공 · 공사비 = 횡단 토공으로 줄 세움 · 종단 반영 켬이면 위 후보만 균형 종단(빠른
# 판) · 고른 하나만 아래에서 바로 셈(419)
outline, metrics = choose(target, pool.candidates, with_profile, quick, cross_score)
# 후보는 이미 측점에 맞춘 선(383-1 · `..._Pool.settle` · 맞추기 전 선은 안 냄 · 419)
snapped = (outline, metrics)
# 추천 = 계산 로직 결과 그대로(합치기 · 측점 맞춤 없음 · 요구 611) · 추적만 평면계획선 길
def finish(route: tuple[Any, dict[str, Any]]) -> dict[str, Any]:
# 정밀 종단은 낼 사본에만 — 묶음 캐시 후보는 빠른 판 그대로 둬 다음 고르기가 같은 잣대끼리
@@ -486,7 +460,7 @@ def generate_initial_route(
"metrics": values,
}
shown = finish(snapped)
shown = finish((outline, metrics))
ends = ends_grade(expected, terrain, shown["metrics"]["length_m"])
return {
"mode": mode,
@@ -12,9 +12,9 @@
· 못 지킨 곳 수(R · 계류 · 측점 — 원지반 기울기는 참고만 · 419 `hard_limits`)가 기준(엄격: 후보
최소 · 풂: 센 벌점 기준 길)보다 많은 후보는 뺀다. 탐색 링크 상한(비탈 곧장 오르는 길 막이)은
남김.
· R = 일반측점이 기본 제약(383-1) — 후보마다 곡선 MC 를 일반측점에 맞춘 선(`snap`)으로 재고
고르고 거름. 못 맞춘 곡선은 「못 지킨 구간」(station) · 고르기는 측점 밖 0 후보를 앞세움.
꺾임만(R 없는 교점)은 안 남김(`..._Radius.curve_every_kink` · 측점 맞춤 · 419).
· 추천 = 평면계획선을 긋는 가이드 선(요구 611 · 769-1) — 계산 로직 결과 그대로(경로 탐색 → 곡선
넣기 → 큰 R) · 같은 쪽 곡선 합치기 · 측점 맞춤은 안 함(평면계획선 길 `..._Follow.merged_follow`
몫). 꺾임만(R 없는 교점)은 안 남김(`..._Radius.curve_every_kink` · 419).
갈래는 **비교값으로** 고른다(`pick` · `GOALS`) — 토공 최소 = 절성토 합 · 기울기 순한 = 최대(같으면
평균) · 영선 노선 = 목표 기울기 편차 · 안전 우선 = 위험 통과 길이(옆경사 35° 이상 · 계류 버퍼 안) ·
@@ -38,7 +38,7 @@ from typing import Any
import numpy as np
from shapely.geometry import LineString
from B05_Profile.B05_Profile_Engine_RouteInitial_Merge import Guard, RouteEarth, merge_same_side
from B05_Profile.B05_Profile_Engine_RouteInitial_Merge import Guard
from B05_Profile.B05_Profile_Engine_RouteInitial_Metrics import (
COST_CROSSING,
COST_EARTH_M3,
@@ -52,7 +52,6 @@ from B05_Profile.B05_Profile_Engine_RouteInitial_Search import (
NoRoute,
RouteGraph,
)
from B05_Profile.B05_Profile_Engine_RouteInitial_Station import with_station_misses
from common_util.common_util_route_geometry import densify_route
from common_util.common_util_route_polyline import build_planned_polyline
@@ -87,8 +86,6 @@ TRACE_DENSIFY_M = 2.0
Measure = Callable[..., dict[str, Any]]
Route = tuple[Any, dict[str, Any]]
# (선, 직접 잼) → 곡선 MC 를 일반측점에 맞춘 선(못 맞추면 그 선 그대로)
Snap = Callable[[Any, Any], Any]
@dataclass
@@ -96,7 +93,6 @@ class Pool:
"""같은 입력의 후보 묶음 — 갈래 고르기 · 비중 탐색이 같이 쓴다."""
graph: RouteGraph
snap: Snap
candidates: list[Route] = field(default_factory=list)
bound: float | None = None
limit: int = 0 # 후보로 남길 못 지킨 곳 수 상한
@@ -111,10 +107,8 @@ class Pool:
follow: Route | None = None
# 제약 안에서 갈 수 있는 띠 폭 어림 — 후보끼리 가장 먼 거리(m)
band_width_m: float = 0.0
# 큰 R · 합치기 · 측점 맞춤이 거리 한도 대신 쓰는 직접 잼(406 · `..._Merge.Guard`)
# 큰 R 이 거리 한도 대신 쓰는 직접 잼(406 · `..._Merge.Guard`)
guard: Guard | None = None
# 합치기 멈춤 토공 = 횡단(노선 + 바뀐 측점 차이 · 420 · 없으면 띠 가늠 · `merge_same_side`)
cross: RouteEarth | None = None
# 후보 줄 세움 횡단 토공 기억(잰 값 id · 종단 반영) — 갈래 사이 다시 안 셈(419)
scores: dict[tuple[int, bool], float] = field(default_factory=dict)
@@ -127,15 +121,6 @@ def hard_limits(metrics: dict[str, Any]) -> dict[str, Any]:
return {**metrics, "over_limit": {**over, "total": over["total"] - over.get("grade", 0)}}
def settle(route: Route, snap: Snap, measure: Measure, **kw: Any) -> Route:
"""측점에 맞춘 선 + 잰 값(측점 밖 곡선 = 못 지킨 구간 · 383-1)."""
outline, metrics = route
moved = snap(outline)
if moved is not outline:
outline, metrics = moved, measure(moved, **kw)
return outline, with_station_misses(metrics, outline)
def _trace(
pool: Pool, nodes: list[tuple[float, float]], min_radius_m: float, measure: Measure, **kw: Any
) -> Route:
@@ -144,9 +129,9 @@ def _trace(
# 경로를 예정노선 규칙으로 단순화 · 묶음 · R 맞춤) 차례로, 위반(R · 계류 · 원지반 기울기 뺌)이
# 최소 R 판(노드마다 곡선 · 옛 규칙)보다 많지 않은 첫 판. 모두 화면이 [확인] 때 같은 선으로 다시
# 짜는 꼴. 같은 노드 · 같은 잼이면 다시 재지 않음(쓸기는 이웃 비중이 같은 길을 자주 냄).
# 고른 판의 꺾임만을 곡선 · 교점 지움으로(419) → 같은 쪽 이웃 곡선 합치기(406 · `..._Merge`) →
# 측점에 맞춤(`settle` · 383-1). 큰 R · 합치기 · 측점 맞춤 모두 거리 한도 대신 직접
# 잼(`pool.guard`) · 비중 계산은 예상노선 이격도(그 갈래 잣대 — 예상 몫을 줄이면 이격이 늘기만).
# 고른 판의 꺾임만을 곡선 · 교점 지움으로(419) — 끝. 같은 쪽 곡선 합치기 · 측점 맞춤은 안 함
# (추천 = 가이드 선 · 요구 611). 큰 R 은 거리 한도 대신 직접 잼(`pool.guard`) · 비중 계산은
# 예상노선 이격도(그 갈래 잣대 — 예상 몫을 줄이면 이격이 늘기만).
traced = pool.traced
guard = pool.guard
if guard is not None and kw.get("weighted"):
@@ -176,10 +161,7 @@ def _trace(
line = curve_every_kink(chosen[0], min_radius_m)
if line is not chosen[0]:
chosen = (line, measure(line, **kw))
chosen = merge_same_side(
chosen, guard, min_radius_m, lambda line: measure(line, **kw), pool.cross
)
traced[key] = settle(chosen, lambda line: pool.snap(line, guard), measure, **kw)
traced[key] = chosen
return traced[key]
@@ -252,13 +234,11 @@ def build_pool(
max_grade: float,
min_radius_m: float,
measure: Measure,
snap: Snap,
cross: RouteEarth | None = None,
) -> Pool:
"""후보 묶음 — 제약 탐색 후보 → 끝 노선 셋으로 비중 크기 맞춤 → 토공 ↔ 기울기 쓸기를 더함.
`measure` 는 원지반 기울기를 뺀 못 지킨 곳 수를 내야 함(`hard_limits`)."""
pool = Pool(graph, snap, follow=follow, guard=Guard(graph.grid, max_grade), cross=cross)
pool = Pool(graph, follow=follow, guard=Guard(graph.grid, max_grade))
pool.candidates, pool.bound, pool.limit = _candidates(pool, max_grade, min_radius_m, measure)
candidates, bound, limit = pool.candidates, pool.bound, pool.limit
ends = [follow[0], pick("earthwork", candidates)[0], pick("gentle", candidates)[0]]
@@ -41,8 +41,9 @@ from common_util.common_util_time import utc_now_text
# 3 = 곡선 MC 일반측점(351-2) · 4 = 다투는 곡선 끔(364-5) · 5 = 후보마다 맞춤(383-1) ·
# 6 = 같은 쪽 곡선 합치기 · 거리 한도 대신 직접 잼(406) · 7 = 측점 밖이 늘면 옛 차례 맞춤(406-4) ·
# 8 = 횡단 토공 · 표준 횡단면 = 저장 패널값 · 엔진 소스 지문(420) ·
# 9 = 측점 맞춤 덜 나쁜 자리에 원지반 기울기 늘어남 · 옮긴 거리(요구 606-2)
ALGO_VERSION = 9
# 9 = 측점 맞춤 덜 나쁜 자리에 원지반 기울기 늘어남 · 옮긴 거리(요구 606-2) ·
# 10 = 추천 4 갈래 합치기 · 측점 맞춤 뺌(가이드 선 · 요구 611)
ALGO_VERSION = 10
# 엔진 소스 지문에 드는 파일(저장소 뿌리 기준) — 추천 결과를 낳는 코드
SOURCE_GLOBS = (
"B05_Profile/B05_Profile_Engine_RouteInitial*.py",
@@ -7,7 +7,7 @@
(화면 묶음 stations · ground 단계와 같은 함수). 파이썬은 머리(schema_version · source ·
coordinate_reference · cad_exchange · options)만 붙임.
DTM · implicit(굳힌 격자 641-1) 지표면만(격자 이진 = `dtm_grid_file`) — 그 밖(tin · nurbs · meshfree) ·
DTM · implicit(굳힌 격자 641-1) 지표면만(격자 이진 = `dtm_grid_file`) — 그 밖 ·
Node 못 돌림 = None → 부르는 쪽이 파이썬 `generate_sections`(짝 · 로그 남김).
측점 열쇠(`key`) = 노선 저장 이월 열쇠(`B05_Profile_Route_Save._carry_cross`).
"""
-35
View File
@@ -82,8 +82,6 @@ def _sample_surface_on_grid(
return dtm_z
models_dir = Path(models_dir)
xx, yy = np.meshgrid(x_coords, y_coords)
query = np.column_stack([xx.ravel(), yy.ravel()])
def _finalize(z_flat: np.ndarray) -> np.ndarray:
z = np.asarray(z_flat, dtype=np.float64).reshape(len(y_coords), len(x_coords))
@@ -93,31 +91,6 @@ def _sample_surface_on_grid(
return z
try:
suffix = "_smooth" if smooth else ""
if method == "tin":
from scipy.interpolate import LinearNDInterpolator
path = models_dir / f"tin_{filter_key}{suffix}.npz"
if not path.exists():
path = models_dir / f"tin_{filter_key}.npz"
d = np.load(path)
verts = np.asarray(d["vertices"], dtype=np.float64)
interp = LinearNDInterpolator(verts[:, :2], verts[:, 2])
return _finalize(interp(query))
if method == "nurbs":
from scipy.interpolate import RectBivariateSpline
d = np.load(models_dir / f"nurbs_{filter_key}.npz")
cx = np.asarray(d["control_x"], dtype=np.float64)
cy = np.asarray(d["control_y"], dtype=np.float64)
cz = np.asarray(d["control_z"], dtype=np.float64)
degree = int(d["degree"][0]) if "degree" in d else 3
spline = RectBivariateSpline(
cy, cx, cz, kx=min(degree, len(cy) - 1), ky=min(degree, len(cx) - 1)
)
return _finalize(spline(y_coords, x_coords).ravel())
if method == "implicit":
# 굳힌 1 m 격자(641-1) 마디 그대로 — B05 종횡단 표본기와 같은 값 · RBF 다시 안 셈(763-1)
grid = DtmGridSampler.open_cached(implicit_grid_npz(models_dir, filter_key))
@@ -127,14 +100,6 @@ def _sample_surface_on_grid(
if off > 1e-3:
raise ValueError(f"굳힌 격자 축이 DTM 과 다름: implicit_{filter_key}_grid.npz")
return _finalize(grid.z[np.ix_(iy, ix)])
if method == "meshfree":
from scipy.interpolate import griddata
d = np.load(models_dir / f"meshfree_{filter_key}.npz")
pts = np.asarray(d["points"], dtype=np.float64)
z = griddata(pts[:, :2], pts[:, 2], query, method="linear")
return _finalize(z)
except FileNotFoundError as error:
name = Path(error.filename).name if error.filename else f"{method}_{filter_key}"
raise FileNotFoundError(f"지표면 모델 파일이 없습니다: {name}") from error
+2 -4
View File
@@ -64,10 +64,8 @@ class RouteSolveRequest(BaseModel):
model_config = ConfigDict(extra="forbid")
filter_key: SurfaceFilter = Field(description="지면 필터 키 (grid_min_z/csf/pmf)")
method: SurfaceMethod = Field(
default="dtm", description="지표면 표현 (dtm/tin/nurbs/implicit/meshfree)"
)
filter_key: SurfaceFilter = Field(description="지면 필터 키 (classification/ptd · _refine)")
method: SurfaceMethod = Field(default="dtm", description="지표면 표현 (dtm/implicit)")
smooth: bool = Field(default=False)
surface_model_id: int | None = Field(default=None, description="기반 지표면 모델 id")
algorithm: str = Field(
+1 -1
View File
@@ -54,7 +54,7 @@ export class TerrainBandSplit {
terrain.updateMatrixWorld(true);
terrain.traverse((child) => {
if (child instanceof THREE.Points) {
// 포인트클라우드(meshfree)는 통째로 near에 둔다 — 클리핑이 점 단위라 가볍다.
// 포인트클라우드는 통째로 near에 둔다 — 클리핑이 점 단위라 가볍다.
this.near.push({
geometry: withOriginalColors(
child.geometry.clone().applyMatrix4(child.matrixWorld),
+1 -7
View File
@@ -3,7 +3,6 @@ import { OrbitControls } from "three/examples/jsm/controls/OrbitControls.js";
import { createTerrainCompass } from "@ui/ui_template_compass";
import { contourMaterials, isMajorContour } from "@ui/ui_template_contour3d";
import { GLTFLoader } from "three/examples/jsm/loaders/GLTFLoader.js";
import { PLYLoader } from "three/examples/jsm/loaders/PLYLoader.js";
import { API_BASE_URL } from "@config/config_frontend";
import { contourAssetUrl, fetchCachedBytes, fetchCachedJson } from "../A00_Common/b_asset_cache";
import { bindCursorPivotControls, TOP_VIEW_TILT } from "../B04_PreProcess/B04_PreProcess_UI_Camera";
@@ -453,12 +452,7 @@ export function createRouteViewer(): RouteViewer {
step("fetch", fetched);
const parsed = performance.now();
terrain = await new Promise<THREE.Object3D>((resolve, reject) => {
if (method === "meshfree") {
const geometry = new PLYLoader().parse(buffer);
resolve(new THREE.Points(geometry, new THREE.PointsMaterial({ size: 0.35 })));
} else {
new GLTFLoader().parse(buffer, "", (gltf) => resolve(gltf.scene), reject);
}
new GLTFLoader().parse(buffer, "", (gltf) => resolve(gltf.scene), reject);
});
terrain.traverse((child) => {
if (child instanceof THREE.Mesh) child.material.side = THREE.DoubleSide;
@@ -256,7 +256,7 @@ def _enforce_stored_designs(
def design_grid_file(project_root: Path, params: dict[str, Any] | None) -> Path | None:
"""설계 지형이 DTM · implicit 면 격자 이진(`dtm_grid_file` · 표본기 `build_surface_sampler` 와
같은 npz · implicit = 굳힌 격자 641-1) · 그 밖 지표면(tin · nurbs · …) · 못 엶 = None
같은 npz · implicit = 굳힌 격자 641-1) · 못 엶 = None
(파이썬 길)."""
from common_util.common_util_surface_sampler import dtm_grid_file, implicit_grid_npz
@@ -14,7 +14,6 @@ export interface DesignDrawingItem {
| "watershed"
| "plan"
| "landuse"
| "plan_lidar"
| "cross_standard"
| "structure"
| "development"
@@ -119,7 +118,6 @@ export interface DesignDrawingResponse {
| "watershed"
| "plan"
| "landuse"
| "plan_lidar"
| "cross_standard"
| "structure"
| "development"
@@ -1,54 +1,15 @@
"""B07 계획평면도(라이다) CAD 조립 — 지표면 격자를 음영기복 그림으로 깔고 노선을 얹는다.
"""B07 지표면 음영기복 그림 — 계획평면도 배경 · D03 지형 층이 같이 씀.
사용자 지시(2026-09-04) — 「라이다 계획평면도는 3D 자료를 탑뷰에서 본 그림이 필요함.
가능한 범위에서 일단 배치해 주면 보고 개선하겠음」.
점구름을 그대로 그리면 수천만 점이라 도면 만들기가 느려진다(용화_LAS 실측 4,900만 점).
이미 만들어 둔 **지표면 격자(DTM)** 로 음영기복 이미지를 서버에서 만들어 배경으로 깐다.
도면 틀이 이미지 요소를 받아 주므로(`Image` 엔티티) PNG 를 그대로 싣는다.
축척·도곽·장 나눔은 계획평면도와 같다(1/1,200 고정) — 같은 자리에 노선이 서야 한다.
이미 만들어 둔 **지표면 격자(DTM)** 로 음영기복 PNG 를 서버에서 만든다 — 점구름을 그대로
그리면 수천만 점이라 느려진다(용화_LAS 실측 4,900만 점). 옛 「계획평면도(라이다)」 짓기는 걷음(608).
"""
import base64
import io
import math
from typing import Any
import numpy as np
from B07_DesignDetail.B07_DesignDetail_Engine_Cad import (
DRAWING_FORMAT,
FRAME_LAYER_ID,
TABLE_LABEL_COLOR,
_layer,
_text_entity,
polyline_entity,
)
from B07_DesignDetail.B07_DesignDetail_Engine_Cad_Plan import (
_COMPASS_MARGIN,
_COMPASS_SIZE,
_FONT_SIZE,
_ROUTE_WIDTH,
_TITLE_FONT_SIZE,
MM,
)
from B07_DesignDetail.B07_DesignDetail_Engine_Template import (
compass_entities,
entities_bbox,
frame_entities,
scale_fields,
)
from config.config_system import DRAWING_SCALE_PLAN
LIDAR_KIND = "plan_lidar"
LIDAR_LABEL = "계획평면도(라이다)"
SHADE_LAYER_ID = "b07-lidar-shade"
ROUTE_LAYER_ID = "b07-lidar-route"
TITLE_LAYER_ID = "b07-lidar-title"
ROUTE_COLOR = "#ffe066"
# 음영기복 광원 — 도면 관행대로 북서(방위각 315°)에서 45° 높이로 비춘다.
_AZIMUTH_DEG = 315.0
_ALTITUDE_DEG = 45.0
@@ -101,120 +62,3 @@ def hillshade_png(z: np.ndarray, valid: np.ndarray, resolution_m: float) -> tupl
image.save(buffer, format="PNG", optimize=True)
data_url = "data:image/png;base64," + base64.b64encode(buffer.getvalue()).decode("ascii")
return (data_url, image.width, image.height)
def build_lidar_plan_drawing(
drawing_id: str,
label: str,
route_xy: list[tuple[float, float]],
shade_image: str | None,
shade_box: tuple[float, float, float, float] | None,
) -> dict[str, Any]:
"""라이다 계획평면도 한 장을 만든다.
`shade_box`는 음영기복 그림이 덮는 실좌표 범위(min_x, min_y, max_x, max_y)다 —
그림 네 모서리를 그 범위 그대로 종이에 놓아야 노선과 좌표가 맞는다.
"""
everything = [*route_xy]
if shade_box:
everything.extend([(shade_box[0], shade_box[1]), (shade_box[2], shade_box[3])])
if not everything:
raise FileNotFoundError(
"라이다 계획평면도에 그릴 자료가 없습니다. B04 전처리에서 지표면을 먼저 만드세요."
)
min_x = min(x for x, _y in everything)
min_y = min(y for _x, y in everything)
def paper(point: tuple[float, float]) -> tuple[float, float]:
return ((point[0] - min_x) * MM, (point[1] - min_y) * MM)
entities: list[dict[str, Any]] = []
if shade_image and shade_box:
left, bottom = paper((shade_box[0], shade_box[1]))
right, top = paper((shade_box[2], shade_box[3]))
entities.append(
{
"id": f"{drawing_id}:shade",
"type": "Image",
"lineColor": "#ffffff",
"lineWidth": 1,
"layerId": SHADE_LAYER_ID,
"shapeData": {
"points": [
{"x": left, "y": bottom},
{"x": right, "y": bottom},
{"x": right, "y": top},
{"x": left, "y": top},
],
"imageData": shade_image,
},
}
)
route = polyline_entity(
drawing_id,
[paper(point) for point in route_xy],
ROUTE_LAYER_ID,
ROUTE_COLOR,
width=_ROUTE_WIDTH,
)
if route:
entities.append(route)
map_bbox = entities_bbox(entities)
if map_bbox:
entities.extend(
compass_entities(
drawing_id,
(
map_bbox[2] + _COMPASS_MARGIN + _COMPASS_SIZE / 2.0,
map_bbox[3] - _COMPASS_SIZE / 2.0,
),
_COMPASS_SIZE,
)
)
bbox = entities_bbox(entities)
if bbox:
min_bx, _min_by, max_bx, max_by = bbox
entities.append(
_text_entity(
f"{drawing_id}:title",
label,
(min_bx + max_bx) / 2.0,
max_by + 12.0,
TITLE_LAYER_ID,
_TITLE_FONT_SIZE,
TABLE_LABEL_COLOR,
)
)
entities.append(
_text_entity(
f"{drawing_id}:scale",
f"S = 1/{DRAWING_SCALE_PLAN:,}",
max_bx,
max_by + 5.0,
TITLE_LAYER_ID,
_FONT_SIZE,
TABLE_LABEL_COLOR,
align="right",
)
)
entities.extend(
frame_entities(
drawing_id,
entities_bbox(entities) or bbox,
fit=False,
fields={"도면명": label, **scale_fields(("", DRAWING_SCALE_PLAN))},
)
)
return {
"format": DRAWING_FORMAT,
"entities": entities,
"layers": [
_layer(SHADE_LAYER_ID, "지표면 음영기복", locked=True),
_layer(ROUTE_LAYER_ID, "계획노선"),
_layer(TITLE_LAYER_ID, "표제"),
_layer(FRAME_LAYER_ID, "도각", locked=True),
],
}
@@ -69,7 +69,6 @@ from B07_DesignDetail.B07_DesignDetail_Engine_Template import add_title_fields
from B07_DesignDetail.B07_DesignDetail_Router_Support_Basin import ( # noqa: F401
_BASIN_MAX_DISTANCE_M,
CONTOUR_FILE,
LIDAR_ID,
PLAN_ID,
STREAM_FILE,
_basins_crs,
@@ -79,7 +78,6 @@ from B07_DesignDetail.B07_DesignDetail_Router_Support_Basin import ( # noqa: F4
_geometry_lines,
_too_far_from_route,
clip_line_to_box,
lidar_source,
plan_source,
plan_stations,
watershed_source,
@@ -18,7 +18,7 @@ from shapely.geometry import box as shapely_box
from B04_PreProcess.B04_PreProcess_Router_Watershed import CONTOUR_FILE, STREAM_FILE
from B07_DesignDetail.B07_DesignDetail_Engine_Cad_Basin import map_area_mm
from B07_DesignDetail.B07_DesignDetail_Engine_Cad_Lidar import LIDAR_LABEL, hillshade_png
from B07_DesignDetail.B07_DesignDetail_Engine_Cad_Lidar import hillshade_png
from B07_DesignDetail.B07_DesignDetail_Engine_Cad_Plan import (
plan_area_mm,
plan_chunks,
@@ -400,13 +400,10 @@ def _clip_rings(
return rings
LIDAR_ID = re.compile(r"^plan_lidar(?:_(\d+))?$")
def _sheet_box(
route_xy: list[tuple[float, float]], scale: int = DRAWING_SCALE_PLAN
) -> tuple[float, float, float, float]:
"""그 장의 도곽 범위(실좌표 m) — 계획평면도·용지도·라이다가 같은 규칙을 쓴다."""
"""그 장의 도곽 범위(실좌표 m) — 계획평면도·용지도가 같은 규칙을 쓴다."""
area_w_mm, area_h_mm = plan_area_mm()
half_w = area_w_mm / 2.0 * scale / 1000.0
half_h = area_h_mm / 2.0 * scale / 1000.0
@@ -438,36 +435,6 @@ def _chunk_route(
return route_xy, chunk, total
def lidar_source(context: Any, longitudinal: dict[str, Any], drawing_id: str) -> dict[str, Any]:
"""라이다 계획평면도 한 장의 입력(노선 + 지표면 음영기복 그림)을 모은다.
지표면은 확정 DTM 격자(`dtm_{필터}[_smooth].npz`)를 도곽 범위로 잘라 쓴다 —
점구름을 그대로 그리면 수천만 점이라 도면 만들기가 느려진다(2026-09-04 사용자 지시).
"""
match = LIDAR_ID.fullmatch(drawing_id)
if not match:
raise ValueError("올바르지 않은 라이다 계획평면도 ID입니다.")
number = int(match.group(1)) if match.group(1) else 1
route_xy, chunk, total = _chunk_route(context, longitudinal, number)
box = _sheet_box(route_xy)
label = LIDAR_LABEL if total <= 1 else f"{LIDAR_LABEL} {chunk['number']}장"
shade_image: str | None = None
shade_box: tuple[float, float, float, float] | None = None
try:
shade_image, shade_box = _hillshade_for_box(context, box)
except (FileNotFoundError, ValueError, OSError) as exc:
# 지표면이 없어도 노선·도각은 그린다 — 빈 화면보다 낫다.
logger.warning("B07 라이다 계획평면도: 음영기복을 만들지 못했습니다 — %s", exc)
return {
"label": label,
"route_xy": route_xy,
"shade_image": shade_image,
"shade_box": shade_box,
}
def _hillshade_for_box(
context: Any, box: tuple[float, float, float, float]
) -> tuple[str, tuple[float, float, float, float]]:
@@ -475,7 +442,7 @@ def _hillshade_for_box(
import numpy as np
params = getattr(context, "surface_params", None) or {}
source_filter = str(params.get("source_filter") or "csf")
source_filter = str(params.get("source_filter") or "")
smooth = bool(params.get("smooth", True)) or params.get("method") == "implicit" # 641-2
models_dir = Path(context.project_root) / "B04_PreProcess" / "models"
candidates = [models_dir / f"dtm_{source_filter}_smooth.npz"] if smooth else []
@@ -18,7 +18,6 @@ class DesignDrawingItem(BaseModel):
"watershed",
"plan",
"landuse",
"plan_lidar",
"cross_standard",
# 구조물 도면(376) — B05 · B06 배치에 연결된 구조물 조합 도면.
"structure",
@@ -64,7 +63,6 @@ class DesignDrawingResponse(BaseModel):
"watershed",
"plan",
"landuse",
"plan_lidar",
"cross_standard",
# 구조물 도면(376) — B05 · B06 배치에 연결된 구조물 조합 도면.
"structure",
+2 -68
View File
@@ -24,25 +24,16 @@ from typing import Literal, Protocol, get_args
import numpy as np
from scipy.interpolate import RegularGridInterpolator
#: 지면 필터 · 지표면 표현 허용 목록 — 파일 이름에 들어감(B04 `_FILTERS` · `MODEL_METHODS`)
#: 지면 필터 · 지표면 표현 허용 목록 — 파일 이름에 들어감(B04 `_FILTERS` · `BUILDERS`)
#: `sheet_*` = 수치지형도 등고선 서피스(B04 `SHEET_METHOD_LABELS` 키 앞에 `sheet_`)
SurfaceFilter = Literal[
"classification",
"grid_min_z",
"csf",
"pmf",
"ptd",
"ransac",
"classification_refine",
"ptd_refine",
"sheet_tin_sheet",
"sheet_tin",
"sheet_biharmonic",
"sheet_anudem",
"sheet_multires",
"sheet_laplace",
]
SurfaceMethod = Literal["dtm", "tin", "nurbs", "implicit", "meshfree"]
SurfaceMethod = Literal["dtm", "implicit"]
def check_surface_keys(filter_key: str, method: str = "dtm") -> None:
@@ -243,31 +234,12 @@ class CallableSurfaceSampler:
return values, valid
@dataclass
class InterpolatedSurfaceSampler:
"""불규칙/곡면 모델 보간기와 DTM footprint 유효성을 결합한다."""
interpolator: Callable[[np.ndarray], np.ndarray]
footprint: SurfaceElevationSampler | None = None
def sample_xy(self, xy: np.ndarray) -> tuple[np.ndarray, np.ndarray]:
xy = np.asarray(xy, dtype=np.float64)
values = np.asarray(self.interpolator(xy), dtype=np.float64).reshape(-1)
valid = np.isfinite(values)
if self.footprint is not None:
_, footprint_valid = self.footprint.sample_xy(xy)
valid &= footprint_valid
values[~valid] = np.nan
return values, valid
def build_surface_sampler(
models_dir: Path | str, source_filter: str, method: str, smooth: bool
) -> SurfaceElevationSampler:
"""1단계 확정 모델을 종·횡단용 일괄 XY 표고 sampler로 연다."""
check_surface_keys(source_filter, method)
models_dir = Path(models_dir)
smooth_suffix = "_smooth" if smooth and method in {"dtm", "tin"} else ""
dtm_smooth = models_dir / f"dtm_{source_filter}_smooth.npz"
dtm_original = models_dir / f"dtm_{source_filter}.npz"
dtm_path = dtm_smooth if smooth and dtm_smooth.exists() else dtm_original
@@ -277,46 +249,8 @@ def build_surface_sampler(
if method == "dtm":
return footprint
if method == "tin":
from scipy.interpolate import LinearNDInterpolator
path = models_dir / f"tin_{source_filter}{smooth_suffix}.npz"
if not path.exists() and smooth_suffix:
path = models_dir / f"tin_{source_filter}.npz"
with np.load(path, allow_pickle=False) as data:
vertices = np.asarray(data["vertices"], dtype=np.float64)
interpolator = LinearNDInterpolator(vertices[:, :2], vertices[:, 2], fill_value=np.nan)
return InterpolatedSurfaceSampler(lambda xy: interpolator(xy), footprint)
if method == "nurbs":
from scipy.interpolate import RectBivariateSpline
path = models_dir / f"nurbs_{source_filter}.npz"
with np.load(path, allow_pickle=False) as data:
control_x = np.asarray(data["control_x"], dtype=np.float64)
control_y = np.asarray(data["control_y"], dtype=np.float64)
control_z = np.asarray(data["control_z"], dtype=np.float64)
degree = int(data["degree"][0]) if "degree" in data else 3
spline = RectBivariateSpline(
control_y,
control_x,
control_z,
kx=min(degree, len(control_y) - 1),
ky=min(degree, len(control_x) - 1),
)
return InterpolatedSurfaceSampler(lambda xy: spline.ev(xy[:, 1], xy[:, 0]), footprint)
if method == "implicit":
# 굳힌 격자(641-1) — 브라우저 · Node 측점 · B06 이 같은 파일을 읽음(`source` 가 붙음)
return DtmGridSampler.open_cached(implicit_grid_npz(models_dir, source_filter))
if method == "meshfree":
from scipy.interpolate import LinearNDInterpolator
path = models_dir / f"meshfree_{source_filter}.npz"
with np.load(path, allow_pickle=False) as data:
points = np.asarray(data["points"], dtype=np.float64)
interpolator = LinearNDInterpolator(points[:, :2], points[:, 2], fill_value=np.nan)
return InterpolatedSurfaceSampler(lambda xy: interpolator(xy), footprint)
raise ValueError(f"지원하지 않는 지표면 모델입니다: {method}")
+7 -64
View File
@@ -8,12 +8,8 @@ import os
# 5. 지형 분석 알고리즘 파라미터
# ─────────────────────────────────────────────────────────────────────────
# Trimesh 메쉬 생성
MESH_GRID_SIZE = float(os.getenv("MESH_GRID_SIZE", "1.0")) # 미터 단위
MESH_SMOOTHING_ITERATIONS = int(os.getenv("MESH_SMOOTHING_ITERATIONS", "0"))
SURFACE_LAS_CHUNK_SIZE = int(os.getenv("SURFACE_LAS_CHUNK_SIZE", "500000"))
SURFACE_DEFAULT_RGB_VALUE = int(os.getenv("SURFACE_DEFAULT_RGB_VALUE", "128"))
SURFACE_GRID_CELL_SIZE_M = float(os.getenv("SURFACE_GRID_CELL_SIZE_M", "2.0"))
# 지형 파일 여러 장 병합 (2026-09-06 사용자 확정)
# ─────────────────────────────────────────────────────────────────────────
@@ -24,45 +20,16 @@ SURFACE_MERGE_MAX_GAP_M = float(os.getenv("SURFACE_MERGE_MAX_GAP_M", "5000"))
# 점이 이 수를 넘으면 아래 칸 크기로 씨닝한다(칸마다 최저점 하나). 넘지 않으면 원본
# 그대로 쓴다 — 작은 자료의 결과는 바뀌지 않는다. 1억점 = 메모리 약 3.2GB.
SURFACE_THIN_TRIGGER_POINTS = int(os.getenv("SURFACE_THIN_TRIGGER_POINTS", "100000000"))
# 씨닝 칸 크기(m). 설계 격자 1m·지면필터 2m·CSF 천 1.5m 보다 촘촘해야 결과가 안 변한다.
# 씨닝 칸 크기(m). 설계 격자 1m·다듬기 칸 1m 보다 촘촘해야 결과가 안 변한다.
SURFACE_THIN_CELL_SIZE_M = float(os.getenv("SURFACE_THIN_CELL_SIZE_M", "0.5"))
# 씨닝 격자가 이보다 많아지면 범위가 비정상이다(먼 파일이 섞였거나 좌표계 불일치).
SURFACE_THIN_MAX_CELLS = int(os.getenv("SURFACE_THIN_MAX_CELLS", "400000000"))
SURFACE_GRID_HEIGHT_THRESHOLD_M = float(os.getenv("SURFACE_GRID_HEIGHT_THRESHOLD_M", "1.5"))
# CSF (Cloth Simulation Filter) 지면 분류 파라미터
SURFACE_CSF_CLOTH_RESOLUTION_M = float(os.getenv("SURFACE_CSF_CLOTH_RESOLUTION_M", "1.5"))
SURFACE_CSF_RIGIDNESS = int(os.getenv("SURFACE_CSF_RIGIDNESS", "1"))
SURFACE_CSF_TIME_STEP = float(os.getenv("SURFACE_CSF_TIME_STEP", "0.65"))
SURFACE_CSF_CLASS_THRESHOLD_M = float(os.getenv("SURFACE_CSF_CLASS_THRESHOLD_M", "0.5"))
# 낙하 후 스프링 정착에 쓰는 여유 반복 수 (실제 낙하 구간은 현장 기복에서 계산한다)
SURFACE_CSF_ITERATIONS = int(os.getenv("SURFACE_CSF_ITERATIONS", "150"))
# 낙하+정착 반복 상한 — 기복 1500m(0.3185m/회)까지 커버한다
SURFACE_CSF_MAX_ITERATIONS = int(os.getenv("SURFACE_CSF_MAX_ITERATIONS", "5000"))
SURFACE_CSF_SLOPE_SMOOTH = os.getenv("SURFACE_CSF_SLOPE_SMOOTH", "True").lower() == "true"
SURFACE_CSF_SLOPE_SMOOTH_THRESHOLD_M = float(
os.getenv("SURFACE_CSF_SLOPE_SMOOTH_THRESHOLD_M", "1.8")
)
# 점진 TIN(PTD) 지면 분류 — 장면 flat · gentle · steep(각 · 거리 묶음) · 시작 간격(씨앗 칸)
# 시작 간격 = 숲 10 m(DJI 기본 20 m 는 급한 볼록 비탈을 못 따라감)
SURFACE_PTD_SCENE = os.getenv("SURFACE_PTD_SCENE", "steep")
SURFACE_PTD_SEED_SPACING_M = float(os.getenv("SURFACE_PTD_SEED_SPACING_M", "10.0"))
# PMF (Progressive Morphological Filter) 지면 분류 파라미터
SURFACE_PMF_CELL_SIZE_M = float(os.getenv("SURFACE_PMF_CELL_SIZE_M", "2.0"))
SURFACE_PMF_MAX_WINDOW_SIZE = int(os.getenv("SURFACE_PMF_MAX_WINDOW_SIZE", "40"))
SURFACE_PMF_INITIAL_WINDOW_SIZE = int(os.getenv("SURFACE_PMF_INITIAL_WINDOW_SIZE", "3"))
SURFACE_PMF_SLOPE = float(os.getenv("SURFACE_PMF_SLOPE", "1.0"))
SURFACE_PMF_MAX_DISTANCE_M = float(os.getenv("SURFACE_PMF_MAX_DISTANCE_M", "2.5"))
# RANSAC (Local plane fitting) 지면 분류 파라미터
SURFACE_RANSAC_DISTANCE_THRESHOLD_M = float(os.getenv("SURFACE_RANSAC_DISTANCE_THRESHOLD_M", "0.3"))
SURFACE_RANSAC_N = int(os.getenv("SURFACE_RANSAC_N", "3"))
SURFACE_RANSAC_ITERATIONS = int(os.getenv("SURFACE_RANSAC_ITERATIONS", "100"))
SURFACE_RANSAC_LOCAL_GRID_SIZE_M = float(os.getenv("SURFACE_RANSAC_LOCAL_GRID_SIZE_M", "10.0"))
SURFACE_RANSAC_SEED = int(os.getenv("SURFACE_RANSAC_SEED", "42"))
# LAS 자체 지면분류(class 2) 필터 — 이 비율 이상 분류돼 있을 때만 필터로 제공한다.
# 미분류 LAS는 classification이 전부 0이라 마스크가 비어 버린다.
SURFACE_CLASSIFIED_GROUND_MIN_RATIO = float(
@@ -83,12 +50,12 @@ SURFACE_GROUND_RATIO_WARN = float(os.getenv("SURFACE_GROUND_RATIO_WARN", "0.01")
SURFACE_ROUTE_EDGE_TRIM_M = float(os.getenv("SURFACE_ROUTE_EDGE_TRIM_M", "3.0"))
# ─────────────────────────────────────────────────────────────────────────
# 5-2. 지표면 모델 생성 파라미터 (TIN/DTM/NURBS/implicit/meshfree)
# 5-2. 지표면 모델 생성 파라미터 (DTM/implicit)
# ─────────────────────────────────────────────────────────────────────────
# 고를 수 있는 지면 필터 전체 — 관리자가 B04 드롭다운에서 요청할 때 쓰인다.
# 자동 전처리가 미리 만드는 조합이 아니다(그건 SURFACE_AUTO_* 가 정한다).
SURFACE_MODEL_SOURCE_FILTERS = tuple(
os.getenv("SURFACE_MODEL_SOURCE_FILTERS", "classification,ptd,grid_min_z,csf,pmf").split(",")
os.getenv("SURFACE_MODEL_SOURCE_FILTERS", "classification,ptd").split(",")
)
# 자동 전처리가 만드는 조합 — 기본 필터 1종 × 아래 표현. 나머지는 관리자 요청 시 만든다.
SURFACE_AUTO_METHODS = tuple(os.getenv("SURFACE_AUTO_METHODS", "dtm").split(","))
@@ -103,12 +70,10 @@ SURFACE_REFINE_WINDOW_M = float(os.getenv("SURFACE_REFINE_WINDOW_M", "18.0"))
SURFACE_REFINE_THRESHOLD_M = float(os.getenv("SURFACE_REFINE_THRESHOLD_M", "0.5"))
SURFACE_REFINE_SCALAR = float(os.getenv("SURFACE_REFINE_SCALAR", "1.25"))
SURFACE_REFINE_LOW_SLOPE = float(os.getenv("SURFACE_REFINE_LOW_SLOPE", "5.0"))
# dtm을 먼저 빌드해야 TIN 등고선 사전 캐시가 dtm footprint를 참조할 수 있다 (PLAN D-6)
SURFACE_MODEL_PRECOMPUTE = tuple(
os.getenv("SURFACE_MODEL_PRECOMPUTE", "dtm,tin,nurbs,implicit,meshfree").split(",")
)
# dtm 먼저 — implicit 이 dtm 메운 칸을 길잡이로 씀(766-2)
SURFACE_MODEL_PRECOMPUTE = tuple(os.getenv("SURFACE_MODEL_PRECOMPUTE", "dtm,implicit").split(","))
SURFACE_MODEL_SMOOTHING_METHODS = tuple(
os.getenv("SURFACE_MODEL_SMOOTHING_METHODS", "dtm,tin").split(",")
os.getenv("SURFACE_MODEL_SMOOTHING_METHODS", "dtm").split(",")
)
SURFACE_MODEL_SYNC_TIMEOUT_SECONDS = int(os.getenv("SURFACE_MODEL_SYNC_TIMEOUT_SECONDS", "0"))
@@ -124,20 +89,12 @@ SURFACE_TILE_SIZE_M = float(os.getenv("SURFACE_TILE_SIZE_M", "50.0"))
SURFACE_MAX_PREVIEW_VERTICES = int(os.getenv("SURFACE_MAX_PREVIEW_VERTICES", "500000"))
# 표현별 파라미터 (old 버전 검증값 기준 — PLAN F)
SURFACE_TIN_MAX_INPUT_POINTS = int(os.getenv("SURFACE_TIN_MAX_INPUT_POINTS", "500000"))
SURFACE_DTM_GRID_RESOLUTION_M = float(os.getenv("SURFACE_DTM_GRID_RESOLUTION_M", "1.0"))
SURFACE_NURBS_DEGREE = int(os.getenv("SURFACE_NURBS_DEGREE", "3"))
SURFACE_NURBS_PATCH_SIZE_M = float(os.getenv("SURFACE_NURBS_PATCH_SIZE_M", "50.0"))
SURFACE_NURBS_CONTROL_POINTS_PER_AXIS = int(
os.getenv("SURFACE_NURBS_CONTROL_POINTS_PER_AXIS", "16")
)
# implicit 중심점 수 = clip(지면점 × 비율, 하한, 상한) · 지면점이 하한보다 적으면 전부(요구 558).
SURFACE_IMPLICIT_POINT_RATIO = float(os.getenv("SURFACE_IMPLICIT_POINT_RATIO", "0.175"))
SURFACE_IMPLICIT_MIN_POINTS = int(os.getenv("SURFACE_IMPLICIT_MIN_POINTS", "50000"))
SURFACE_IMPLICIT_MAX_POINTS = int(os.getenv("SURFACE_IMPLICIT_MAX_POINTS", "2000000"))
SURFACE_IMPLICIT_SMOOTHING = float(os.getenv("SURFACE_IMPLICIT_SMOOTHING", "0.1"))
SURFACE_MESHFREE_MAX_MODEL_POINTS = int(os.getenv("SURFACE_MESHFREE_MAX_MODEL_POINTS", "500000"))
SURFACE_MESHFREE_POINT_RADIUS_M = float(os.getenv("SURFACE_MESHFREE_POINT_RADIUS_M", "0.15"))
# 스무딩 파라미터
SURFACE_SMOOTHING_DTM_SIGMA_M = float(os.getenv("SURFACE_SMOOTHING_DTM_SIGMA_M", "0.5"))
@@ -145,11 +102,6 @@ SURFACE_SMOOTHING_DTM_SPLINE_SMOOTH = float(os.getenv("SURFACE_SMOOTHING_DTM_SPL
SURFACE_SMOOTHING_DTM_PREVIEW_RESOLUTION_M = float(
os.getenv("SURFACE_SMOOTHING_DTM_PREVIEW_RESOLUTION_M", "0.5")
)
SURFACE_SMOOTHING_TIN_TAUBIN_ITERATIONS = int(
os.getenv("SURFACE_SMOOTHING_TIN_TAUBIN_ITERATIONS", "10")
)
SURFACE_SMOOTHING_TIN_TAUBIN_LAMBDA = float(os.getenv("SURFACE_SMOOTHING_TIN_TAUBIN_LAMBDA", "0.5"))
SURFACE_SMOOTHING_TIN_TAUBIN_MU = float(os.getenv("SURFACE_SMOOTHING_TIN_TAUBIN_MU", "-0.53"))
# 등고선 파라미터
SURFACE_CONTOUR_INTERVAL_M = float(os.getenv("SURFACE_CONTOUR_INTERVAL_M", "1.0"))
@@ -170,7 +122,7 @@ SHEET_SURFACE_METHODS = [
method.strip()
for method in os.getenv(
"SHEET_SURFACE_METHODS",
"tin_sheet,tin,biharmonic,anudem,multires,laplace",
"laplace",
).split(",")
if method.strip()
]
@@ -209,23 +161,14 @@ def build_surface_model_config() -> dict:
"hole_fill": "laplace",
"tile_size_meters": SURFACE_TILE_SIZE_M,
"max_preview_vertices": SURFACE_MAX_PREVIEW_VERTICES,
"tin_max_input_points": SURFACE_TIN_MAX_INPUT_POINTS,
"dtm_grid_resolution_meters": SURFACE_DTM_GRID_RESOLUTION_M,
"nurbs_degree": SURFACE_NURBS_DEGREE,
"nurbs_patch_size_meters": SURFACE_NURBS_PATCH_SIZE_M,
"nurbs_control_points_per_axis": SURFACE_NURBS_CONTROL_POINTS_PER_AXIS,
"implicit_point_ratio": SURFACE_IMPLICIT_POINT_RATIO,
"implicit_min_points": SURFACE_IMPLICIT_MIN_POINTS,
"implicit_max_points": SURFACE_IMPLICIT_MAX_POINTS,
"implicit_smoothing": SURFACE_IMPLICIT_SMOOTHING,
"meshfree_max_model_points": SURFACE_MESHFREE_MAX_MODEL_POINTS,
"meshfree_point_radius_meters": SURFACE_MESHFREE_POINT_RADIUS_M,
"smoothing_dtm_sigma_meters": SURFACE_SMOOTHING_DTM_SIGMA_M,
"smoothing_dtm_spline_smooth": SURFACE_SMOOTHING_DTM_SPLINE_SMOOTH,
"smoothing_dtm_preview_resolution_meters": SURFACE_SMOOTHING_DTM_PREVIEW_RESOLUTION_M,
"smoothing_tin_taubin_iterations": SURFACE_SMOOTHING_TIN_TAUBIN_ITERATIONS,
"smoothing_tin_taubin_lambda": SURFACE_SMOOTHING_TIN_TAUBIN_LAMBDA,
"smoothing_tin_taubin_mu": SURFACE_SMOOTHING_TIN_TAUBIN_MU,
"contour_interval_meters": SURFACE_CONTOUR_INTERVAL_M,
"contour_grid_resolution_meters": SURFACE_CONTOUR_GRID_RESOLUTION_M,
}
+1 -1
View File
@@ -75,7 +75,7 @@ const LATEST = {
{ x: 0, y: 0 },
{ x: 100, y: 0 },
],
surface_params: { source_filter: "csf", method: "tin", smooth: true, contour_interval_m: 1 },
surface_params: { source_filter: "ptd", method: "implicit", smooth: true, contour_interval_m: 1 },
route_params: {},
};
@@ -24,7 +24,7 @@ from common_util.common_util_surface_confirmation import (
from common_util.common_util_surface_sampler import build_surface_sampler, check_surface_keys
ROOT = Path(__file__).resolve().parents[2]
LAS_PARAMS = {"source_filter": "csf", "method": "dtm", "smooth": True, "contour_interval_m": 1.0}
LAS_PARAMS = {"source_filter": "ptd", "method": "dtm", "smooth": True, "contour_interval_m": 1.0}
SHEET_PARAMS = {
"source_filter": SHEET_ANALYSIS_FILTER,
"method": "dtm",
@@ -51,7 +51,7 @@ def _project(tmp_path: Path, *, las: bool, sheet: bool) -> Path:
models = tmp_path / "B04_PreProcess" / "models"
models.mkdir(parents=True)
if las:
_dtm(models / "dtm_csf_smooth.npz", LAS_Z)
_dtm(models / "dtm_ptd_smooth.npz", LAS_Z)
if sheet: # B04 처럼 원본 · 평활 둘 다(평활 쪽을 읽음)
_dtm(models / f"dtm_{SHEET_ANALYSIS_FILTER}.npz", SHEET_Z - 1.0)
_dtm(models / f"dtm_{SHEET_ANALYSIS_FILTER}_smooth.npz", SHEET_Z)
@@ -40,9 +40,9 @@ def runs() -> dict:
caught: dict = {}
real = engine.build_pool
def spy(graph, follow, max_grade, min_radius_m, measure, snap, **kw):
def spy(graph, follow, max_grade, min_radius_m, measure, **kw):
caught["measure"] = measure
caught["pool"] = real(graph, follow, max_grade, min_radius_m, measure, snap, **kw)
caught["pool"] = real(graph, follow, max_grade, min_radius_m, measure, **kw)
return caught["pool"]
engine._POOLS.clear()
@@ -102,8 +102,8 @@ def _plant_like_b06(root: Path) -> None:
stations = sections.generate_irregular_sections(
root,
ROUTE,
"csf",
"tin",
"ptd",
"implicit",
False,
extra_stations=sections.resolve_extra_stations(root, pipes),
)
@@ -113,8 +113,8 @@ def _plant_like_b06(root: Path) -> None:
def _save(root: Path, screen: list[tuple[float, str]]) -> None:
"""B05 [저장] — 화면은 관 · A군 목록만 실어 보냄."""
request = RouteConfirmRequest(
filter_key="csf",
method="tin",
filter_key="ptd",
method="implicit",
surface_model_id=1,
irregular_stations=[{"chainage_m": c, "structure": s} for c, s in screen],
)
+2 -2
View File
@@ -53,8 +53,8 @@ def test_종단_보기를_열든_안_열든_같은_요청(got):
assert got["closed"]["lines"] == got["open"]["lines"]
body = _confirm_body(got["closed"]["lines"])
assert body == {
"filter_key": "csf",
"method": "tin",
"filter_key": "ptd",
"method": "implicit",
"smooth": True,
"surface_model_id": 7,
"irregular_stations": [],
@@ -60,12 +60,12 @@ def test_미리_뜨기는_안_쓰고_저장과_같은_횡단(tmp_path: Path, mon
monkeypatch.setattr(sections, "build_surface_sampler", lambda *_a, **_k: _Slope())
before = _files(tmp_path)
stations, crosses = sections.sample_structure_sections(
tmp_path, ROUTE, "csf", "tin", False, extra_stations=((130.0, ""),)
tmp_path, ROUTE, "ptd", "implicit", False, extra_stations=((130.0, ""),)
)
assert _files(tmp_path) == before # 파일 · 폴더 하나도 안 생김
assert [c["chainage_m"] for c in crosses] == [130.0] and len(stations) == 1
sections.generate_irregular_sections(
tmp_path, ROUTE, "csf", "tin", False, extra_stations=((130.0, "파형강관 D800"),)
tmp_path, ROUTE, "ptd", "implicit", False, extra_stations=((130.0, "파형강관 D800"),)
)
written = json.loads(
(tmp_path / "B06_Section/cross_sections/cross_00130m.json").read_text(encoding="utf-8")
@@ -235,7 +235,7 @@ def test_같은_입력_두_번이면_Node_0_쓰기_0(node, monkeypatch, tmp_path
return copy.deepcopy(stored)
if fn is prebuild.get_project_storage_relative_path:
return "p"
return {"source_filter": "csf", "method": "dtm", "smooth": False}
return {"source_filter": "ptd", "method": "dtm", "smooth": False}
async def write_haul(connection, *, route_id, data_patch):
written_haul.append(data_patch)
@@ -154,8 +154,8 @@ class _NoXyz(dict):
def test_저장본이_맞으면_점군을_풀지_않음() -> None:
source = _NoXyz(classification=np.zeros(3, np.uint8), bounds=np.zeros((3, 2)))
data = _LazyStructured(source, _LazyXyz(source, 3))
summary = summarize_masks(data, {"csf": np.array([True, False, True])})
assert summary["csf"]["total_point_count"] == 3
summary = summarize_masks(data, {"ptd": np.array([True, False, True])})
assert summary["ptd"]["total_point_count"] == 3
assert "classification" in data and "rgb" not in data
assert len(data["classification"]) == 3
+21 -19
View File
@@ -54,30 +54,30 @@ def sheet_project(tmp_path: Path, monkeypatch: pytest.MonkeyPatch) -> dict:
encoding="utf-8",
)
calls: list[str] = []
real = SHEET_METHOD_BUILDERS["tin_sheet"]
real = SHEET_METHOD_BUILDERS["laplace"]
def counting(*args, **kwargs): # type: ignore[no-untyped-def]
calls.append("tin_sheet")
calls.append("laplace")
return real(*args, **kwargs)
monkeypatch.setitem(SHEET_METHOD_BUILDERS, "tin_sheet", counting)
monkeypatch.setitem(SHEET_METHOD_BUILDERS, "laplace", counting)
return {"root": tmp_path, "processed": processed, "models": models, "calls": calls}
def _build(project: dict, route: np.ndarray = ROUTE) -> list[dict]:
return sheet.build_sheet_surface_model(
project["root"], project["processed"], project["models"], route, EPSG, ["tin_sheet"]
project["root"], project["processed"], project["models"], route, EPSG, ["laplace"]
)
def test_도엽_3D_면_같은_입력_두번째는_안_셈(sheet_project: dict) -> None:
first = _build(sheet_project)
assert sheet_project["calls"] == ["tin_sheet"]
assert (sheet_project["models"] / "dtm_sheet_tin_sheet.key.json").is_file()
assert sheet_project["calls"] == ["laplace"]
assert (sheet_project["models"] / "dtm_sheet_laplace.key.json").is_file()
started = time.perf_counter()
second = _build(sheet_project)
assert time.perf_counter() - started < 0.5
assert sheet_project["calls"] == ["tin_sheet"]
assert sheet_project["calls"] == ["laplace"]
assert second == first
assert first[0]["generation_params"]["file_stamp"]
@@ -97,11 +97,11 @@ def test_도엽_3D_면_열쇠가_바뀌면_다시_셈(
elif change == "algo":
monkeypatch.setattr(sheet, "SHEET_SURFACE_ALGO", sheet.SHEET_SURFACE_ALGO + 1)
elif change == "key_gone":
(sheet_project["models"] / "dtm_sheet_tin_sheet.key.json").unlink()
(sheet_project["models"] / "dtm_sheet_laplace.key.json").unlink()
else:
(sheet_project["models"] / "dtm_sheet_tin_sheet_preview.glb").unlink()
(sheet_project["models"] / "dtm_sheet_laplace_preview.glb").unlink()
_build(sheet_project, route)
assert sheet_project["calls"] == ["tin_sheet", "tin_sheet"]
assert sheet_project["calls"] == ["laplace", "laplace"]
# ── 지표면 모델 id 유지(R2) ─────────────────────────────────────────────────────
@@ -173,12 +173,12 @@ def _model(path: str, stamp: str) -> dict:
"model_type": "dtm",
"model_file_path": path,
"resolution_m": 1.0,
"generation_params": {"source_filter": "csf", "file_stamp": stamp},
"generation_params": {"source_filter": "ptd", "file_stamp": stamp},
"layers": [
{
"layer_name": "dtm_csf_preview",
"layer_name": "dtm_ptd_preview",
"geometry_type": "MESH",
"file_path": "B04_PreProcess/models/dtm_csf_preview.glb",
"file_path": "B04_PreProcess/models/dtm_ptd_preview.glb",
"file_format": "glb",
}
],
@@ -212,7 +212,7 @@ def _save(db: _FakeDb, models: list[dict]) -> list[int]:
},
"models": models,
},
source_filters=["csf"],
source_filters=["ptd"],
)
)
finally:
@@ -221,8 +221,10 @@ def _save(db: _FakeDb, models: list[dict]) -> list[int]:
def test_파일_지문이_같으면_모델_id_유지() -> None:
db = _FakeDb()
path = "B04_PreProcess/models/dtm_csf.npz"
first = _save(db, [_model(path, "1:2"), _model("B04_PreProcess/models/tin_csf.npz", "3:4")])
path = "B04_PreProcess/models/dtm_ptd.npz"
first = _save(
db, [_model(path, "1:2"), _model("B04_PreProcess/models/implicit_ptd.npz", "3:4")]
)
second = _save(db, [_model(path, "1:2")])
assert second == [first[0]]
# 이번 결과에 없는 옛 행은 지움 · 같은 파일 행은 하나
@@ -235,8 +237,8 @@ def test_도엽_한_방식_요청도_id_유지() -> None:
from B04_PreProcess import B04_PreProcess_Repository as repo
db = _FakeDb()
path = "B04_PreProcess/models/dtm_sheet_tin_sheet.npz"
entry = {**_model(path, "5:6"), "generation_params": {"source_filter": "sheet_tin_sheet"}}
path = "B04_PreProcess/models/dtm_sheet_laplace.npz"
entry = {**_model(path, "5:6"), "generation_params": {"source_filter": "sheet_laplace"}}
entry["generation_params"]["file_stamp"] = "5:6"
first = asyncio.run(repo.save_sheet_surface_models(db, "p", [entry])) # type: ignore[arg-type]
second = asyncio.run(repo.save_sheet_surface_models(db, "p", [entry])) # type: ignore[arg-type]
@@ -265,7 +267,7 @@ def test_비용면_계류_메모리_기억(tmp_path: Path, monkeypatch: pytest.M
monkeypatch.setattr(solver, "_cost_surface_signature", lambda *_a: signature["value"])
monkeypatch.setattr(watershed, "_load_features", lambda *_a: [])
monkeypatch.setattr(crs, "resolve_project_crs", lambda _root: "EPSG:5186")
selection = {"source_filter": "csf", "method": "dtm", "smooth": False}
selection = {"source_filter": "ptd", "method": "dtm", "smooth": False}
terrain, _ = engine.load_terrain_and_streams(tmp_path, selection)
again, _ = engine.load_terrain_and_streams(tmp_path, selection)
assert loads == [1] and again is terrain and terrain.signature
+5 -5
View File
@@ -62,14 +62,14 @@ def test_고정_표본은_JSON_길과_같은_점(tmp_path: Path, monkeypatch: py
def test_원천이_그대로면_다시_안_만듦(tmp_path: Path) -> None:
source = tmp_path / "ground_points_csf.npz"
source = tmp_path / "ground_points_ptd.npz"
_npz(source, 30, rgb=False)
target = preview.ensure_point_preview(source, "csf")
assert target.name == "preview_points_csf.bin"
target = preview.ensure_point_preview(source, "ptd")
assert target.name == "preview_points_ptd.bin"
stamp = target.stat().st_mtime_ns
assert preview.ensure_point_preview(source, "csf").stat().st_mtime_ns == stamp
assert preview.ensure_point_preview(source, "ptd").stat().st_mtime_ns == stamp
_npz(source, 40, rgb=False)
head, _, rgb = _decode(preview.ensure_point_preview(source, "csf"))
head, _, rgb = _decode(preview.ensure_point_preview(source, "ptd"))
assert head["sampled_count"] == 40 and rgb is None
@@ -18,7 +18,7 @@ from B05_Profile import B05_Profile_Router_Terrain as terrain_router
from common_util.common_util_surface_sampler import dtm_grid_bytes
PID = "1bb0f876-983f-4ca5-9cc7-f59ee1f77e0a"
PARAMS = {"source_filter": "csf", "method": "dtm", "smooth": True, "contour_interval_m": 1.0}
PARAMS = {"source_filter": "ptd", "method": "dtm", "smooth": True, "contour_interval_m": 1.0}
def _dtm(path: Path, z: float) -> None:
@@ -35,11 +35,11 @@ def _dtm(path: Path, z: float) -> None:
def _project(tmp_path: Path) -> Path:
models = tmp_path / "B04_PreProcess" / "models"
models.mkdir(parents=True)
_dtm(models / "dtm_csf.npz", 100.0)
_dtm(models / "dtm_csf_smooth.npz", 101.0)
_dtm(models / "dtm_ptd.npz", 100.0)
_dtm(models / "dtm_ptd_smooth.npz", 101.0)
manifest = {
"config_signature": "cfg1",
"source_filters": {"csf": {"methods": {"dtm": {"smooth": {"smoothing_signature": "sm1"}}}}},
"source_filters": {"ptd": {"methods": {"dtm": {"smooth": {"smoothing_signature": "sm1"}}}}},
}
(models / "manifest.json").write_text(json.dumps(manifest), encoding="utf-8")
return tmp_path
@@ -83,19 +83,19 @@ def _grid(smooth: bool, etag: str | None = None):
def test_평활_격자를_저장본으로_내고_같은_ETag_면_304(tmp_path, monkeypatch):
root = _project(tmp_path)
models = root / "B04_PreProcess" / "models"
_patch(monkeypatch, root, ("dtm", "B04_PreProcess/models/dtm_csf.npz"))
_patch(monkeypatch, root, ("dtm", "B04_PreProcess/models/dtm_ptd.npz"))
response = _grid(True)
assert response.status_code == 200
saved = models / "dtm_csf_smooth_grid.bin"
saved = models / "dtm_ptd_smooth_grid.bin"
assert Path(response.path) == saved
assert saved.read_bytes() == dtm_grid_bytes(models / "dtm_csf_smooth.npz")
assert saved.read_bytes() == dtm_grid_bytes(models / "dtm_ptd_smooth.npz")
assert _grid(True, response.headers["etag"]).status_code == 304
assert Path(_grid(False).path).name == "dtm_csf_grid.bin"
assert Path(_grid(False).path).name == "dtm_ptd_grid.bin"
def test_DTM_아닌_표현은_404(tmp_path, monkeypatch):
root = _project(tmp_path)
_patch(monkeypatch, root, ("tin", "B04_PreProcess/models/tin_csf.npz"))
_patch(monkeypatch, root, ("implicit", "B04_PreProcess/models/implicit_ptd.npz"))
assert _grid(True).status_code == 404
_patch(monkeypatch, root, None)
assert _grid(True).status_code == 404
@@ -103,16 +103,19 @@ def test_DTM_아닌_표현은_404(tmp_path, monkeypatch):
def test_받을_자리_서명은_manifest_구성과_평활(tmp_path):
models = _project(tmp_path) / "B04_PreProcess" / "models"
confirmed = (7, "dtm", "B04_PreProcess/models/dtm_csf.npz")
confirmed = (7, "dtm", "B04_PreProcess/models/dtm_ptd.npz")
grid = terrain_router._design_grid(PARAMS, models, confirmed)
assert grid == {"model_id": 7, "smooth": True, "signature": "cfg1.sm1"}
raw = terrain_router._design_grid({**PARAMS, "smooth": False}, models, confirmed)
assert raw["signature"] == "cfg1."
# 서버 표본기와 다른 파일 · DTM 아닌 표현 · 확정 없음 = 격자 안 냄(서버 표본 길)
assert (
terrain_router._design_grid({**PARAMS, "source_filter": "pmf"}, models, confirmed) is None
terrain_router._design_grid(
{**PARAMS, "source_filter": "classification"}, models, confirmed
)
is None
)
assert terrain_router._design_grid({**PARAMS, "method": "tin"}, models, confirmed) is None
assert terrain_router._design_grid({**PARAMS, "method": "implicit"}, models, confirmed) is None
assert terrain_router._design_grid(PARAMS, models, None) is None
@@ -122,7 +125,7 @@ def test_지반고_자릿수_기본_3_묶음은_6(tmp_path, monkeypatch):
models.mkdir(parents=True)
x = np.arange(0.0, 11.0)
np.savez(
models / "dtm_csf_smooth.npz",
models / "dtm_ptd_smooth.npz",
x=x,
y=x,
z=np.add.outer(x * 0.0, x * 0.1234567),
@@ -144,10 +147,10 @@ def test_지반고_자릿수_기본_3_묶음은_6(tmp_path, monkeypatch):
def test_격자_모델은_params_파일_확정_행이_어긋나도(monkeypatch):
"""578-1 — 확정 행(sheet)이 params(csf)와 어긋나도 서버 표본기와 같은 `dtm_{filter}.npz` 행을 고름."""
"""578-1 — 확정 행(sheet)이 params(ptd)와 어긋나도 서버 표본기와 같은 `dtm_{filter}.npz` 행을 고름."""
rows = [
(9, "dtm", "B04_PreProcess/models/dtm_sheet_laplace.npz"), # 확정(어긋남)
(7, "dtm", "B04_PreProcess/models/dtm_csf.npz"),
(7, "dtm", "B04_PreProcess/models/dtm_ptd.npz"),
]
class _Cursor:
@@ -163,6 +166,6 @@ def test_격자_모델은_params_파일_확정_행이_어긋나도(monkeypatch):
yield _Cursor()
pick = asyncio.run(terrain_router._design_model(_Connection(), PID, PARAMS))
assert pick == (7, "dtm", "B04_PreProcess/models/dtm_csf.npz")
other = {**PARAMS, "source_filter": "pmf"}
assert pick == (7, "dtm", "B04_PreProcess/models/dtm_ptd.npz")
other = {**PARAMS, "source_filter": "classification"}
assert asyncio.run(terrain_router._design_model(_Connection(), PID, other)) is None
+1 -1
View File
@@ -331,7 +331,7 @@ def chain_world(monkeypatch, tmp_path):
chain,
"redesign_solve_request",
lambda *a: b05_router.RouteSolveRequest(
filter_key="csf",
filter_key="ptd",
method="dtm",
smooth=False,
algorithm="as_planned",
@@ -125,7 +125,7 @@ def _source_project(root: Path) -> None:
def test_원본에서_다시_읽기(tmp_path: Path) -> None:
_source_project(tmp_path)
surface = {"source_filter": "csf", "method": "dtm", "smooth": False}
surface = {"source_filter": "ptd", "method": "dtm", "smooth": False}
stored = load_design_route(tmp_path, surface, (20.0, 80.0))
assert [v.y for v in stored.vertices] == [5.0, 5.0] # 정본이 먼저
fresh = load_design_route(tmp_path, surface, (20.0, 80.0), from_source=True)
@@ -175,7 +175,7 @@ def test_사용_범위_바꾸면_예상노선_다시_쓰고_미리_셈(
return "x"
async def surface_params(*_):
return {"source_filter": "csf", "method": "dtm", "smooth": False}
return {"source_filter": "ptd", "method": "dtm", "smooth": False}
async def precompute(_pid):
calls.append("precompute")
@@ -122,8 +122,8 @@ def hill_runs() -> dict:
@pytest.mark.parametrize("mode", ["follow", "gentle", "earthwork", "weighted"])
def test_recommend_result_on_stations(hill_runs: dict, mode: str) -> None:
"""추적(평면계획선 길)만 측점 맞춤 · 추천 갈래는 계산 로직 결과 그대로(769-1 · 요구 611)."""
result = hill_runs[mode]
assert result["off_station_curves"] == 0
outline = build_planned_polyline(
[tuple(n[k] for k in ("x", "y")) for n in result["nodes"]],
min_radius_m=12.0,
@@ -133,16 +133,18 @@ def test_recommend_result_on_stations(hill_runs: dict, mode: str) -> None:
)
# 화면이 펴서 다시 짠 선 = 서버 선
assert np.allclose(outline.vertices, result["planned"], atol=1e-3)
assert station_hits(outline) == len(outline.curves) > 0
if mode != "follow": # 갈래는 기울기 상한 안(언덕은 비켜 갈 수 있음)
if mode == "follow":
assert result["off_station_curves"] == 0
assert station_hits(outline) == len(outline.curves) > 0
# 갈래는 기울기 상한 안(언덕은 비켜 갈 수 있음) · 비중은 제약을 벌점으로만(넘으면 목록)
if mode in ("gentle", "earthwork"):
assert result["metrics"]["over_limit"]["total"] == 0
assert result["metrics"]["max_grade_pct"] <= 14.0 + 1e-6
# 측점 맞춤 필수 — 맞추기 전 선은 안 냄(419)
assert "unsnapped" not in result
def test_합성_지형_modes_before_after(hill_runs: dict) -> None:
"""합성 언덕(14 % · R12) — 갈래별 최대 기울기와 측점 밖 곡선을 함께 지킨다."""
"""합성 언덕(14 % · R12) — 추적만 측점 위 · 추천 갈래는 측점 맞춤 없음(769-1)."""
on_total, curves_total = 0, 0
for mode, result in hill_runs.items():
after = result["metrics"]
@@ -155,7 +157,8 @@ def test_합성_지형_modes_before_after(hill_runs: dict) -> None:
f" · 곡선 {curves}(측점 밖 {result['off_station_curves']})"
)
print(f"측점 위 곡선 {on_total}/{curves_total} = {100 * on_total / curves_total:.0f}%")
assert curves_total and on_total == curves_total
assert curves_total and hill_runs["follow"]["off_station_curves"] == 0
assert on_total < curves_total # 추천 갈래는 맞추지 않음
def test_합성_예상노선_eight_outlines() -> None:
@@ -203,13 +203,12 @@ def hill_runs() -> dict:
@pytest.mark.parametrize("mode", ["gradeline", "earthwork", "gentle", "weighted"])
def test_hill_branches_rules(hill_runs: dict, mode: str) -> None:
"""언덕 — 꺾임만 0 · 측점 밖 0 · 못 지킨 곳 수 = R · 계류 · 측점 · 맞추기 전 선 없음."""
"""언덕 — 꺾임만 0 · 못 지킨 곳 수 = R · 계류 · 맞추기 전 선 없음(측점 맞춤 없음 · 769-1)."""
result = hill_runs[mode]
interior = result["nodes"][1:-1]
assert interior and all(node["radius_m"] is not None for node in interior) # 삼각 0
assert result["off_station_curves"] == 0
over = result["metrics"]["over_limit"]
assert over["total"] == over["radius"] + over["stream"] + over.get("station", 0)
assert over["total"] == over["radius"] + over["stream"]
assert "unsnapped" not in result
assert result["ends_over_limit"] is False and result["ends_grade_pct"] < 14.0
+1 -1
View File
@@ -211,7 +211,7 @@ def test_E1_재설계_요청이_설정을_옮김():
"enforce_pipe_clearance": True,
"station_display": {"station": 12, "cumulative": 240},
}
request = redesign_solve_request(params, {"source_filter": "csf"}, 5)
request = redesign_solve_request(params, {"source_filter": "ptd"}, 5)
assert request.design_speed_kph == 30
assert request.enforce_pipe_clearance is True
assert request.station_display == {"station": 12, "cumulative": 240}
+2 -2
View File
@@ -167,8 +167,8 @@ async def _aboom(*_a, **_k):
_REQUEST = RouteConfirmRequest(
filter_key="csf",
method="tin",
filter_key="ptd",
method="implicit",
surface_model_id=1,
uphill_overrides=[{"chainage_m": 40.0, "side": "right"}],
)
+6 -6
View File
@@ -21,7 +21,7 @@ from common_util.common_util_surface_sampler import build_surface_sampler, check
from M02_MasterTemplete import M02_Item_Snapshot as snap
BAD_IDS = ["../../..", "..", "a/b", "a\\b", "C:/Windows", "/etc", ".hidden", "", "x" * 65]
BAD_KEYS = ["../other/dtm_csf", "csf/../../x", "C:/x", "", "CSF"]
BAD_KEYS = ["../other/dtm_ptd", "ptd/../../x", "C:/x", "", "PTD", "csf"]
def _structure(structure_id):
@@ -80,7 +80,7 @@ def test_H2_스키마_허용_목록(bad):
with pytest.raises(ValidationError):
RouteSolveRequest(filter_key=bad, **body)
with pytest.raises(ValidationError):
RouteSolveRequest(filter_key="csf", method=bad, **body)
RouteSolveRequest(filter_key="ptd", method=bad, **body)
with pytest.raises(ValidationError):
RouteConfirmRequest(filter_key=bad)
with pytest.raises(ValidationError):
@@ -88,7 +88,7 @@ def test_H2_스키마_허용_목록(bad):
def test_H2_스키마_정상값_통과():
assert RouteConfirmRequest(filter_key="grid_min_z", method="tin").method == "tin"
assert RouteConfirmRequest(filter_key="ptd_refine", method="implicit").method == "implicit"
assert RouteConfirmRequest().filter_key is None
@@ -97,13 +97,13 @@ def test_H2_파일_이름_만드는_곳(bad, tmp_path):
with pytest.raises(ValueError):
check_surface_keys(bad)
with pytest.raises(ValueError):
check_surface_keys("csf", bad)
check_surface_keys("ptd", bad)
with pytest.raises(ValueError):
build_surface_sampler(tmp_path, bad, "dtm", False)
with pytest.raises(ValueError):
_load_dtm_grid(tmp_path, bad, False)
with pytest.raises(ValueError):
_source_npz_paths(tmp_path, "csf", bad, False)
_source_npz_paths(tmp_path, "ptd", bad, False)
with pytest.raises(ValueError):
load_or_build_skeleton(tmp_path, bad, "dtm", False)
check_surface_keys("classification", "meshfree")
check_surface_keys("classification", "implicit")
@@ -58,7 +58,7 @@ class FakeWfs:
@pytest.fixture
def world(tmp_path, monkeypatch):
"""D 프로젝트(수치지형도) · B 프로젝트(라이다 csf 확정 + 수치지형도 서피스 · B04 지적도)."""
"""D 프로젝트(수치지형도) · B 프로젝트(라이다 ptd 확정 + 수치지형도 서피스 · B04 지적도)."""
import common_util.common_util_crs as crs
import common_util.common_util_surface_confirmation as confirm
import config.config_db as config_db
@@ -66,8 +66,8 @@ def world(tmp_path, monkeypatch):
d_root, b_root = tmp_path / "D", tmp_path / "B"
_dem(d_root / "B04_PreProcess" / "models_plan" / "dtm_sheet_laplace.npz", 0.3)
b_models = b_root / "B04_PreProcess" / "models"
_dem(b_models / "dtm_csf.npz", 0.6)
_dem(b_models / "dtm_csf_smooth.npz", 0.5)
_dem(b_models / "dtm_ptd.npz", 0.6)
_dem(b_models / "dtm_ptd_smooth.npz", 0.5)
_dem(b_models / "dtm_sheet_laplace_smooth.npz", 0.2)
processed = b_root / "B04_PreProcess" / "processed"
processed.mkdir(parents=True)
@@ -76,7 +76,7 @@ def world(tmp_path, monkeypatch):
(processed / "vworld_white_meta.json").write_text(json.dumps(meta), "utf-8")
kinds = {"D": "PLAN", "B": "ROAD"}
state = {"pid": "D"}
params = {"source_filter": "csf", "method": "dtm", "smooth": True}
params = {"source_filter": "ptd", "method": "dtm", "smooth": True}
async def project_row(project_id):
return {"project_kind": kinds[state["pid"]], "storage_path": state["pid"]}
@@ -183,7 +183,7 @@ def test_lidar_zero_order_whole_grid(world) -> None:
n, cell = 200, 1.0
x, y = X0 + np.arange(n) * cell, Y0 + np.arange(n) * cell
valley = np.abs(x - x.mean())[None, :] * 0.3 + (y - Y0)[:, None] * 0.2
npz = world["B"] / "B04_PreProcess" / "models" / "dtm_csf_smooth.npz"
npz = world["B"] / "B04_PreProcess" / "models" / "dtm_ptd_smooth.npz"
np.savez(npz, x=x, y=y, z=valley, valid_mask=np.ones(valley.shape, bool), resolution=[cell])
res = world["get"]("B", "zero-order", terrain="lidar")
got = res.json()
+10 -9
View File
@@ -183,23 +183,24 @@ def test_갈음_안_지은_필터_파일_지움(tmp_path: Path) -> None:
"dtm_sheet_laplace.key.json",
]
stale = [
"dtm_csf.npz",
"dtm_csf_preview.glb",
"dtm_csf_smooth.npz",
"dtm_csf_smooth_grid.bin",
"tin_csf.npz",
"contour_csf_dtm_smooth_1.0m.json",
"dtm_ptd.npz",
"dtm_ptd_preview.glb",
"dtm_ptd_smooth.npz",
"dtm_ptd_smooth_grid.bin",
"implicit_ptd.npz",
"contour_ptd_dtm_smooth_1.0m.json",
]
for name in keep + stale:
(models / name).write_bytes(b"x")
for name in ["mask_csf.npy", "ground_points_csf.npz", "preview_points_csf.bin"]:
for name in ["mask_ptd.npy", "ground_points_ptd.npz", "preview_points_ptd.bin"]:
(processed / name).write_bytes(b"x")
(processed / "mask_classification.npy").write_bytes(b"x")
manifest = {"source_filters": {"classification": {}, "csf": {}}, "x": 1}
manifest = {"source_filters": {"classification": {}, "ptd": {}}, "x": 1}
(models / "manifest.json").write_text(json.dumps(manifest), encoding="utf-8")
gone = lidar_upload.drop_stale_filters(tmp_path, ["classification"])
assert gone == sorted(
f"B04_PreProcess/models/{n}" for n in ["dtm_csf.npz", "dtm_csf_smooth.npz", "tin_csf.npz"]
f"B04_PreProcess/models/{n}"
for n in ["dtm_ptd.npz", "dtm_ptd_smooth.npz", "implicit_ptd.npz"]
)
assert sorted(p.name for p in models.iterdir()) == sorted([*keep, "manifest.json"])
assert [p.name for p in processed.iterdir()] == ["mask_classification.npy"]
+14 -14
View File
@@ -52,31 +52,31 @@ def _band_contour(path: Path) -> bool:
def test_implicit_만_지어도_낡은_dtm_같이_다시(tmp_path: Path) -> None:
data, mask = _banded_plane()
masks = {"csf": mask}
build_all_terrain_models(data, masks, tmp_path, _config("csf", ["dtm", "implicit"], 1.0))
narrow = _valid_area(tmp_path / "dtm_csf.npz")
contour = tmp_path / "contour_csf_implicit_1.0m.json"
masks = {"ptd": mask}
build_all_terrain_models(data, masks, tmp_path, _config("ptd", ["dtm", "implicit"], 1.0))
narrow = _valid_area(tmp_path / "dtm_ptd.npz")
contour = tmp_path / "contour_ptd_implicit_1.0m.json"
assert not _band_contour(contour) # 1 m 닫기 = 띠 빈칸 · 106 m 줄 잘림
# 화면 [서피스 생성] = 고른 표현 하나만 — 외곽 설정(지문)이 바뀌었으니 dtm 도 다시
manifest = build_all_terrain_models(data, masks, tmp_path, _config("csf", ["implicit"], 3.0))
methods = manifest["source_filters"]["csf"]["methods"]
manifest = build_all_terrain_models(data, masks, tmp_path, _config("ptd", ["implicit"], 3.0))
methods = manifest["source_filters"]["ptd"]["methods"]
assert methods["dtm"]["status"] == "completed"
wide = _valid_area(tmp_path / "dtm_csf.npz")
grid_area = _valid_area(implicit_grid_npz(tmp_path, "csf"))
wide = _valid_area(tmp_path / "dtm_ptd.npz")
grid_area = _valid_area(implicit_grid_npz(tmp_path, "ptd"))
glb_area = methods["implicit"]["footprint_area_m2"]
assert wide > 1.05 * narrow
assert grid_area == wide # 격자 외곽 = dtm 외곽
assert abs(glb_area - wide) / wide < 0.05 # 3D 외곽 = 격자 외곽
assert _band_contour(contour) # 등고선이 넓힌 띠까지
assert not contour_cache_is_stale(contour, tmp_path / "implicit_csf.npz")
assert not contour_cache_is_stale(contour, tmp_path / "implicit_ptd.npz")
def test_dtm_이_바뀌면_등고선_캐시_낡음(tmp_path: Path) -> None:
model = tmp_path / "implicit_csf.npz"
dtm = tmp_path / "dtm_csf.npz"
grid = tmp_path / "implicit_csf_grid.npz"
contour = tmp_path / "contour_csf_implicit_1.0m.json"
model = tmp_path / "implicit_ptd.npz"
dtm = tmp_path / "dtm_ptd.npz"
grid = tmp_path / "implicit_ptd_grid.npz"
contour = tmp_path / "contour_ptd_implicit_1.0m.json"
assert contour_cache_is_stale(contour, model) # 없음
for index, path in enumerate((model, dtm, grid, contour)):
path.write_bytes(b"x")
@@ -86,7 +86,7 @@ def test_dtm_이_바뀌면_등고선_캐시_낡음(tmp_path: Path) -> None:
os.utime(newer, ns=(2 * 10**18, 2 * 10**18))
assert contour_cache_is_stale(contour, model)
os.utime(newer, ns=(10**18, 10**18))
tin_smooth = tmp_path / "tin_csf_smooth.npz"
tin_smooth = tmp_path / "tin_ptd_smooth.npz"
tin_smooth.write_bytes(b"x")
os.utime(tin_smooth, ns=(10**18, 10**18))
os.utime(dtm, ns=(2 * 10**18, 2 * 10**18))
@@ -135,9 +135,10 @@ def test_implicit_파일_없으면_까닭_그대로(tmp_path) -> None:
def test_깨진_모델은_DTM_대체_없이_오류(tmp_path) -> None:
models = _models(tmp_path)
(models / "tin_classification.npz").write_bytes(b"broken")
(models / "implicit_classification_grid.npz").unlink(missing_ok=True)
(models / "implicit_classification.npz").write_bytes(b"broken")
with pytest.raises(ValueError, match="지표면 모델을 읽지 못함"):
solver._build_cost_surface(models, "classification", "tin", False)
solver._build_cost_surface(models, "classification", "implicit", False)
source = (Path(solver.__file__)).read_text(encoding="utf-8")
assert "return dtm_z\n\n return dtm_z" not in source
assert "except Exception:\n return dtm_z" not in source
@@ -128,7 +128,7 @@ def b04_world(tmp_path, monkeypatch):
def test_분석은_보기만_단계_확정_그대로(b04_world) -> None:
request = SurfaceAnalyzeRequest(input_file_id=3, source_filters=["csf"], methods=["implicit"])
request = SurfaceAnalyzeRequest(input_file_id=3, source_filters=["ptd"], methods=["implicit"])
answer = asyncio.run(b04.analyze_surface(uuid4(), request))
assert answer.surface_model_ids == [1, 2]
sql = " | ".join(b04_world["log"])
@@ -1,180 +0,0 @@
"""실데이터 등가성 — 래스터화가 기존 griddata 경로와 같은 표고를 내는지.
PLAN 2026-08-17 완료 조건 검증.
- `meshfree`: Delaunay를 그대로 쓰므로 표고 완전 일치
- `tin`: 저장 faces와 재Delaunay는 다른 삼각망이라 차이가 남는다. 그 차이가
**보간 버그가 아니라 삼각망 선택 차이**임을 대조로 분리한다.
**주의: griddata가 호출당 3분 걸려 이 파일 전체 실행에 약 8분이 든다.**
실제 프로젝트 산출물이 없으면 skip.
"""
import time
from pathlib import Path
import numpy as np
import pytest
from scipy.interpolate import griddata
from scipy.spatial import Delaunay
from B04_PreProcess.B04_PreProcess_Engine_Contour import (
_grid_axes,
extract_contours,
rasterize_triangle_mesh,
)
MODELS = Path(
"C:/Program_coding/임도설계 및 견적자동화 프로그램 개발/storage/1/3/"
"f9543035-4a91-4df0-b467-cdc515723507/B04_PreProcess/models"
)
def _load(name):
path = MODELS / name
if not path.exists():
pytest.skip(f"실데이터 없음: {path}")
return np.load(path)
def _axes_for(xy):
return _grid_axes(
float(np.min(xy[:, 0])),
float(np.max(xy[:, 0])),
float(np.min(xy[:, 1])),
float(np.max(xy[:, 1])),
1.0,
)
def _compare(old, new):
return {
"only_old": int((np.isfinite(old) & ~np.isfinite(new)).sum()),
"only_new": int((np.isfinite(new) & ~np.isfinite(old)).sum()),
"both": np.isfinite(old) & np.isfinite(new),
}
def test_meshfree_matches_griddata_exactly():
"""meshfree — 같은 Delaunay를 쓰므로 표고가 완전히 일치해야 한다."""
points = _load("meshfree_csf.npz")["points"]
x_coords, y_coords = _axes_for(points)
xx, yy = np.meshgrid(x_coords, y_coords)
triangulation = Delaunay(np.asarray(points[:, :2], dtype=np.float64))
new = rasterize_triangle_mesh(points, triangulation.simplices, x_coords, y_coords)
old = griddata(points[:, :2], points[:, 2], (xx, yy), method="linear")
result = _compare(old, new)
assert result["both"].sum() > 10_000, "비교 대상 셀이 너무 적다"
assert result["only_old"] == 0, "기존이 채우던 셀이 비었다"
assert result["only_new"] == 0, "없던 셀이 생겼다"
assert np.abs(old[result["both"]] - new[result["both"]]).max() < 1e-9
def test_rasterize_matches_griddata_on_same_triangulation():
"""래스터화 정확성 — 같은 재Delaunay를 넣으면 griddata와 완전 일치.
이 대조가 통과하면 tin에서 남는 차이는 전부 삼각망 선택 차이다.
"""
vertices = _load("tin_csf.npz")["vertices"]
x_coords, y_coords = _axes_for(vertices)
xx, yy = np.meshgrid(x_coords, y_coords)
triangulation = Delaunay(np.asarray(vertices[:, :2], dtype=np.float64))
control = rasterize_triangle_mesh(vertices, triangulation.simplices, x_coords, y_coords)
old = griddata(vertices[:, :2], vertices[:, 2], (xx, yy), method="linear")
result = _compare(old, control)
assert result["only_old"] == 0 and result["only_new"] == 0
assert np.abs(old[result["both"]] - control[result["both"]]).max() < 1e-9
def test_tin_difference_comes_only_from_triangulation():
"""tin 저장 faces — griddata 대비 차이가 재Delaunay 대비 차이와 같아야 한다.
두 대조가 일치하면 차이의 원인은 삼각망 하나뿐이다. 저장 faces는 긴 변 제거·
외곽 클리핑이 반영된 실제 TIN이므로(`ModelBuild.py:38-51`) 재Delaunay보다 덮는
면적이 작고, 그 몫이 `only_old`로 나온다.
"""
data = _load("tin_csf.npz")
vertices, faces = data["vertices"], data["faces"]
x_coords, y_coords = _axes_for(vertices)
xx, yy = np.meshgrid(x_coords, y_coords)
stored = rasterize_triangle_mesh(vertices, faces, x_coords, y_coords)
triangulation = Delaunay(np.asarray(vertices[:, :2], dtype=np.float64))
control = rasterize_triangle_mesh(vertices, triangulation.simplices, x_coords, y_coords)
old = griddata(vertices[:, :2], vertices[:, 2], (xx, yy), method="linear")
vs_griddata = _compare(old, stored)
vs_control = _compare(control, stored)
assert vs_griddata["only_old"] == vs_control["only_old"]
assert vs_griddata["only_new"] == vs_control["only_new"] == 0
diff = np.abs(old[vs_griddata["both"]] - stored[vs_griddata["both"]])
changed = diff > 1e-6
print(
f"\n[tin] 공통셀={int(vs_griddata['both'].sum()):,} "
f"저장 faces가 안 덮는 셀={vs_griddata['only_old']:,}\n"
f" 값 다른 셀={int(changed.sum()):,} ({changed.mean():.2%}) "
f"최대={diff.max():.4f}m 차이셀 평균={diff[changed].mean():.4f}m"
)
# 대각선 선택이 다른 사각형에서만 어긋나므로 국소 지형 기복을 넘지 않아야 한다.
assert diff.max() < 1.0, "표고차 1m 초과 — 삼각망 차이로 설명 불가"
@pytest.mark.parametrize(
"npz_name,representation,cached_name,exact",
[
("tin_csf.npz", "triangular_mesh", "contour_csf_tin_1.0m.json", False),
("meshfree_csf.npz", "meshfree_surfels", "contour_csf_meshfree_1.0m.json", True),
("dtm_csf.npz", "regular_grid", "contour_csf_dtm_1.0m.json", True),
],
)
def test_contour_output_matches_cached(npz_name, representation, cached_name, exact):
"""최종 산출물 대조 — 캐시된 등고선 JSON과 레벨·정점 수를 맞춰 본다."""
import json
npz_path = MODELS / npz_name
cached_path = MODELS / cached_name
if not npz_path.exists() or not cached_path.exists():
pytest.skip("실데이터 없음")
started = time.perf_counter()
lines = extract_contours(npz_path, representation, 1.0)
elapsed = time.perf_counter() - started
with open(cached_path, encoding="utf-8") as handle:
cached = json.load(handle)["contours"]
new_levels = {round(float(item["level"]), 3) for item in lines}
old_levels = {round(float(item["level"]), 3) for item in cached}
new_points = sum(len(item["coordinates"]) for item in lines)
old_points = sum(len(item["coordinates"]) for item in cached)
print(
f"\n[{npz_name}] {elapsed:.2f}초 세그먼트 {len(lines):,}/{len(cached):,} "
f"정점 {new_points:,}/{old_points:,}"
)
assert new_levels == old_levels, "등고선 레벨 집합이 바뀌었다"
if exact:
assert len(lines) == len(cached)
assert new_points == old_points
else:
# tin은 삼각망 차이만큼 경계 세그먼트가 미세하게 달라진다.
assert abs(new_points - old_points) / old_points < 0.01
def test_contour_extraction_is_fast():
"""가속이 목적이므로 속도를 회귀로 못박는다 (변경 전 tin 119초·meshfree 302초)."""
for npz_name, representation in (
("tin_csf.npz", "triangular_mesh"),
("meshfree_csf.npz", "meshfree_surfels"),
):
path = MODELS / npz_name
if not path.exists():
pytest.skip("실데이터 없음")
started = time.perf_counter()
extract_contours(path, representation, 1.0)
elapsed = time.perf_counter() - started
print(f"\n[속도] {npz_name} {elapsed:.2f}초")
assert elapsed < 30.0, f"{npz_name} {elapsed:.1f}초 — 가속 실패"
-195
View File
@@ -1,195 +0,0 @@
"""B04 삼각망 래스터화 — 위치탐색 없이 격자 표고를 만드는 함수의 정확도 검증.
기존 경로는 scipy `griddata`가 질의점마다 삼각형 98만개 망을 탐색해 등고선 1회에
2~5분이 걸렸다(PLAN 2026-08-17). 삼각형을 격자에 직접 래스터화하면 탐색이 사라진다.
정밀도 유지가 조건이므로 무게중심 선형보간의 해석해 일치를 먼저 못박는다.
"""
import numpy as np
import pytest
from B04_PreProcess.B04_PreProcess_Engine_Contour import rasterize_triangle_mesh
def _axes(x_min, x_max, y_min, y_max, step=1.0):
cols = int(round((x_max - x_min) / step)) + 1
rows = int(round((y_max - y_min) / step)) + 1
return (
np.linspace(x_min, x_max, cols, dtype=np.float64),
np.linspace(y_min, y_max, rows, dtype=np.float64),
)
def _plane(x, y):
"""검증용 기준 평면 — 선형보간이 정확하면 오차 0이어야 한다."""
return 2.0 * x + 3.0 * y + 5.0
def test_single_triangle_matches_plane_solution():
"""평면 위 삼각형 하나 — 내부 격자점이 해석해와 일치."""
verts = np.array(
[
[0.0, 0.0, _plane(0.0, 0.0)],
[10.0, 0.0, _plane(10.0, 0.0)],
[0.0, 10.0, _plane(0.0, 10.0)],
]
)
tris = np.array([[0, 1, 2]])
x_coords, y_coords = _axes(0.0, 10.0, 0.0, 10.0)
z_grid = rasterize_triangle_mesh(verts, tris, x_coords, y_coords)
xx, yy = np.meshgrid(x_coords, y_coords)
filled = np.isfinite(z_grid)
assert filled.any()
assert np.allclose(z_grid[filled], _plane(xx[filled], yy[filled]), atol=1e-9)
def test_cells_outside_mesh_stay_nan():
"""삼각망 밖은 NaN — footprint 마스크가 기존과 같이 동작해야 한다."""
verts = np.array([[0.0, 0.0, 1.0], [2.0, 0.0, 1.0], [0.0, 2.0, 1.0]])
tris = np.array([[0, 1, 2]])
x_coords, y_coords = _axes(0.0, 10.0, 0.0, 10.0)
z_grid = rasterize_triangle_mesh(verts, tris, x_coords, y_coords)
# 빗변 바깥쪽 (9, 9)는 삼각형 밖
assert np.isnan(z_grid[9, 9])
# 꼭짓점 (0, 0)은 삼각형 안
assert z_grid[0, 0] == pytest.approx(1.0)
def test_two_triangles_fill_full_square():
"""사각형을 이루는 두 삼각형 — 격자 전체가 빈칸 없이 채워진다."""
corners = [(0.0, 0.0), (4.0, 0.0), (4.0, 4.0), (0.0, 4.0)]
verts = np.array([[x, y, _plane(x, y)] for x, y in corners])
tris = np.array([[0, 1, 2], [0, 2, 3]])
x_coords, y_coords = _axes(0.0, 4.0, 0.0, 4.0)
z_grid = rasterize_triangle_mesh(verts, tris, x_coords, y_coords)
assert np.isfinite(z_grid).all()
xx, yy = np.meshgrid(x_coords, y_coords)
assert np.allclose(z_grid, _plane(xx, yy), atol=1e-9)
def test_shared_edge_is_consistent():
"""두 삼각형이 공유하는 변 위의 격자점은 어느 쪽으로 계산해도 같은 값."""
verts = np.array(
[[0.0, 0.0, 0.0], [4.0, 0.0, 8.0], [4.0, 4.0, 20.0], [0.0, 4.0, 12.0]]
) # z = 2x + 3y
tris = np.array([[0, 1, 2], [0, 2, 3]])
x_coords, y_coords = _axes(0.0, 4.0, 0.0, 4.0)
z_grid = rasterize_triangle_mesh(verts, tris, x_coords, y_coords)
xx, yy = np.meshgrid(x_coords, y_coords)
assert np.allclose(z_grid, 2.0 * xx + 3.0 * yy, atol=1e-9)
def test_degenerate_triangle_is_skipped():
"""면적 0 삼각형이 섞여 있어도 죽지 않고 나머지를 채운다."""
verts = np.array(
[
[0.0, 0.0, 1.0],
[4.0, 0.0, 1.0],
[0.0, 4.0, 1.0],
[1.0, 1.0, 99.0],
[2.0, 2.0, 99.0],
[3.0, 3.0, 99.0],
]
)
tris = np.array([[0, 1, 2], [3, 4, 5]]) # 두 번째는 일직선
x_coords, y_coords = _axes(0.0, 4.0, 0.0, 4.0)
z_grid = rasterize_triangle_mesh(verts, tris, x_coords, y_coords)
filled = np.isfinite(z_grid)
assert filled.any()
assert np.allclose(z_grid[filled], 1.0)
def test_triangle_smaller_than_cell_still_lands():
"""격자 간격보다 작은 삼각형도 자기가 덮는 격자점을 채운다 (실데이터는 m²당 13.7면)."""
verts = np.array([[1.6, 1.6, 7.0], [2.4, 1.6, 7.0], [2.0, 2.4, 7.0]])
tris = np.array([[0, 1, 2]])
x_coords, y_coords = _axes(0.0, 4.0, 0.0, 4.0)
z_grid = rasterize_triangle_mesh(verts, tris, x_coords, y_coords)
assert z_grid[2, 2] == pytest.approx(7.0)
assert np.isnan(z_grid[0, 0])
def test_float32_utm_axis_does_not_drop_cells():
"""float32 축(UTM 좌표) — 등간격 환산으로 인덱스를 내면 셀을 놓친다.
`_grid_axes`가 float32를 쓰는데 18만대 좌표에서 해상도가 0.015625m라 간격이
0.984~1.0으로 흔들린다. 등간격 가정 시 최대 0.81셀 어긋나 실측 8,819셀이
비었다(2026-08-17). 축을 직접 탐색해야 한다.
"""
x_coords = np.linspace(183433.6, 183805.8, 374, dtype=np.float32)
y_coords = np.linspace(489168.8, 489570.2, 403, dtype=np.float32)
# 격자 전체를 덮는 큰 삼각형 2개 → 모든 셀이 채워져야 한다.
x0, x1 = float(x_coords[0]) - 5.0, float(x_coords[-1]) + 5.0
y0, y1 = float(y_coords[0]) - 5.0, float(y_coords[-1]) + 5.0
corners = [(x0, y0), (x1, y0), (x1, y1), (x0, y1)]
verts = np.array([[x, y, _plane(x, y)] for x, y in corners])
tris = np.array([[0, 1, 2], [0, 2, 3]])
z_grid = rasterize_triangle_mesh(verts, tris, x_coords, y_coords)
assert np.isfinite(z_grid).all(), f"빈 셀 {int((~np.isfinite(z_grid)).sum()):,}개"
xx, yy = np.meshgrid(np.asarray(x_coords, np.float64), np.asarray(y_coords, np.float64))
assert np.abs(z_grid - _plane(xx, yy)).max() < 1e-6
def test_small_triangles_on_float32_axis_cover_every_cell():
"""실데이터처럼 격자보다 작은 삼각형이 촘촘할 때도 누락이 없어야 한다."""
x_coords = np.linspace(183433.6, 183463.6, 31, dtype=np.float32)
y_coords = np.linspace(489168.8, 489198.8, 31, dtype=np.float32)
# 0.5m 간격 정점망을 사각형→삼각형 2개로 쪼개 촘촘한 TIN을 흉내낸다.
gx = np.arange(float(x_coords[0]) - 1.0, float(x_coords[-1]) + 1.5, 0.5)
gy = np.arange(float(y_coords[0]) - 1.0, float(y_coords[-1]) + 1.5, 0.5)
mx, my = np.meshgrid(gx, gy)
verts = np.column_stack([mx.ravel(), my.ravel(), _plane(mx.ravel(), my.ravel())])
n_col = len(gx)
tris = []
for r in range(len(gy) - 1):
for c in range(n_col - 1):
p00, p10 = r * n_col + c, r * n_col + c + 1
p01, p11 = (r + 1) * n_col + c, (r + 1) * n_col + c + 1
tris.append([p00, p10, p11])
tris.append([p00, p11, p01])
z_grid = rasterize_triangle_mesh(verts, np.array(tris), x_coords, y_coords)
assert np.isfinite(z_grid).all(), f"빈 셀 {int((~np.isfinite(z_grid)).sum()):,}개"
xx, yy = np.meshgrid(np.asarray(x_coords, np.float64), np.asarray(y_coords, np.float64))
assert np.abs(z_grid - _plane(xx, yy)).max() < 1e-6
def test_descending_axis_is_supported():
"""내림차순 축(예: 북쪽이 위인 y축)도 같은 결과를 낸다."""
corners = [(0.0, 0.0), (4.0, 0.0), (4.0, 4.0), (0.0, 4.0)]
verts = np.array([[x, y, _plane(x, y)] for x, y in corners])
tris = np.array([[0, 1, 2], [0, 2, 3]])
x_coords = np.linspace(0.0, 4.0, 5)
y_desc = np.linspace(4.0, 0.0, 5)
z_grid = rasterize_triangle_mesh(verts, tris, x_coords, y_desc)
xx, yy = np.meshgrid(x_coords, y_desc)
assert np.isfinite(z_grid).all()
assert np.allclose(z_grid, _plane(xx, yy), atol=1e-9)
def test_empty_mesh_returns_all_nan():
x_coords, y_coords = _axes(0.0, 3.0, 0.0, 3.0)
z_grid = rasterize_triangle_mesh(
np.zeros((0, 3)), np.zeros((0, 3), dtype=np.int64), x_coords, y_coords
)
assert z_grid.shape == (len(y_coords), len(x_coords))
assert np.isnan(z_grid).all()
+2 -60
View File
@@ -1,31 +1,18 @@
"""B04 지면 필터 회귀 테스트 — 고정 상수 때문에 지면점이 사라지던 결함 방지.
원 결함(2026-09-01): CSF 하강 예산이 0.3185m x 150회 = 47.8m로 고정돼 기복이
그보다 큰 산악지에서 천이 지면에 닿지 못했고(용화.las 지면점 0.06%),
grid_min_z의 3x3 minimum_filter가 급경사에서 기준면을 경사만큼 파고들었다.
"""
"""B04 지면 필터 회귀 테스트 — 업체 분류 필터 · 자동 고름(CSF · 격자 최저는 608 에서 걷음)."""
import numpy as np
import pytest
from B04_PreProcess.B04_PreProcess_Engine_Filter_Classification import (
filter_classification,
has_classified_ground,
)
from B04_PreProcess.B04_PreProcess_Engine_Filter_CSF import filter_csf
from B04_PreProcess.B04_PreProcess_Engine_Filter_Grid import filter_grid_min_z
from B04_PreProcess.B04_PreProcess_Engine_Ground import resolve_auto_source_filters
# CSF 고정 하강 예산 (옛 코드): 9.8 * 0.05 * 0.65 * 150
LEGACY_CSF_DROP_BUDGET_M = 47.775
def _sloped_ground(relief_m: float, slope: float = 0.3, spacing_m: float = 0.5) -> np.ndarray:
"""지정한 경사로 relief_m 만큼 오르내리는 지면 점군을 만든다.
경사와 기복을 따로 잡는다 — 결함은 경사가 아니라 **절대 기복**이 CSF 하강
예산을 넘을 때 터졌다. 점 간격은 격자·천 셀보다 촘촘하게 둔다(실제 LAS와 동일).
셀 안 표고 편차가 임계값을 넘으면 필터와 무관하게 점이 떨어지기 때문이다.
점 간격은 격자 칸보다 촘촘하게 둔다(실제 LAS와 동일).
"""
axis = np.arange(0.0, relief_m / slope, spacing_m, dtype=np.float64)
x, y = np.meshgrid(axis, axis)
@@ -49,51 +36,6 @@ def _structured(points: np.ndarray, classification: np.ndarray | None = None) ->
return data
@pytest.mark.parametrize("relief_m", [20.0, 120.0])
def test_csf_recovers_ground_regardless_of_relief(relief_m: float) -> None:
"""기복이 옛 고정 예산(47.8m)을 넘어도 지면이 살아남는다."""
ground = _sloped_ground(relief_m)
mask = filter_csf(_structured(ground))
assert mask.mean() > 0.5, f"기복 {relief_m}m 에서 지면 회수율 {mask.mean():.1%}"
# 옛 결함은 낮은 쪽 띠만 남기고 높은 쪽을 통째로 버렸다.
high_side = ground[:, 2] > ground[:, 2].min() + LEGACY_CSF_DROP_BUDGET_M
if high_side.any():
assert mask[high_side].mean() > 0.5, "옛 하강 예산 위쪽 지면이 잡히지 않는다"
def test_csf_rejects_canopy_points() -> None:
"""6단계 수목 필터가 실제로 동작한다 (피연산자 역전 시 항상 참이 되던 자리)."""
ground = _sloped_ground(60.0, slope=0.3)
canopy = ground.copy()
canopy[:, 2] += 12.0
points = np.vstack([ground, canopy])
mask = filter_csf(_structured(points))
is_canopy = np.zeros(len(points), dtype=bool)
is_canopy[len(ground) :] = True
assert not mask[is_canopy].any(), "수목 반사점이 지면으로 분류됐다"
assert mask[~is_canopy].mean() > 0.5
def test_grid_min_z_keeps_ground_on_steep_slope() -> None:
"""급경사에서 3x3 최소필터가 기준면을 파고들던 결함 방지."""
# 경사 100%. 셀 1m 안 표고 편차는 1.0m로 임계값(1.5m) 이내지만, 옛 코드는
# 3x3 최소필터가 2m 밖 셀 값을 끌어와 기준면을 2.0m 낮췄다 -> 임계값 초과.
ground = _sloped_ground(120.0, slope=1.0)
mask = filter_grid_min_z(_structured(ground), cell_size=1.0, height_threshold=1.5)
assert mask.mean() > 0.9, f"급경사 지면 회수율 {mask.mean():.1%}"
def test_grid_min_z_rejects_points_above_threshold() -> None:
ground = _sloped_ground(10.0)
above = ground.copy()
above[:, 2] += 5.0 # 임계값 1.5m 초과
points = np.vstack([ground, above])
mask = filter_grid_min_z(_structured(points))
assert not mask[len(ground) :].any()
def test_classification_filter_reads_class_two() -> None:
points = _sloped_ground(10.0)
classification = np.ones(len(points), dtype=np.uint8)
@@ -296,7 +296,7 @@ def _generate(project_root, options=None):
return run_section_generation(
project_root,
"B05_Profile/route/route_main.geojson",
"csf",
"ptd",
"dtm",
True,
options=options or SectionGenerationOptions(station_interval_m=20.0),
@@ -71,7 +71,7 @@ def test_request_wins_over_stored():
def test_schema_defaults_to_none_and_ships_the_flag():
"""스키마 기본은 None(미지정) — 저장값을 덮어쓰지 않는다."""
request = RouteSolveRequest(
filter_key="csf",
filter_key="ptd",
bp={"x": 0.0, "y": 0.0},
ep={"x": 100.0, "y": 0.0},
)
@@ -60,7 +60,7 @@ def test_샘플링_조건이_없으면_None(tmp_path: Path) -> None:
path.write_text(json.dumps({"filter_key": "", "method": "dtm"}), encoding="utf-8")
assert load_sampling_snapshot(tmp_path) is None
path.write_text(
json.dumps({"filter_key": "csf", "method": "tin", "smooth": True}), encoding="utf-8"
json.dumps({"filter_key": "ptd", "method": "implicit", "smooth": True}), encoding="utf-8"
)
snapshot = load_sampling_snapshot(tmp_path)
assert snapshot is not None and snapshot["method"] == "tin" and snapshot["smooth"] is True
assert snapshot is not None and snapshot["method"] == "implicit" and snapshot["smooth"] is True
+4 -1
View File
@@ -104,7 +104,10 @@ def test_every_mode_keeps_grade_limit_on_hill(hill: Terrain, mode: str) -> None:
_kept(result)
metrics = result["metrics"]
assert metrics["over_limit"]["total"] == 0
assert metrics["max_grade_pct"] <= GRADE * 100.0 + 1e-6
# 영선 = 상한 기울기를 걷는 선 — 곡선 넣기로 조금 넘을 수 있음(원지반 기울기 참고만 · 419 ·
# 측점 맞춤 없는 가이드 선 · 769-1)
if mode != "gradeline":
assert metrics["max_grade_pct"] <= GRADE * 100.0 + 1e-6
# 기점 · 종점은 옮기지 않고 복도(100 m) 밖으로 안 나감
assert result["planned"][0] == [50.0, 200.0]
assert result["planned"][-1] == [750.0, 200.0]
+6 -35
View File
@@ -72,14 +72,14 @@ def test_sheet_surface_npz_and_sampler_roundtrip(sheet_project: dict) -> None:
sheet_project["models"],
route_xy,
EPSG,
methods=["tin_sheet"],
methods=["laplace"],
)
assert len(models) == 1
model = models[0]
assert model["model_type"] == "dtm"
assert model["generation_params"]["source_filter"] == "sheet_tin_sheet"
assert model["generation_params"]["source_filter"] == "sheet_laplace"
npz_path = sheet_project["models"] / "dtm_sheet_tin_sheet.npz"
npz_path = sheet_project["models"] / "dtm_sheet_laplace.npz"
assert npz_path.is_file()
with np.load(npz_path) as data:
for key in ("x", "y", "z", "valid_mask"):
@@ -90,17 +90,17 @@ def test_sheet_surface_npz_and_sampler_roundtrip(sheet_project: dict) -> None:
assert data["z"].shape == (len(y), len(x))
# 프리뷰 glb
assert (sheet_project["models"] / "dtm_sheet_tin_sheet_preview.glb").stat().st_size > 0
assert (sheet_project["models"] / "dtm_sheet_laplace_preview.glb").stat().st_size > 0
# 종·횡단 공용 sampler로 왕복 — 등고선 안쪽은 유효 표고.
# 최고 등고선(120) 안쪽은 마루 캡(경사 연장)으로 최대 +간격-0.5m까지 올라갈 수 있다.
sampler = build_surface_sampler(sheet_project["models"], "sheet_tin_sheet", "dtm", smooth=False)
sampler = build_surface_sampler(sheet_project["models"], "sheet_laplace", "dtm", smooth=False)
inner = np.array([[X0, Y0 + 200.0], [X0 + 100.0, Y0 - 100.0]], dtype=np.float64)
z, valid = sampler.sample_xy(inner)
assert valid.all()
assert np.all((z >= 100.0 - 1e-6) & (z <= 124.5 + 1e-3))
# TIN 바깥(등고선 밖 + margin 밖)은 무효
# 바깥(등고선 밖 + margin 밖)은 무효
outer = np.array([[X0 + 10000.0, Y0]], dtype=np.float64)
_, outer_valid = sampler.sample_xy(outer)
assert not outer_valid.any()
@@ -185,32 +185,3 @@ def test_cone_reproduction_beats_over_relaxation() -> None:
relax_laplace(heavy, np.isfinite(burned), 200)
# 수렴시키면 harmonic 해로 끌려가 오차가 눈에 띄게 커진다
assert float(np.abs(heavy - truth)[band].mean()) > base_error
def test_tin_sheet_follows_source_contour_lines(sheet_project: dict) -> None:
"""TIN(도엽선)은 원본 등고선 자리에서 그 표고를 그대로 낸다 — 격자 계단을 타지 않는다.
링을 셀 격자와 어긋난 위치(+0.37m)에 두어, 격자에 구운 라인 셀이 아니라 벡터
정점을 쓰는지 가른다.
"""
route_xy = np.array([[X0 - 200.0, Y0 - 200.0], [X0 + 200.0, Y0 + 200.0]], dtype=np.float64)
models = build_sheet_surface_model(
sheet_project["root"],
sheet_project["processed"],
sheet_project["models"],
route_xy,
EPSG,
methods=["tin_sheet"],
)
assert len(models) == 1
sampler = build_surface_sampler(sheet_project["models"], "sheet_tin_sheet", "dtm", smooth=False)
# 각 링 변 위의 점 — 표고는 그 링의 값이어야 한다
probes, expected = [], []
for half, elev in ((300.0, 110.0), (150.0, 120.0)):
for offset in (-100.0, 0.0, 100.0):
probes.append([X0 + offset + 0.37, Y0 - half + 0.37])
expected.append(elev)
z, valid = sampler.sample_xy(np.asarray(probes))
assert valid.all()
assert np.abs(z - np.asarray(expected)).max() < 1.0