|
1 | | -""" pyplots.ai |
| 1 | +"""pyplots.ai |
2 | 2 | stereonet-equal-area: Structural Geology Stereonet (Equal-Area Projection) |
3 | 3 | Library: altair 6.0.0 | Python 3.14.3 |
4 | 4 | Quality: 79/100 | Created: 2026-03-15 |
|
95 | 95 | } |
96 | 96 | ) |
97 | 97 |
|
98 | | -# Density heatmap using simple Gaussian KDE with numpy |
99 | | -n_grid = 60 |
| 98 | +# Density contours using Gaussian KDE - extract iso-lines via marching squares |
| 99 | +n_grid = 100 |
100 | 100 | gx = np.linspace(-r_prim, r_prim, n_grid) |
101 | 101 | gy = np.linspace(-r_prim, r_prim, n_grid) |
102 | 102 | gxx, gyy = np.meshgrid(gx, gy) |
103 | | -bw = 0.12 |
| 103 | +bw = 0.10 |
104 | 104 | density = np.zeros_like(gxx) |
105 | 105 | for px, py in zip(pole_x, pole_y, strict=True): |
106 | 106 | density += np.exp(-((gxx - px) ** 2 + (gyy - py) ** 2) / (2 * bw**2)) |
107 | 107 | density /= len(pole_x) * 2 * np.pi * bw**2 |
108 | | - |
109 | 108 | circ_mask = gxx**2 + gyy**2 <= r_prim**2 |
110 | | -density[~circ_mask] = np.nan |
111 | | -cell_w = gx[1] - gx[0] |
112 | | -cell_h = gy[1] - gy[0] |
| 109 | +density[~circ_mask] = 0.0 |
| 110 | + |
| 111 | +# Extract contour lines using simple marching squares |
| 112 | +contour_levels = np.linspace(np.nanmax(density) * 0.15, np.nanmax(density) * 0.85, 5) |
| 113 | +contour_rows = [] |
| 114 | +for lev_idx, level in enumerate(contour_levels): |
| 115 | + for ri in range(n_grid - 1): |
| 116 | + for ci in range(n_grid - 1): |
| 117 | + corners = [density[ri, ci], density[ri, ci + 1], density[ri + 1, ci + 1], density[ri + 1, ci]] |
| 118 | + above = [c >= level for c in corners] |
| 119 | + n_above = sum(above) |
| 120 | + if n_above == 0 or n_above == 4: |
| 121 | + continue |
| 122 | + edges = [] |
| 123 | + edge_pairs = [ |
| 124 | + (0, 1, ri, ci, ri, ci + 1), |
| 125 | + (1, 2, ri, ci + 1, ri + 1, ci + 1), |
| 126 | + (2, 3, ri + 1, ci + 1, ri + 1, ci), |
| 127 | + (3, 0, ri + 1, ci, ri, ci), |
| 128 | + ] |
| 129 | + for c1, c2, r1, col1, r2, col2 in edge_pairs: |
| 130 | + if above[c1] != above[c2]: |
| 131 | + t = (level - corners[c1]) / (corners[c2] - corners[c1]) if corners[c2] != corners[c1] else 0.5 |
| 132 | + ex = gxx[r1, col1] + t * (gxx[r2, col2] - gxx[r1, col1]) |
| 133 | + ey = gyy[r1, col1] + t * (gyy[r2, col2] - gyy[r1, col1]) |
| 134 | + edges.append((ex, ey)) |
| 135 | + if len(edges) >= 2: |
| 136 | + seg_id = f"c{lev_idx}_{ri}_{ci}" |
| 137 | + contour_rows.append( |
| 138 | + {"x": edges[0][0], "y": edges[0][1], "seg_id": seg_id, "level": float(lev_idx), "order": 0} |
| 139 | + ) |
| 140 | + contour_rows.append( |
| 141 | + {"x": edges[1][0], "y": edges[1][1], "seg_id": seg_id, "level": float(lev_idx), "order": 1} |
| 142 | + ) |
| 143 | + |
| 144 | +df_contours = ( |
| 145 | + pd.DataFrame(contour_rows) if contour_rows else pd.DataFrame(columns=["x", "y", "seg_id", "level", "order"]) |
| 146 | +) |
113 | 147 |
|
114 | | -density_rows = [] |
115 | | -for ri in range(n_grid): |
116 | | - for ci in range(n_grid): |
117 | | - if circ_mask[ri, ci] and density[ri, ci] > np.nanmax(density) * 0.15: |
118 | | - density_rows.append({"x": gxx[ri, ci], "y": gyy[ri, ci], "density": density[ri, ci]}) |
119 | | -df_density = pd.DataFrame(density_rows) if density_rows else pd.DataFrame(columns=["x", "y", "density"]) |
| 148 | +# Mean orientation for each feature type (for annotations) |
| 149 | +mean_rows = [] |
| 150 | +for ft in ["Bedding", "Fault", "Joint"]: |
| 151 | + mask = np.array(feature_types) == ft |
| 152 | + mx = np.mean(pole_x[mask]) |
| 153 | + my = np.mean(pole_y[mask]) |
| 154 | + ms = np.mean(strikes[mask]) |
| 155 | + md = np.mean(dips[mask]) |
| 156 | + mean_rows.append({"x": mx, "y": my, "feature_type": ft, "label": f"μ {ms:.0f}°/{md:.0f}°"}) |
| 157 | +df_means = pd.DataFrame(mean_rows) |
120 | 158 |
|
121 | 159 | # Plot |
122 | 160 | color_map = {"Bedding": "#306998", "Fault": "#C1432E", "Joint": "#2A9D8F"} |
123 | 161 | color_scale = alt.Scale(domain=list(color_map.keys()), range=list(color_map.values())) |
124 | 162 | x_enc = alt.X("x:Q", axis=None, scale=alt.Scale(domain=[-1.25, 1.25])) |
125 | 163 | y_enc = alt.Y("y:Q", axis=None, scale=alt.Scale(domain=[-1.25, 1.25])) |
126 | 164 |
|
127 | | -density_layer = ( |
| 165 | +# Interactive selection for highlighting by feature type |
| 166 | +selection = alt.selection_point(fields=["feature_type"], bind="legend") |
| 167 | + |
| 168 | +# Density contour lines (smooth, not blocky) |
| 169 | +contour_layer = ( |
128 | 170 | ( |
129 | | - alt.Chart(df_density) |
130 | | - .mark_rect(opacity=0.25) |
| 171 | + alt.Chart(df_contours) |
| 172 | + .mark_line(strokeWidth=1.0) |
131 | 173 | .encode( |
132 | | - x=alt.X("x:Q", axis=None, scale=alt.Scale(domain=[-1.25, 1.25]), bin=alt.Bin(step=cell_w)), |
133 | | - y=alt.Y("y:Q", axis=None, scale=alt.Scale(domain=[-1.25, 1.25]), bin=alt.Bin(step=cell_h)), |
134 | | - color=alt.Color("max(density):Q", scale=alt.Scale(scheme="greys"), legend=None), |
| 174 | + x=x_enc, |
| 175 | + y=y_enc, |
| 176 | + detail="seg_id:N", |
| 177 | + order="order:O", |
| 178 | + color=alt.Color("level:Q", scale=alt.Scale(scheme="greys", domain=[0, len(contour_levels)]), legend=None), |
| 179 | + opacity=alt.value(0.35), |
135 | 180 | ) |
136 | 181 | ) |
137 | | - if len(df_density) > 0 |
| 182 | + if len(df_contours) > 0 |
138 | 183 | else alt.Chart(pd.DataFrame({"x": [0], "y": [0]})).mark_point(size=0).encode(x="x:Q", y="y:Q") |
139 | 184 | ) |
140 | 185 |
|
|
152 | 197 |
|
153 | 198 | great_circles = ( |
154 | 199 | alt.Chart(df_gc) |
155 | | - .mark_line(strokeWidth=1.0, opacity=0.35) |
| 200 | + .mark_line(strokeWidth=1.0) |
156 | 201 | .encode( |
157 | 202 | x=x_enc, |
158 | 203 | y=y_enc, |
159 | 204 | detail="gc_id:N", |
160 | 205 | order="order:O", |
161 | 206 | color=alt.Color("feature_type:N", scale=color_scale, legend=None), |
| 207 | + opacity=alt.condition(selection, alt.value(0.35), alt.value(0.05)), |
162 | 208 | ) |
| 209 | + .add_params(selection) |
163 | 210 | ) |
164 | 211 |
|
165 | 212 | prim_circle = alt.Chart(df_circle).mark_line(strokeWidth=2.5, color="#333333").encode(x=x_enc, y=y_enc, order="order:O") |
|
178 | 225 |
|
179 | 226 | poles_layer = ( |
180 | 227 | alt.Chart(df_poles) |
181 | | - .mark_point(filled=True, size=180, stroke="white", strokeWidth=1.2, opacity=0.9) |
| 228 | + .mark_point(filled=True, size=180, stroke="white", strokeWidth=1.2) |
182 | 229 | .encode( |
183 | 230 | x=x_enc, |
184 | 231 | y=y_enc, |
|
190 | 237 | titleFontSize=20, |
191 | 238 | labelFontSize=18, |
192 | 239 | symbolSize=300, |
193 | | - orient="right", |
| 240 | + orient="top-right", |
194 | 241 | titleColor="#333333", |
195 | 242 | labelColor="#333333", |
196 | 243 | ), |
197 | 244 | ), |
| 245 | + opacity=alt.condition(selection, alt.value(0.9), alt.value(0.15)), |
198 | 246 | tooltip=[ |
199 | 247 | alt.Tooltip("feature_type:N", title="Type"), |
200 | | - alt.Tooltip("strike:Q", title="Strike"), |
201 | | - alt.Tooltip("dip:Q", title="Dip"), |
| 248 | + alt.Tooltip("strike:Q", title="Strike (°)"), |
| 249 | + alt.Tooltip("dip:Q", title="Dip (°)"), |
202 | 250 | ], |
203 | 251 | ) |
| 252 | + .add_params(selection) |
| 253 | +) |
| 254 | + |
| 255 | +# Mean orientation markers with labels |
| 256 | +mean_markers = ( |
| 257 | + alt.Chart(df_means) |
| 258 | + .mark_point(shape="cross", size=400, strokeWidth=3) |
| 259 | + .encode( |
| 260 | + x=x_enc, |
| 261 | + y=y_enc, |
| 262 | + color=alt.Color("feature_type:N", scale=color_scale, legend=None), |
| 263 | + opacity=alt.condition(selection, alt.value(1.0), alt.value(0.15)), |
| 264 | + ) |
| 265 | + .add_params(selection) |
| 266 | +) |
| 267 | + |
| 268 | +mean_labels = ( |
| 269 | + alt.Chart(df_means) |
| 270 | + .mark_text(fontSize=16, fontWeight="bold", dy=-18, color="#333333") |
| 271 | + .encode(x=x_enc, y=y_enc, text="label:N", opacity=alt.condition(selection, alt.value(0.9), alt.value(0.1))) |
| 272 | + .add_params(selection) |
204 | 273 | ) |
205 | 274 |
|
206 | 275 | chart = ( |
207 | | - alt.layer(density_layer, grid_circles, cross_lines, great_circles, prim_circle, tick_marks, dir_labels, poles_layer) |
| 276 | + alt.layer( |
| 277 | + contour_layer, |
| 278 | + grid_circles, |
| 279 | + cross_lines, |
| 280 | + great_circles, |
| 281 | + prim_circle, |
| 282 | + tick_marks, |
| 283 | + dir_labels, |
| 284 | + poles_layer, |
| 285 | + mean_markers, |
| 286 | + mean_labels, |
| 287 | + ) |
208 | 288 | .properties( |
209 | 289 | width=1200, |
210 | 290 | height=1200, |
|
0 commit comments