|
| 1 | +""" pyplots.ai |
| 2 | +stereonet-equal-area: Structural Geology Stereonet (Equal-Area Projection) |
| 3 | +Library: seaborn 0.13.2 | Python 3.14.3 |
| 4 | +Quality: 89/100 | Created: 2026-03-15 |
| 5 | +""" |
| 6 | + |
| 7 | +import matplotlib.pyplot as plt |
| 8 | +import numpy as np |
| 9 | +import pandas as pd |
| 10 | +import seaborn as sns |
| 11 | + |
| 12 | + |
| 13 | +# Data - field measurements from a geological mapping campaign |
| 14 | +np.random.seed(42) |
| 15 | + |
| 16 | +bedding_strike = np.random.normal(45, 8, 40) % 360 |
| 17 | +bedding_dip = np.clip(np.random.normal(35, 5, 40), 1, 89) |
| 18 | + |
| 19 | +joint1_strike = np.random.normal(315, 10, 30) % 360 |
| 20 | +joint1_dip = np.clip(np.random.normal(75, 8, 30), 1, 89) |
| 21 | + |
| 22 | +joint2_strike = np.random.normal(90, 12, 25) % 360 |
| 23 | +joint2_dip = np.clip(np.random.normal(60, 10, 25), 1, 89) |
| 24 | + |
| 25 | +fault_strike = np.random.normal(180, 15, 15) % 360 |
| 26 | +fault_dip = np.clip(np.random.normal(70, 10, 15), 1, 89) |
| 27 | + |
| 28 | +strikes = np.concatenate([bedding_strike, joint1_strike, joint2_strike, fault_strike]) |
| 29 | +dips = np.concatenate([bedding_dip, joint1_dip, joint2_dip, fault_dip]) |
| 30 | +feature_types = ["Bedding"] * 40 + ["Joint Set 1"] * 30 + ["Joint Set 2"] * 25 + ["Fault"] * 15 |
| 31 | + |
| 32 | +df = pd.DataFrame({"strike": strikes, "dip": dips, "feature_type": feature_types}) |
| 33 | + |
| 34 | +# Equal-area projection: convert pole trend/plunge to x, y |
| 35 | +pole_trend = (df["strike"].values + 270) % 360 |
| 36 | +pole_plunge = 90 - df["dip"].values |
| 37 | +pole_trend_rad = np.radians(pole_trend) |
| 38 | +pole_plunge_rad = np.radians(pole_plunge) |
| 39 | +pole_r = np.sqrt(2) * np.sin((np.pi / 2 - pole_plunge_rad) / 2) |
| 40 | +df["pole_x"] = pole_r * np.sin(pole_trend_rad) |
| 41 | +df["pole_y"] = pole_r * np.cos(pole_trend_rad) |
| 42 | + |
| 43 | +# Scale marker size by dip angle for visual hierarchy (steeper dips = larger markers) |
| 44 | +df["marker_size"] = 80 + (df["dip"] / 90) * 120 |
| 45 | + |
| 46 | +# Colorblind-safe palette: blue, orange, purple, teal (avoids green-red pairing) |
| 47 | +sns.set_theme(style="white", font="DejaVu Sans", font_scale=1.2) |
| 48 | +sns.set_context("talk", rc={"font.family": "DejaVu Sans"}) |
| 49 | +feature_colors = ["#306998", "#E88D39", "#8B5FC7", "#17A2B8"] |
| 50 | +palette = sns.color_palette(feature_colors) |
| 51 | +feature_names = ["Bedding", "Joint Set 1", "Joint Set 2", "Fault"] |
| 52 | +palette_dict = dict(zip(feature_names, palette, strict=True)) |
| 53 | + |
| 54 | +fig, ax = plt.subplots(figsize=(12, 12), facecolor="#FAFAFA") |
| 55 | +ax.set_facecolor("#FAFAFA") |
| 56 | +ax.set_aspect("equal") |
| 57 | + |
| 58 | +# Subtle outer glow ring for depth |
| 59 | +for ring_r, ring_alpha in [(1.02, 0.04), (1.04, 0.02)]: |
| 60 | + glow = plt.Circle((0, 0), ring_r, fill=False, color="#888888", linewidth=0.6, alpha=ring_alpha) |
| 61 | + ax.add_patch(glow) |
| 62 | + |
| 63 | +# Primitive circle with refined styling |
| 64 | +theta_circle = np.linspace(0, 2 * np.pi, 300) |
| 65 | +ax.fill(np.sin(theta_circle), np.cos(theta_circle), color="white", zorder=0) |
| 66 | +ax.plot(np.sin(theta_circle), np.cos(theta_circle), color="#2C2C2C", linewidth=1.8, zorder=4) |
| 67 | + |
| 68 | +# Grid: small circles (constant plunge) |
| 69 | +for plunge_deg in range(10, 90, 10): |
| 70 | + r = np.sqrt(2) * np.sin((np.pi / 2 - np.radians(plunge_deg)) / 2) |
| 71 | + grid_circle = plt.Circle((0, 0), r, fill=False, color="#CCCCCC", linewidth=0.3, linestyle=(0, (5, 5)), zorder=1) |
| 72 | + ax.add_patch(grid_circle) |
| 73 | + |
| 74 | +# Grid: great circles for N-S and E-W reference planes |
| 75 | +for ref_strike in [0, 90]: |
| 76 | + ref_strike_rad = np.radians(ref_strike) |
| 77 | + for ref_dip in range(10, 90, 10): |
| 78 | + ref_dip_rad = np.radians(ref_dip) |
| 79 | + alpha = np.linspace(0.01, np.pi - 0.01, 150) |
| 80 | + vx = np.cos(alpha) * np.sin(ref_strike_rad) + np.sin(alpha) * np.cos(ref_dip_rad) * np.cos(ref_strike_rad) |
| 81 | + vy = np.cos(alpha) * np.cos(ref_strike_rad) - np.sin(alpha) * np.cos(ref_dip_rad) * np.sin(ref_strike_rad) |
| 82 | + vz = -np.sin(alpha) * np.sin(ref_dip_rad) |
| 83 | + gc_plunge = np.arcsin(np.clip(-vz, -1, 1)) |
| 84 | + gc_trend = np.arctan2(vx, vy) |
| 85 | + gc_r = np.sqrt(2) * np.sin((np.pi / 2 - gc_plunge) / 2) |
| 86 | + gc_x = gc_r * np.sin(gc_trend) |
| 87 | + gc_y = gc_r * np.cos(gc_trend) |
| 88 | + ax.plot(gc_x, gc_y, color="#CCCCCC", linewidth=0.3, linestyle=(0, (5, 5)), zorder=1) |
| 89 | + |
| 90 | +# Tick marks every 10 degrees with graduated lengths |
| 91 | +for azimuth in range(0, 360, 10): |
| 92 | + az_rad = np.radians(azimuth) |
| 93 | + if azimuth % 90 == 0: |
| 94 | + tick_inner, tick_lw = 0.93, 1.2 |
| 95 | + elif azimuth % 30 == 0: |
| 96 | + tick_inner, tick_lw = 0.95, 0.9 |
| 97 | + else: |
| 98 | + tick_inner, tick_lw = 0.97, 0.6 |
| 99 | + ax.plot( |
| 100 | + [tick_inner * np.sin(az_rad), np.sin(az_rad)], |
| 101 | + [tick_inner * np.cos(az_rad), np.cos(az_rad)], |
| 102 | + color="#2C2C2C", |
| 103 | + linewidth=tick_lw, |
| 104 | + zorder=4, |
| 105 | + ) |
| 106 | + |
| 107 | +# Cardinal direction labels with refined typography |
| 108 | +for az, label in [(0, "N"), (90, "E"), (180, "S"), (270, "W")]: |
| 109 | + az_rad = np.radians(az) |
| 110 | + fontsize = 26 if label == "N" else 22 |
| 111 | + ax.text( |
| 112 | + 1.09 * np.sin(az_rad), |
| 113 | + 1.09 * np.cos(az_rad), |
| 114 | + label, |
| 115 | + ha="center", |
| 116 | + va="center", |
| 117 | + fontsize=fontsize, |
| 118 | + fontweight="bold", |
| 119 | + color="#1A1A1A", |
| 120 | + zorder=6, |
| 121 | + ) |
| 122 | + |
| 123 | +# Degree labels every 30 degrees (skip cardinals) |
| 124 | +for az in range(30, 360, 30): |
| 125 | + if az in [90, 180, 270]: |
| 126 | + continue |
| 127 | + az_rad = np.radians(az) |
| 128 | + ax.text( |
| 129 | + 1.08 * np.sin(az_rad), |
| 130 | + 1.08 * np.cos(az_rad), |
| 131 | + f"{az}°", |
| 132 | + ha="center", |
| 133 | + va="center", |
| 134 | + fontsize=16, |
| 135 | + color="#777777", |
| 136 | + zorder=6, |
| 137 | + ) |
| 138 | + |
| 139 | +# Density contours per feature type using seaborn kdeplot with hue |
| 140 | +n_collections_before = len(ax.collections) |
| 141 | +n_lines_before = len(ax.lines) |
| 142 | + |
| 143 | +sns.kdeplot( |
| 144 | + data=df, |
| 145 | + x="pole_x", |
| 146 | + y="pole_y", |
| 147 | + hue="feature_type", |
| 148 | + hue_order=feature_names, |
| 149 | + palette=[sns.desaturate(c, 0.3) for c in feature_colors], |
| 150 | + fill=True, |
| 151 | + levels=5, |
| 152 | + thresh=0.2, |
| 153 | + alpha=0.15, |
| 154 | + bw_adjust=0.7, |
| 155 | + ax=ax, |
| 156 | + zorder=2, |
| 157 | + legend=False, |
| 158 | +) |
| 159 | +sns.kdeplot( |
| 160 | + data=df, |
| 161 | + x="pole_x", |
| 162 | + y="pole_y", |
| 163 | + hue="feature_type", |
| 164 | + hue_order=feature_names, |
| 165 | + palette=[sns.desaturate(c, 0.5) for c in feature_colors], |
| 166 | + levels=5, |
| 167 | + thresh=0.2, |
| 168 | + linewidths=0.6, |
| 169 | + bw_adjust=0.7, |
| 170 | + ax=ax, |
| 171 | + zorder=2, |
| 172 | + legend=False, |
| 173 | +) |
| 174 | + |
| 175 | +# Clip density contours to the primitive circle |
| 176 | +clip_circle = plt.Circle((0, 0), 1.0, transform=ax.transData, fill=False) |
| 177 | +ax.add_patch(clip_circle) |
| 178 | +clip_circle.set_visible(False) |
| 179 | +for c in ax.collections[n_collections_before:]: |
| 180 | + c.set_clip_path(clip_circle) |
| 181 | +for line in ax.lines[n_lines_before:]: |
| 182 | + line.set_clip_path(clip_circle) |
| 183 | + |
| 184 | +# Great circles for representative planes (3 per feature type) |
| 185 | +for feature_type in df["feature_type"].unique(): |
| 186 | + subset = df[df["feature_type"] == feature_type] |
| 187 | + indices = np.linspace(0, len(subset) - 1, min(3, len(subset)), dtype=int) |
| 188 | + for idx in indices: |
| 189 | + row = subset.iloc[idx] |
| 190 | + s_rad = np.radians(row["strike"]) |
| 191 | + d_rad = np.radians(row["dip"]) |
| 192 | + alpha = np.linspace(0.01, np.pi - 0.01, 200) |
| 193 | + vx = np.cos(alpha) * np.sin(s_rad) + np.sin(alpha) * np.cos(d_rad) * np.cos(s_rad) |
| 194 | + vy = np.cos(alpha) * np.cos(s_rad) - np.sin(alpha) * np.cos(d_rad) * np.sin(s_rad) |
| 195 | + vz = -np.sin(alpha) * np.sin(d_rad) |
| 196 | + gc_plunge = np.arcsin(np.clip(-vz, -1, 1)) |
| 197 | + gc_trend = np.arctan2(vx, vy) |
| 198 | + gc_r = np.sqrt(2) * np.sin((np.pi / 2 - gc_plunge) / 2) |
| 199 | + gc_x = gc_r * np.sin(gc_trend) |
| 200 | + gc_y = gc_r * np.cos(gc_trend) |
| 201 | + ax.plot(gc_x, gc_y, color=palette_dict[feature_type], linewidth=2.0, alpha=0.45, zorder=3) |
| 202 | + |
| 203 | +# Poles using seaborn scatterplot with hue and size encoding |
| 204 | +sns.scatterplot( |
| 205 | + data=df, |
| 206 | + x="pole_x", |
| 207 | + y="pole_y", |
| 208 | + hue="feature_type", |
| 209 | + hue_order=feature_names, |
| 210 | + style="feature_type", |
| 211 | + style_order=feature_names, |
| 212 | + markers={"Bedding": "o", "Joint Set 1": "s", "Joint Set 2": "D", "Fault": "^"}, |
| 213 | + size="marker_size", |
| 214 | + sizes=(80, 200), |
| 215 | + palette=palette_dict, |
| 216 | + edgecolor="white", |
| 217 | + linewidth=0.9, |
| 218 | + alpha=0.88, |
| 219 | + ax=ax, |
| 220 | + zorder=5, |
| 221 | + legend=False, |
| 222 | +) |
| 223 | + |
| 224 | +# Style |
| 225 | +ax.set_xlim(-1.22, 1.22) |
| 226 | +ax.set_ylim(-1.22, 1.22) |
| 227 | +ax.set_xlabel("") |
| 228 | +ax.set_ylabel("") |
| 229 | +ax.set_xticks([]) |
| 230 | +ax.set_yticks([]) |
| 231 | +sns.despine(ax=ax, left=True, bottom=True) |
| 232 | + |
| 233 | +# Build a custom legend with marker shapes matching the plot |
| 234 | +marker_map = {"Bedding": "o", "Joint Set 1": "s", "Joint Set 2": "D", "Fault": "^"} |
| 235 | +handles = [] |
| 236 | +for ft in feature_names: |
| 237 | + handles.append( |
| 238 | + plt.Line2D( |
| 239 | + [0], |
| 240 | + [0], |
| 241 | + marker=marker_map[ft], |
| 242 | + color="none", |
| 243 | + markerfacecolor=palette_dict[ft], |
| 244 | + markersize=11, |
| 245 | + markeredgecolor="white", |
| 246 | + markeredgewidth=0.6, |
| 247 | + label=ft, |
| 248 | + ) |
| 249 | + ) |
| 250 | + |
| 251 | +legend = ax.legend( |
| 252 | + handles=handles, |
| 253 | + title="Feature Type", |
| 254 | + loc="upper left", |
| 255 | + fontsize=15, |
| 256 | + title_fontsize=17, |
| 257 | + framealpha=0.92, |
| 258 | + edgecolor="#BBBBBB", |
| 259 | + fancybox=True, |
| 260 | + shadow=False, |
| 261 | + bbox_to_anchor=(-0.02, 1.02), |
| 262 | +) |
| 263 | +legend.get_title().set_fontweight("semibold") |
| 264 | + |
| 265 | +# Measurement count annotation |
| 266 | +n_total = len(df) |
| 267 | +ax.text( |
| 268 | + 0.98, |
| 269 | + -0.02, |
| 270 | + f"n = {n_total} measurements", |
| 271 | + transform=ax.transAxes, |
| 272 | + fontsize=13, |
| 273 | + color="#888888", |
| 274 | + ha="right", |
| 275 | + va="top", |
| 276 | + style="italic", |
| 277 | +) |
| 278 | + |
| 279 | +ax.set_title("stereonet-equal-area · seaborn · pyplots.ai", fontsize=24, fontweight="medium", color="#2C2C2C", pad=30) |
| 280 | + |
| 281 | +# Save |
| 282 | +plt.tight_layout() |
| 283 | +plt.savefig("plot.png", dpi=300, bbox_inches="tight", facecolor="#FAFAFA") |
0 commit comments