Coverage for graphica/core/grid_data.py: 100.0%

1"""(x, y, z) の点から 2D マップの格子を作る。GUI には依存しない。 

2 

3(x, y) が完全な格子なら pivot でそのまま並べ(is_regular=True、値を変えない)、そうでなければ griddata で補間する。 

4""" 

5import numpy as np 

6import pandas as pd 

7from scipy.interpolate import griddata, RegularGridInterpolator 

8from typing import Any 

9 

10# ドラッグでぴったり水平・垂直には引けないので、各軸の範囲に対するこの割合までは水平・垂直とみなす 

11SLICE_AXIS_ALIGNMENT_TOLERANCE = 0.01 

12 

13GRID_INTERP_METHODS = ('linear', 'cubic', 'nearest') 

14 

15 

16class GridDataError(ValueError): 

17 """2D の格子を作れない。""" 

18 

19 

20def is_regular_grid(x: Any, y: Any) -> bool: 

21 """x と y の全部の組み合わせが、ちょうど1回ずつあるか。""" 

22 unique_x = np.unique(x) 

23 unique_y = np.unique(y) 

24 if len(unique_x) * len(unique_y) != len(x): 

25 return False 

26 pairs = set(zip(x.tolist(), y.tolist())) 

27 return len(pairs) == len(x) 

28 

29 

30def compute_z_grid(x: Any, y: Any, z: Any, interp_method: str = 'linear', 

31 resolution: tuple[int, int] | list[int] | None = None) -> dict[str, Any]: 

32 """{'x_grid', 'y_grid'(ソート済みの1次元), 'z_grid'(shape=(len(y), len(x))、データの無い所は nan), 'is_regular'}。 

33 

34 resolution=(nx, ny) は補間するときだけ使う(None なら点の数から決める)。interp_method も補間のときだけ。 

35 """ 

36 if interp_method not in GRID_INTERP_METHODS: 

37 raise GridDataError(f"不明な補間方法です: {interp_method}") 

38 

39 x = np.asarray(x, dtype=float) 

40 y = np.asarray(y, dtype=float) 

41 z = np.asarray(z, dtype=float) 

42 

43 valid = ~(np.isnan(x) | np.isnan(y) | np.isnan(z)) 

44 x, y, z = x[valid], y[valid], z[valid] 

45 

46 if len(x) == 0: 

47 raise GridDataError("有効なデータ点がありません(すべて欠損値です)。") 

48 

49 if is_regular_grid(x, y): 

50 x_grid = np.unique(x) 

51 y_grid = np.unique(y) 

52 pivot = pd.DataFrame({'x': x, 'y': y, 'z': z}).pivot(index='y', columns='x', values='z') 

53 pivot = pivot.reindex(index=y_grid, columns=x_grid) 

54 return { 

55 'x_grid': x_grid, 

56 'y_grid': y_grid, 

57 'z_grid': pivot.values, 

58 'is_regular': True, 

59 } 

60 

61 if resolution is None: 

62 n = max(int(np.sqrt(len(x)) * 2), 10) 

63 nx, ny = n, n 

64 else: 

65 nx, ny = resolution 

66 if nx < 2 or ny < 2: 

67 raise GridDataError("補間グリッドの解像度は2以上である必要があります。") 

68 

69 if len(x) < 3: 

70 raise GridDataError( 

71 "散在データの補間には最低3点が必要です(規則格子ではない点が" 

72 f"{len(x)}点しかありません)。" 

73 ) 

74 

75 x_grid = np.linspace(x.min(), x.max(), nx) 

76 y_grid = np.linspace(y.min(), y.max(), ny) 

77 xx, yy = np.meshgrid(x_grid, y_grid) 

78 z_grid = griddata((x, y), z, (xx, yy), method=interp_method) 

79 

80 return { 

81 'x_grid': x_grid, 

82 'y_grid': y_grid, 

83 'z_grid': z_grid, 

84 'is_regular': False, 

85 } 

86 

87 

88def extract_slice(x_grid: Any, y_grid: Any, z_grid: Any, start: tuple[float, float], end: tuple[float, float], 

89 n_points: int = 200) -> dict[str, Any]: 

90 """格子上の線分に沿った断面を返す({'axis_values', 'axis_kind', 'z_values'})。 

91 

92 axis_kind はほぼ水平なら 'x'、ほぼ垂直なら 'y'、斜めなら 'distance'(始点からの距離)。格子の外は nan。 

93 """ 

94 x0, y0 = float(start[0]), float(start[1]) 

95 x1, y1 = float(end[0]), float(end[1]) 

96 if x0 == x1 and y0 == y1: 

97 raise GridDataError("始点と終点が同じ位置です(長さ0の線分は抽出できません)。") 

98 

99 x_grid = np.asarray(x_grid, dtype=float) 

100 y_grid = np.asarray(y_grid, dtype=float) 

101 z_grid = np.asarray(z_grid, dtype=float) 

102 

103 x_range = x_grid.max() - x_grid.min() if len(x_grid) > 1 else 1.0 

104 y_range = y_grid.max() - y_grid.min() if len(y_grid) > 1 else 1.0 

105 

106 is_horizontal = abs(y1 - y0) <= SLICE_AXIS_ALIGNMENT_TOLERANCE * (y_range or 1.0) 

107 is_vertical = abs(x1 - x0) <= SLICE_AXIS_ALIGNMENT_TOLERANCE * (x_range or 1.0) 

108 

109 if is_horizontal and not is_vertical: 

110 axis_kind = 'x' 

111 axis_values = np.linspace(x0, x1, n_points) 

112 sample_x, sample_y = axis_values, np.full(n_points, y0) 

113 elif is_vertical and not is_horizontal: 

114 axis_kind = 'y' 

115 axis_values = np.linspace(y0, y1, n_points) 

116 sample_x, sample_y = np.full(n_points, x0), axis_values 

117 else: 

118 axis_kind = 'distance' 

119 sample_x = np.linspace(x0, x1, n_points) 

120 sample_y = np.linspace(y0, y1, n_points) 

121 axis_values = np.sqrt((sample_x - x0) ** 2 + (sample_y - y0) ** 2) 

122 

123 interpolator = RegularGridInterpolator( 

124 (y_grid, x_grid), z_grid, method='linear', bounds_error=False, fill_value=np.nan, 

125 ) 

126 z_values = interpolator(np.column_stack([sample_y, sample_x])) 

127 

128 return {'axis_values': axis_values, 'axis_kind': axis_kind, 'z_values': z_values}