"""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()