|
1 | | -""" pyplots.ai |
| 1 | +"""pyplots.ai |
2 | 2 | stereonet-equal-area: Structural Geology Stereonet (Equal-Area Projection) |
3 | 3 | Library: letsplot 4.9.0 | Python 3.14.3 |
4 | 4 | Quality: 84/100 | Created: 2026-03-15 |
|
16 | 16 | geom_density2d, |
17 | 17 | geom_path, |
18 | 18 | geom_point, |
| 19 | + geom_polygon, |
19 | 20 | geom_segment, |
20 | 21 | geom_text, |
21 | 22 | ggplot, |
22 | 23 | ggsize, |
23 | 24 | labs, |
24 | 25 | layer_tooltips, |
25 | 26 | scale_color_manual, |
| 27 | + scale_fill_manual, |
26 | 28 | scale_x_continuous, |
27 | 29 | scale_y_continuous, |
28 | 30 | theme, |
|
104 | 106 | sc = 1.0 / np.sqrt(1.0 - pz) |
105 | 107 | gx = px * sc |
106 | 108 | gy = py * sc |
107 | | - for j in range(len(alphas)): |
108 | | - gc_rows.append({"x": gx[j], "y": gy[j], "group": i, "feature_type": feature_types[i]}) |
| 109 | + # Clip to unit circle boundary |
| 110 | + r2 = gx**2 + gy**2 |
| 111 | + mask_inside = r2 <= 1.0 |
| 112 | + gx_clip = gx[mask_inside] |
| 113 | + gy_clip = gy[mask_inside] |
| 114 | + for j in range(len(gx_clip)): |
| 115 | + gc_rows.append({"x": gx_clip[j], "y": gy_clip[j], "group": i, "feature_type": feature_types[i]}) |
109 | 116 |
|
110 | 117 | gc_df = pd.DataFrame(gc_rows) |
111 | 118 |
|
|
129 | 136 | # Cardinal direction labels |
130 | 137 | label_df = pd.DataFrame({"x": [0, 1.09, 0, -1.09], "y": [1.09, 0, -1.09, 0], "label": ["N", "E", "S", "W"]}) |
131 | 138 |
|
132 | | -# Mean pole annotations per feature type (circular mean for strike) |
| 139 | +# Mean pole annotations and confidence ellipses per feature type |
133 | 140 | mean_annotations = [] |
134 | | -for ft in ["Bedding", "Joint", "Fault", "Foliation"]: |
| 141 | +ellipse_rows = [] |
| 142 | +for idx, ft in enumerate(["Bedding", "Joint", "Fault", "Foliation"]): |
135 | 143 | mask = np.array(feature_types) == ft |
136 | 144 | mx = pole_x[mask].mean() |
137 | 145 | my = pole_y[mask].mean() |
|
140 | 148 | ms = np.degrees(np.arctan2(np.sin(s_rad).mean(), np.cos(s_rad).mean())) % 360 |
141 | 149 | md = dips[mask].mean() |
142 | 150 | mean_annotations.append({"x": mx, "y": my, "label": f"{ft}\n{ms:.0f}/{md:.0f}"}) |
| 151 | + # Confidence ellipse (1-sigma) around each cluster |
| 152 | + px_ft = pole_x[mask] |
| 153 | + py_ft = pole_y[mask] |
| 154 | + cov = np.cov(px_ft, py_ft) |
| 155 | + eigvals, eigvecs = np.linalg.eigh(cov) |
| 156 | + angle = np.arctan2(eigvecs[1, 1], eigvecs[0, 1]) |
| 157 | + t = np.linspace(0, 2 * np.pi, 40) |
| 158 | + ex = np.sqrt(eigvals[1]) * np.cos(t) |
| 159 | + ey = np.sqrt(eigvals[0]) * np.sin(t) |
| 160 | + rx = mx + ex * np.cos(angle) - ey * np.sin(angle) |
| 161 | + ry = my + ex * np.sin(angle) + ey * np.cos(angle) |
| 162 | + for j in range(len(t)): |
| 163 | + ellipse_rows.append({"x": rx[j], "y": ry[j], "group": idx, "feature_type": ft}) |
| 164 | + |
143 | 165 | mean_df = pd.DataFrame(mean_annotations) |
| 166 | +ellipse_df = pd.DataFrame(ellipse_rows) |
144 | 167 |
|
145 | 168 | # Colorblind-safe palette (4 feature types) |
146 | 169 | color_values = ["#306998", "#CC79A7", "#E69F00", "#009E73"] |
|
169 | 192 | ) |
170 | 193 | # Great circles |
171 | 194 | + geom_path(aes(x="x", y="y", group="group", color="feature_type"), data=gc_df, size=0.4, alpha=0.3) |
| 195 | + # Confidence ellipses around each cluster |
| 196 | + + geom_polygon( |
| 197 | + aes(x="x", y="y", group="group", fill="feature_type"), data=ellipse_df, alpha=0.12, size=0, show_legend=False |
| 198 | + ) |
172 | 199 | # Poles to planes with tooltips |
173 | 200 | + geom_point(aes(x="x", y="y", color="feature_type"), data=poles_df, size=5, alpha=0.85, tooltips=pole_tooltips) |
174 | 201 | # Mean pole markers (larger, outlined) |
|
177 | 204 | + geom_text( |
178 | 205 | aes(x="x", y="y", label="label"), |
179 | 206 | data=mean_df, |
180 | | - size=11, |
181 | | - color="#333333", |
182 | | - nudge_y=-0.08, |
| 207 | + size=14, |
| 208 | + color="#222222", |
| 209 | + nudge_y=-0.12, |
183 | 210 | fontface="italic", |
184 | 211 | show_legend=False, |
185 | 212 | ) |
|
191 | 218 | + geom_text(aes(x="x", y="y", label="label"), data=label_df, size=18, color="#333333", fontface="bold") |
192 | 219 | # Color scale |
193 | 220 | + scale_color_manual(values=color_values, name="Feature Type") |
| 221 | + + scale_fill_manual(values=color_values) |
194 | 222 | + coord_fixed() |
195 | 223 | + scale_x_continuous(limits=[-1.25, 1.25]) |
196 | 224 | + scale_y_continuous(limits=[-1.25, 1.25]) |
|
0 commit comments