PySDM_examples.Luettmer_et_al_2026.plot
shared plot functions for homogeneous freezing notebooks
1"""shared plot functions for homogeneous freezing notebooks""" 2 3from matplotlib import pyplot, ticker 4import numpy as np 5import seaborn as sns 6from cycler import cycler 7from PySDM import Formulae 8 9formulae = Formulae(particle_shape_and_density="MixedPhaseSpheres") 10 11ax_title_size = 18 12ax_lab_fsize = 15 13ax_lab_fsize2 = 19 14tick_fsize = 15 15T_frz_bins = np.linspace(-40, -34, num=60, endpoint=True) 16T_frz_bins_kelvin = np.linspace(230, 240, num=100, endpoint=True) 17 18 19def cumulative_histogram(data, bins, reverse=False, density=True): 20 hist, bin_edges = np.histogram(data, bins=bins, density=False) 21 22 if reverse: 23 cum_hist = np.cumsum(hist[::-1])[::-1] 24 cum_hist_0 = cum_hist[0] 25 else: 26 cum_hist = np.cumsum(hist) 27 cum_hist_0 = cum_hist[-1] 28 29 if density: 30 cum_hist = cum_hist / cum_hist_0 31 32 bin_centers = 0.5 * (bin_edges[:-1] + bin_edges[1:]) 33 return cum_hist, bin_centers 34 35 36def plot_thermodynamics_and_bulk( 37 simulation, 38 title_add=None, 39 show_conc=False, 40 show_jhom=True, 41 show_tf=True, 42 t_lim=None, 43): 44 if title_add is None: 45 title_add = ["", "", "", ""] 46 plot_daw = False 47 output = simulation["ensemble_member_outputs"][0] 48 time = output["t"] 49 T = np.asarray(output["T"]) 50 RH = np.asarray(output["RH"]) / 100 51 RHi = np.asarray(output["RH_ice"]) / 100 52 qc = np.asarray(output["LWC"]) 53 qi = np.asarray(output["IWC"]) 54 qv = np.asarray(output["qv"]) 55 T_frz = np.asarray(output["T_frz"]) 56 if show_conc: 57 nc, ni = np.asarray(output["ns"]), np.asarray(output["ni"]) 58 59 if t_lim is None: 60 t_lim = np.amax(time) 61 62 first_ice_idx = np.where(qi > 1e-7)[0][0] 63 first_ice_time = time[first_ice_idx] 64 first_T_frz = T[first_ice_idx] 65 66 if not show_tf: 67 rc, ri = np.asarray(output["rs"]), np.asarray(output["ri"]) 68 69 if show_jhom: 70 svp = Formulae( 71 saturation_vapour_pressure="FlatauWalkoCotton" 72 ).saturation_vapour_pressure 73 a_w_ice = svp.pvs_ice(T) / svp.pvs_water(T) 74 d_a_w_ice = (RHi - 1) * a_w_ice 75 76 j_hom_rate = Formulae( 77 homogeneous_ice_nucleation_rate="KoopMurray2016" 78 ).homogeneous_ice_nucleation_rate 79 koop_murray_2016 = j_hom_rate.j_hom(T, d_a_w_ice) 80 j_hom_rate = Formulae( 81 homogeneous_ice_nucleation_rate="KoopMurray2016_DWA" 82 ).homogeneous_ice_nucleation_rate 83 KoopMurray2016_DWA = j_hom_rate.j_hom(T, d_a_w_ice) 84 abs_diff_j_hom = (KoopMurray2016_DWA - koop_murray_2016) / koop_murray_2016 85 else: 86 radius = np.asarray(output["radius"]) 87 multiplicity = np.asarray(output["multiplicity"]) 88 _, axs = pyplot.subplots(1, 4, figsize=(20, 5), constrained_layout=True) 89 90 # Temperture profile 91 iax = 0 92 ax = axs[iax] 93 ax.plot(time, RH, color="red", linestyle="dashdot", label=r"$S_\text{w}$") 94 ax.plot(time, RHi, color="blue", linestyle="--", label=r"$S_\text{i}$") 95 ax.set_xlabel("time [s]", fontsize=ax_lab_fsize2) 96 ax.set_ylabel("saturation ratio", fontsize=ax_lab_fsize2) 97 ax.legend(loc="center left", fontsize=ax_lab_fsize2) 98 ax.set_xlim(time[0], t_lim) 99 ax.tick_params(labelsize=tick_fsize) 100 ax.set_title(title_add[iax] + r"ambient thermodynamics", fontsize=ax_lab_fsize) 101 ax.grid(visible=True) 102 ax.axvline(x=first_ice_time, color="black", linestyle=":") 103 104 twin = ax.twinx() 105 twin.plot(time, T, color="black", linestyle="-", label="T") 106 twin.set_xlabel("time [s]", fontsize=ax_lab_fsize2) 107 twin.set_ylabel("temperature [K]", fontsize=ax_lab_fsize2) 108 twin.legend(loc="upper left", fontsize=ax_lab_fsize2) 109 twin.tick_params(labelsize=tick_fsize) 110 111 # mixing ratio and number concentration 112 iax = 1 113 ax = axs[iax] 114 ax.plot(time, qc, color="red", linestyle="dashdot", label=r"$q_\text{w}$") 115 ax.plot(time, qi, color="blue", linestyle="--", label=r"$q_\text{i}$") 116 ax.plot(time, qv, color="black", linestyle="-", label=r"$q_\text{v}$") 117 ax.set_yscale("log") 118 ax.set_ylim(1e-5, 1e-2) 119 ax.set_xlabel("time [s]", fontsize=ax_lab_fsize2) 120 ax.set_ylabel(r"mixing ratio [$\mathrm{kg \, kg^{-1}}$]", fontsize=ax_lab_fsize2) 121 ax.legend(fontsize=ax_lab_fsize2) 122 ax.tick_params(labelsize=tick_fsize) 123 ax.set_xlim(time[0], t_lim) 124 ax.grid(visible=True) 125 ax.axvline(x=first_ice_time, color="black", linestyle=":") 126 ax.set_title(title_add[iax] + r"bulk quantities", fontsize=ax_lab_fsize) 127 if show_conc: 128 twin = ax.twinx() 129 twin.plot(time, nc, color="red", linestyle="densly dashdot", label="water") 130 twin.plot(time, ni, color="blue", linestyle="densly dashed", label="ice") 131 twin.set_yscale("log") 132 twin.set_xlabel("time [s]", fontsize=ax_lab_fsize2) 133 twin.set_ylabel( 134 r"number concentration [$\mathrm{kg^{-1}}$]", fontsize=ax_lab_fsize2 135 ) 136 twin.tick_params(labelsize=tick_fsize) 137 # tfrz histogram or mean radius 138 iax = 2 139 ax = axs[iax] 140 if show_tf: 141 ax.hist( 142 T_frz, 143 bins=T_frz_bins_kelvin, 144 density=True, 145 cumulative=-1, 146 alpha=1.0, 147 histtype="step", 148 linewidth=1.5, 149 ) 150 ax.set_xlim(left=234, right=239) 151 ax.set_xlabel("freezing temperature [K]", fontsize=ax_lab_fsize2) 152 ax.set_ylabel("frozen fraction", fontsize=ax_lab_fsize2) 153 ax.tick_params(labelsize=tick_fsize) 154 ax.grid(visible=True) 155 ax.axvline(x=first_T_frz, color="black", linestyle=":") 156 ax.set_title(title_add[iax] + r"$T_{frz}$ histogram", fontsize=ax_lab_fsize) 157 else: 158 ax.plot(time, rc * 1e6, color="red", linestyle="dashdot", label="water") 159 ax.plot(time, ri * 1e6, color="blue", linestyle="--", label="ice") 160 ax.set_yscale("log") 161 ax.set_ylim(1e-2, 1e2) 162 ax.set_xlabel("time [s]", fontsize=ax_lab_fsize2) 163 ax.set_ylabel("mean radius [µm]", fontsize=ax_lab_fsize2) 164 ax.legend(fontsize=ax_lab_fsize2) 165 ax.set_xlim(time[0], t_lim) 166 ax.tick_params(labelsize=tick_fsize) 167 ax.grid(visible=True) 168 ax.axvline(x=first_ice_time, color="black", linestyle=":") 169 ax.set_title(title_add[iax] + r" mean radius", fontsize=ax_lab_fsize) 170 # Water activity difference profile 171 iax = 3 172 ax = axs[iax] 173 if show_jhom: 174 lin_s_SP2023 = "--" 175 lin_s_KM2016 = "-" 176 if simulation["settings"]["hom_freezing"] == "KoopMurray2016_DWA": 177 lin_s_SP2023 = "-" 178 lin_s_KM2016 = "--" 179 180 ax.plot( 181 time, 182 koop_murray_2016, 183 color="blue", 184 linestyle=lin_s_KM2016, 185 label="JHOM-T", 186 ) 187 ax.plot( 188 time, 189 KoopMurray2016_DWA, 190 color="red", 191 linestyle=lin_s_SP2023, 192 label="JHOM-DWA", 193 ) 194 ax.set_ylabel( 195 r"nucleation rate [$\mathrm{m^{-3} \, s^{-1}}$]", fontsize=ax_lab_fsize2 196 ) 197 ax.set_ylim(1e-30, 1e30) 198 ax.set_title(title_add[iax] + r"nucleation rates", fontsize=ax_lab_fsize) 199 ax.legend(loc="upper left", fontsize=ax_lab_fsize2) 200 ax.set_yscale("log") 201 ax.set_xlim(time[0], t_lim) 202 ax.axvline(x=first_ice_time, color="black", linestyle=":") 203 ax.set_xlabel("time [s]", fontsize=ax_lab_fsize2) 204 if plot_daw: 205 twin = ax.twinx() 206 twin.plot( 207 time, 208 d_a_w_ice, 209 color="gray", 210 linestyle="dashdot", 211 label=r"$\Delta a_{w}$", 212 ) 213 twin.set_ylim(0.2, 0.35) 214 twin.set_ylabel("water activity difference", fontsize=ax_lab_fsize2) 215 twin.tick_params(labelsize=tick_fsize) 216 twin.legend(loc="lower left", fontsize=ax_lab_fsize2) 217 else: 218 twin = ax.twinx() 219 twin.plot( 220 time, 221 abs_diff_j_hom, 222 label=r"$\Delta J_{\mathrm{hom}}$", 223 color="gray", 224 linestyle="dashdot", 225 ) 226 twin.set_ylim(-10, 10) 227 twin.set_ylabel("relative error", fontsize=ax_lab_fsize2) 228 twin.tick_params(labelsize=tick_fsize) 229 twin.legend(loc="lower left", fontsize=ax_lab_fsize2) 230 else: 231 ax.scatter(radius * 1e6, multiplicity) 232 ax.set_yscale("log") 233 ax.set_xscale("log") 234 ax.set_xlim(1e-3, 5e-0) 235 ax.set_xlabel("initial radius [µm]", fontsize=ax_lab_fsize2) 236 ax.set_ylabel("multiplicity", fontsize=ax_lab_fsize2) 237 ax.set_title(title_add[iax] + r"CCN size distribution", fontsize=ax_lab_fsize) 238 239 ax.tick_params(labelsize=tick_fsize) 240 ax.grid(visible=True) 241 242 243def plot_freezing_temperatures_histogram(ax, simulation, plot_rhi=False): 244 number_of_ensemble_runs = simulation["settings"]["number_of_ensemble_runs"] 245 246 for i in range(number_of_ensemble_runs): 247 output = simulation["ensemble_member_outputs"][i] 248 if plot_rhi: 249 var = np.asarray(output["RHi_frz"]) 250 bins = np.linspace(1, 1.6, num=60, endpoint=True) 251 cumulative = 1 252 else: 253 var = np.asarray(output["T_frz"]) 254 bins = T_frz_bins_kelvin 255 cumulative = -1 256 title = "Nucleation rate=" + simulation["settings"]["hom_freezing"] 257 258 ax.hist( 259 var, 260 bins=bins, 261 density=True, 262 cumulative=cumulative, 263 alpha=1.0, 264 histtype="step", 265 linewidth=1.5, 266 ) 267 268 if plot_rhi: 269 ax.set_xlim(left=1.0, right=1.6) 270 ax.set_xlabel("freezing supersaturation wrt ice", fontsize=ax_lab_fsize) 271 else: 272 ax.set_xlim(left=234, right=239) 273 ax.axvline(x=235, color="k", linestyle="--") 274 ax.set_xlabel("freezing temperature [K]", fontsize=ax_lab_fsize) 275 ax.set_title(title, fontsize=ax_lab_fsize) 276 ax.set_ylabel("frozen fraction", fontsize=ax_lab_fsize) 277 ax.tick_params(labelsize=tick_fsize) 278 return ax 279 280 281def plot_freezing_temperatures_histogram_allinone( 282 ax, simulations, title, lloc="upper right" 283): 284 285 colors = ["black", "blue", "red"] 286 linestyles = ["-", "--", ":"] 287 288 for k, simulation in enumerate(simulations): 289 290 number_of_ensemble_runs = simulation["settings"]["number_of_ensemble_runs"] 291 n_sd = simulation["settings"]["n_sd"] 292 histogram_list = np.zeros((number_of_ensemble_runs, len(T_frz_bins_kelvin) - 1)) 293 for i in range(number_of_ensemble_runs): 294 output = simulation["ensemble_member_outputs"][i] 295 T_frz = np.asarray(output["T_frz"]) 296 297 hist, T_frz_bins_center = cumulative_histogram( 298 T_frz, T_frz_bins_kelvin, reverse=True 299 ) 300 histogram_list[i, :] = hist 301 302 max_line = np.max(histogram_list, axis=0) 303 mean_line = np.mean(histogram_list, axis=0) 304 min_line = np.min(histogram_list, axis=0) 305 306 ax.plot( 307 T_frz_bins_center, 308 mean_line, 309 color=colors[k], 310 linestyle=linestyles[k], 311 label=r"$n_\text{sd}$: " + f"{int(n_sd):5.0f}", 312 ) 313 ax.fill_between( 314 T_frz_bins_center, min_line, max_line, color=colors[k], alpha=0.2 315 ) 316 ax.set_xlim(left=234.5, right=239) 317 ax.set_title(title, fontsize=ax_lab_fsize) 318 ax.set_xlabel("freezing temperature [K]", fontsize=ax_lab_fsize) 319 ax.set_ylabel("frozen fraction", fontsize=ax_lab_fsize) 320 ax.tick_params(labelsize=tick_fsize) 321 ax.xaxis.set_major_locator(ticker.MultipleLocator(1)) 322 ax.xaxis.set_minor_locator(ticker.MultipleLocator(0.5)) 323 ax.grid(visible=True) 324 ax.legend(loc=lloc, fontsize=ax_lab_fsize) 325 return ax 326 327 328def plot_freezing_temperatures_2d_histogram_seaborn( 329 ensemble_simulations, 330 hom_freezing_type, 331 title="", 332 height=4, 333 width=5, 334): 335 sns.set_theme(style="ticks") 336 337 second_axis = True 338 y_log = True 339 340 ens_variable_name = ensemble_simulations["ens_variable_name"] 341 simulations = ensemble_simulations[hom_freezing_type] 342 343 ens_variable = [] 344 ens_variable_sec = [] 345 T_frz_hist = [] 346 multiplicity_hist = [] 347 348 for simulation in simulations: 349 if ens_variable_name == "sig": 350 ens_variable_value = simulation["settings"]["sigma_droplet_distribution"] 351 else: 352 ens_variable_value = simulation["settings"][ens_variable_name] 353 354 output = simulation["ensemble_member_outputs"][0] 355 qi = np.asarray(output["IWC"]) 356 T_frz = np.asarray(output["T_frz"]) 357 multiplicity = np.asarray(output["multiplicity"]) 358 first_ice_idx = np.where(qi > 1e-7)[0][0] 359 360 if ens_variable_name == "n_ccn": 361 rs = np.asarray(output["rs"]) 362 ens_variable_sec_value = rs[first_ice_idx - 1] 363 else: 364 time = np.asarray(output["t"]) 365 T = np.asarray(output["T"]) 366 ens_variable_sec_value = abs(np.mean(np.diff(T) / np.diff(time))) 367 368 T_frz_hist.extend(T_frz) 369 multiplicity_hist.extend(multiplicity) 370 ens_variable.extend(np.full_like(T_frz, ens_variable_value)) 371 ens_variable_sec.append(ens_variable_sec_value) 372 373 ens_variable = np.array(ens_variable) 374 ens_variable_sec = np.array(ens_variable_sec) 375 y_label, y_label_sec = "", "" 376 ylim = None 377 if ens_variable_name == "w_updraft": 378 y_label = r"vertical updraft [$\mathrm{m \, s^{-1}}$]" 379 y_label_sec = r"cooling rate [$\mathrm{mK \, s^{-1}}$]" 380 ens_variable_sec = ens_variable_sec * 1e3 381 binwidth = 0.25 382 elif ens_variable_name == "n_ccn": 383 y_label = r"ccn concentration at STP [$\mathrm{cm^{-3}}$]" 384 y_label_sec = r"radius [$\mathrm{\mu m}$]" 385 ens_variable = ens_variable / 1.0e6 386 ens_variable_sec = ens_variable_sec * 1.0e6 387 binwidth = 0.25 388 elif ens_variable_name == "sig": 389 y_label = r"standard deviation DSD" 390 second_axis = False 391 y_log = False 392 binwidth = 0.1 393 ylim = (1.5, 2.0) 394 ens_variable_label = np.unique(np.sort(ens_variable)) 395 396 xlim = (232.5, 240) 397 h = sns.JointGrid( 398 x=T_frz_hist, 399 y=ens_variable, 400 xlim=xlim, 401 ylim=ylim, 402 ) 403 if y_log: 404 h.ax_joint.set(yscale="log") 405 406 h.plot_joint( 407 sns.histplot, 408 stat="probability", 409 binwidth=binwidth, 410 discrete=(False, False), 411 weights=multiplicity_hist, 412 pmax=0.8, 413 ) 414 415 h.plot_marginals( 416 sns.histplot, 417 element="step", 418 ) 419 h.set_axis_labels("freezing temperature [K]", y_label, fontsize=ax_lab_fsize) 420 h.ax_joint.set_title( 421 title, 422 pad=50, 423 fontsize=ax_title_size, 424 ) 425 h.ax_joint.tick_params(labelsize=tick_fsize) 426 h.ax_joint.grid(True, axis="x", which="both", linestyle="-", linewidth=0.6) 427 h.ax_marg_x.tick_params(labelsize=tick_fsize) 428 h.ax_marg_y.tick_params(labelsize=tick_fsize) 429 430 h.ax_joint.xaxis.set_major_locator(ticker.MultipleLocator(2, offset=1)) 431 h.ax_joint.xaxis.set_minor_locator(ticker.MultipleLocator(1, offset=1)) 432 h.ax_marg_y.remove() 433 434 if second_axis: 435 ax2 = h.ax_joint.secondary_yaxis("right", functions=(lambda y: y, lambda y: y)) 436 ax2.set_yticks(ens_variable_label) 437 ax2.set_yticklabels([f"{v:.1f}" for v in ens_variable_sec]) 438 ax2.minorticks_off() 439 ax2.tick_params( 440 axis="y", 441 which="both", 442 direction="out", 443 length=h.ax_joint.yaxis.get_ticklines()[0].get_markersize(), 444 width=h.ax_joint.yaxis.get_ticklines()[0].get_markeredgewidth(), 445 labelsize=tick_fsize, 446 ) 447 ax2.set_ylabel(y_label_sec, fontsize=ax_lab_fsize) 448 449 h.fig.set_size_inches(width, height) 450 451 452def plot_ensemble_bulk( 453 ax, ensemble_simulations, var_name, title_add="" 454): # pylint: disable=too-many-nested-blocks 455 colors = ["blue", "red", "cyan"] 456 linestyles = ["-", "--", ":"] 457 pyplot.rcParams["axes.prop_cycle"] = cycler(color=colors) + cycler( 458 linestyle=linestyles 459 ) 460 461 for ensemble_simulation in ensemble_simulations: 462 ens_var = np.asarray(ensemble_simulation["ens_variable"]) 463 ens_var_name = ensemble_simulation["ens_variable_name"] 464 hom_freezing_types = ensemble_simulation["hom_freezing_types"] 465 hom_freezing_labels = ensemble_simulation["hom_freezing_labels"] 466 len_ens_var = len(ens_var) 467 468 if ens_var_name == "n_ccn": 469 ens_var_scale = 1.0 / 1e6 470 else: 471 ens_var_scale = 1.0 472 if ens_var_name == "sig": 473 ens_var_name = "sigma_droplet_distribution" 474 475 for j, hom_freezing_type in enumerate(hom_freezing_types): 476 simulations = ensemble_simulation[hom_freezing_type] 477 number_of_ensemble_runs = simulations[0]["settings"][ 478 "number_of_ensemble_runs" 479 ] 480 var = np.zeros((len_ens_var, number_of_ensemble_runs)) 481 for i in range(len_ens_var): 482 for simulation in simulations: 483 if simulation["settings"][ens_var_name] == ens_var[i]: 484 for h in range(number_of_ensemble_runs): 485 output = simulation["ensemble_member_outputs"][h] 486 if var_name == "freezing_fraction": 487 ni = np.asarray(output["ni"])[-1] 488 nc = np.asarray(output["ns"])[0] 489 var[i, h] = (1 - (nc - ni) / nc) * 100 490 else: 491 var[i, h] = np.asarray(output[var_name])[-1] 492 493 if number_of_ensemble_runs > 1: 494 ax.fill_betweenx( 495 ens_var * ens_var_scale, 496 np.min(var, axis=1), 497 np.max(var, axis=1), 498 alpha=0.2, 499 ) 500 else: 501 ax.scatter( 502 var[:, 0], 503 ens_var * ens_var_scale, 504 ) 505 ax.plot( 506 var[:, 0], 507 ens_var * ens_var_scale, 508 label=hom_freezing_labels[j], 509 ) 510 511 title, x_label, y_label, ens_label = "", "", "", "" 512 if var_name == "ni": 513 ax.set_xscale("log") 514 x_label = r"ice number concentration [$\mathrm{kg^{-1}}$]" 515 title = r"$n_\text{i}$" 516 ax.set_xlim(1e6, 1e10) 517 elif var_name == "IWC": 518 ax.set_xscale("linear") 519 x_label = r"mixing ratio [$\mathrm{kg \, kg^{-1}}$]" 520 title = "ice mixing ratio" 521 ax.set_xlim(3e-4, 1e-3) 522 elif var_name == "freezing_fraction": 523 title = r"$n_\text{frz}$" 524 x_label = r"frozen fraction [$\mathrm{\%}$]" 525 ax.set_xlim(0, 20) 526 527 if ens_var_name == "n_ccn": 528 ax.set_yscale("log") 529 y_label = r"ccn concentration [$\mathrm{cm^{-3}}$]" 530 ens_label = r"$n_\text{ccn}$ ensemble" 531 elif ens_var_name == "w_updraft": 532 ax.set_yscale("linear") 533 y_label = r"vertical updraft [$\mathrm{m \, s^{-1}}$]" 534 ens_label = "w ensemble" 535 elif ens_var_name in ("sigma_droplet_distribution", "sig"): 536 y_label = r"standard deviation DSD" 537 ens_label = r"$\sigma$ ensemble" 538 elif ens_var_name == "n_sd": 539 ax.set_yscale("log") 540 y_label = "number of super-particles" 541 ens_label = r"$n_\text{sd}$ ensemble" 542 543 ax.set_title(title_add + " " + title + " for " + ens_label, fontsize=ax_lab_fsize) 544 ax.set_xlabel(x_label, fontsize=ax_lab_fsize2) 545 ax.set_ylabel(y_label, fontsize=ax_lab_fsize2) 546 ax.grid(visible=True) 547 ax.legend(fontsize=ax_lab_fsize) 548 ax.tick_params(labelsize=tick_fsize + 2) 549 return ax
formulae =
<PySDM.formulae.Formulae object>
ax_title_size =
18
ax_lab_fsize =
15
ax_lab_fsize2 =
19
tick_fsize =
15
T_frz_bins =
array([-40. , -39.89830508, -39.79661017, -39.69491525,
-39.59322034, -39.49152542, -39.38983051, -39.28813559,
-39.18644068, -39.08474576, -38.98305085, -38.88135593,
-38.77966102, -38.6779661 , -38.57627119, -38.47457627,
-38.37288136, -38.27118644, -38.16949153, -38.06779661,
-37.96610169, -37.86440678, -37.76271186, -37.66101695,
-37.55932203, -37.45762712, -37.3559322 , -37.25423729,
-37.15254237, -37.05084746, -36.94915254, -36.84745763,
-36.74576271, -36.6440678 , -36.54237288, -36.44067797,
-36.33898305, -36.23728814, -36.13559322, -36.03389831,
-35.93220339, -35.83050847, -35.72881356, -35.62711864,
-35.52542373, -35.42372881, -35.3220339 , -35.22033898,
-35.11864407, -35.01694915, -34.91525424, -34.81355932,
-34.71186441, -34.61016949, -34.50847458, -34.40677966,
-34.30508475, -34.20338983, -34.10169492, -34. ])
T_frz_bins_kelvin =
array([230. , 230.1010101 , 230.2020202 , 230.3030303 ,
230.4040404 , 230.50505051, 230.60606061, 230.70707071,
230.80808081, 230.90909091, 231.01010101, 231.11111111,
231.21212121, 231.31313131, 231.41414141, 231.51515152,
231.61616162, 231.71717172, 231.81818182, 231.91919192,
232.02020202, 232.12121212, 232.22222222, 232.32323232,
232.42424242, 232.52525253, 232.62626263, 232.72727273,
232.82828283, 232.92929293, 233.03030303, 233.13131313,
233.23232323, 233.33333333, 233.43434343, 233.53535354,
233.63636364, 233.73737374, 233.83838384, 233.93939394,
234.04040404, 234.14141414, 234.24242424, 234.34343434,
234.44444444, 234.54545455, 234.64646465, 234.74747475,
234.84848485, 234.94949495, 235.05050505, 235.15151515,
235.25252525, 235.35353535, 235.45454545, 235.55555556,
235.65656566, 235.75757576, 235.85858586, 235.95959596,
236.06060606, 236.16161616, 236.26262626, 236.36363636,
236.46464646, 236.56565657, 236.66666667, 236.76767677,
236.86868687, 236.96969697, 237.07070707, 237.17171717,
237.27272727, 237.37373737, 237.47474747, 237.57575758,
237.67676768, 237.77777778, 237.87878788, 237.97979798,
238.08080808, 238.18181818, 238.28282828, 238.38383838,
238.48484848, 238.58585859, 238.68686869, 238.78787879,
238.88888889, 238.98989899, 239.09090909, 239.19191919,
239.29292929, 239.39393939, 239.49494949, 239.5959596 ,
239.6969697 , 239.7979798 , 239.8989899 , 240. ])
def
cumulative_histogram(data, bins, reverse=False, density=True):
20def cumulative_histogram(data, bins, reverse=False, density=True): 21 hist, bin_edges = np.histogram(data, bins=bins, density=False) 22 23 if reverse: 24 cum_hist = np.cumsum(hist[::-1])[::-1] 25 cum_hist_0 = cum_hist[0] 26 else: 27 cum_hist = np.cumsum(hist) 28 cum_hist_0 = cum_hist[-1] 29 30 if density: 31 cum_hist = cum_hist / cum_hist_0 32 33 bin_centers = 0.5 * (bin_edges[:-1] + bin_edges[1:]) 34 return cum_hist, bin_centers
def
plot_thermodynamics_and_bulk( simulation, title_add=None, show_conc=False, show_jhom=True, show_tf=True, t_lim=None):
37def plot_thermodynamics_and_bulk( 38 simulation, 39 title_add=None, 40 show_conc=False, 41 show_jhom=True, 42 show_tf=True, 43 t_lim=None, 44): 45 if title_add is None: 46 title_add = ["", "", "", ""] 47 plot_daw = False 48 output = simulation["ensemble_member_outputs"][0] 49 time = output["t"] 50 T = np.asarray(output["T"]) 51 RH = np.asarray(output["RH"]) / 100 52 RHi = np.asarray(output["RH_ice"]) / 100 53 qc = np.asarray(output["LWC"]) 54 qi = np.asarray(output["IWC"]) 55 qv = np.asarray(output["qv"]) 56 T_frz = np.asarray(output["T_frz"]) 57 if show_conc: 58 nc, ni = np.asarray(output["ns"]), np.asarray(output["ni"]) 59 60 if t_lim is None: 61 t_lim = np.amax(time) 62 63 first_ice_idx = np.where(qi > 1e-7)[0][0] 64 first_ice_time = time[first_ice_idx] 65 first_T_frz = T[first_ice_idx] 66 67 if not show_tf: 68 rc, ri = np.asarray(output["rs"]), np.asarray(output["ri"]) 69 70 if show_jhom: 71 svp = Formulae( 72 saturation_vapour_pressure="FlatauWalkoCotton" 73 ).saturation_vapour_pressure 74 a_w_ice = svp.pvs_ice(T) / svp.pvs_water(T) 75 d_a_w_ice = (RHi - 1) * a_w_ice 76 77 j_hom_rate = Formulae( 78 homogeneous_ice_nucleation_rate="KoopMurray2016" 79 ).homogeneous_ice_nucleation_rate 80 koop_murray_2016 = j_hom_rate.j_hom(T, d_a_w_ice) 81 j_hom_rate = Formulae( 82 homogeneous_ice_nucleation_rate="KoopMurray2016_DWA" 83 ).homogeneous_ice_nucleation_rate 84 KoopMurray2016_DWA = j_hom_rate.j_hom(T, d_a_w_ice) 85 abs_diff_j_hom = (KoopMurray2016_DWA - koop_murray_2016) / koop_murray_2016 86 else: 87 radius = np.asarray(output["radius"]) 88 multiplicity = np.asarray(output["multiplicity"]) 89 _, axs = pyplot.subplots(1, 4, figsize=(20, 5), constrained_layout=True) 90 91 # Temperture profile 92 iax = 0 93 ax = axs[iax] 94 ax.plot(time, RH, color="red", linestyle="dashdot", label=r"$S_\text{w}$") 95 ax.plot(time, RHi, color="blue", linestyle="--", label=r"$S_\text{i}$") 96 ax.set_xlabel("time [s]", fontsize=ax_lab_fsize2) 97 ax.set_ylabel("saturation ratio", fontsize=ax_lab_fsize2) 98 ax.legend(loc="center left", fontsize=ax_lab_fsize2) 99 ax.set_xlim(time[0], t_lim) 100 ax.tick_params(labelsize=tick_fsize) 101 ax.set_title(title_add[iax] + r"ambient thermodynamics", fontsize=ax_lab_fsize) 102 ax.grid(visible=True) 103 ax.axvline(x=first_ice_time, color="black", linestyle=":") 104 105 twin = ax.twinx() 106 twin.plot(time, T, color="black", linestyle="-", label="T") 107 twin.set_xlabel("time [s]", fontsize=ax_lab_fsize2) 108 twin.set_ylabel("temperature [K]", fontsize=ax_lab_fsize2) 109 twin.legend(loc="upper left", fontsize=ax_lab_fsize2) 110 twin.tick_params(labelsize=tick_fsize) 111 112 # mixing ratio and number concentration 113 iax = 1 114 ax = axs[iax] 115 ax.plot(time, qc, color="red", linestyle="dashdot", label=r"$q_\text{w}$") 116 ax.plot(time, qi, color="blue", linestyle="--", label=r"$q_\text{i}$") 117 ax.plot(time, qv, color="black", linestyle="-", label=r"$q_\text{v}$") 118 ax.set_yscale("log") 119 ax.set_ylim(1e-5, 1e-2) 120 ax.set_xlabel("time [s]", fontsize=ax_lab_fsize2) 121 ax.set_ylabel(r"mixing ratio [$\mathrm{kg \, kg^{-1}}$]", fontsize=ax_lab_fsize2) 122 ax.legend(fontsize=ax_lab_fsize2) 123 ax.tick_params(labelsize=tick_fsize) 124 ax.set_xlim(time[0], t_lim) 125 ax.grid(visible=True) 126 ax.axvline(x=first_ice_time, color="black", linestyle=":") 127 ax.set_title(title_add[iax] + r"bulk quantities", fontsize=ax_lab_fsize) 128 if show_conc: 129 twin = ax.twinx() 130 twin.plot(time, nc, color="red", linestyle="densly dashdot", label="water") 131 twin.plot(time, ni, color="blue", linestyle="densly dashed", label="ice") 132 twin.set_yscale("log") 133 twin.set_xlabel("time [s]", fontsize=ax_lab_fsize2) 134 twin.set_ylabel( 135 r"number concentration [$\mathrm{kg^{-1}}$]", fontsize=ax_lab_fsize2 136 ) 137 twin.tick_params(labelsize=tick_fsize) 138 # tfrz histogram or mean radius 139 iax = 2 140 ax = axs[iax] 141 if show_tf: 142 ax.hist( 143 T_frz, 144 bins=T_frz_bins_kelvin, 145 density=True, 146 cumulative=-1, 147 alpha=1.0, 148 histtype="step", 149 linewidth=1.5, 150 ) 151 ax.set_xlim(left=234, right=239) 152 ax.set_xlabel("freezing temperature [K]", fontsize=ax_lab_fsize2) 153 ax.set_ylabel("frozen fraction", fontsize=ax_lab_fsize2) 154 ax.tick_params(labelsize=tick_fsize) 155 ax.grid(visible=True) 156 ax.axvline(x=first_T_frz, color="black", linestyle=":") 157 ax.set_title(title_add[iax] + r"$T_{frz}$ histogram", fontsize=ax_lab_fsize) 158 else: 159 ax.plot(time, rc * 1e6, color="red", linestyle="dashdot", label="water") 160 ax.plot(time, ri * 1e6, color="blue", linestyle="--", label="ice") 161 ax.set_yscale("log") 162 ax.set_ylim(1e-2, 1e2) 163 ax.set_xlabel("time [s]", fontsize=ax_lab_fsize2) 164 ax.set_ylabel("mean radius [µm]", fontsize=ax_lab_fsize2) 165 ax.legend(fontsize=ax_lab_fsize2) 166 ax.set_xlim(time[0], t_lim) 167 ax.tick_params(labelsize=tick_fsize) 168 ax.grid(visible=True) 169 ax.axvline(x=first_ice_time, color="black", linestyle=":") 170 ax.set_title(title_add[iax] + r" mean radius", fontsize=ax_lab_fsize) 171 # Water activity difference profile 172 iax = 3 173 ax = axs[iax] 174 if show_jhom: 175 lin_s_SP2023 = "--" 176 lin_s_KM2016 = "-" 177 if simulation["settings"]["hom_freezing"] == "KoopMurray2016_DWA": 178 lin_s_SP2023 = "-" 179 lin_s_KM2016 = "--" 180 181 ax.plot( 182 time, 183 koop_murray_2016, 184 color="blue", 185 linestyle=lin_s_KM2016, 186 label="JHOM-T", 187 ) 188 ax.plot( 189 time, 190 KoopMurray2016_DWA, 191 color="red", 192 linestyle=lin_s_SP2023, 193 label="JHOM-DWA", 194 ) 195 ax.set_ylabel( 196 r"nucleation rate [$\mathrm{m^{-3} \, s^{-1}}$]", fontsize=ax_lab_fsize2 197 ) 198 ax.set_ylim(1e-30, 1e30) 199 ax.set_title(title_add[iax] + r"nucleation rates", fontsize=ax_lab_fsize) 200 ax.legend(loc="upper left", fontsize=ax_lab_fsize2) 201 ax.set_yscale("log") 202 ax.set_xlim(time[0], t_lim) 203 ax.axvline(x=first_ice_time, color="black", linestyle=":") 204 ax.set_xlabel("time [s]", fontsize=ax_lab_fsize2) 205 if plot_daw: 206 twin = ax.twinx() 207 twin.plot( 208 time, 209 d_a_w_ice, 210 color="gray", 211 linestyle="dashdot", 212 label=r"$\Delta a_{w}$", 213 ) 214 twin.set_ylim(0.2, 0.35) 215 twin.set_ylabel("water activity difference", fontsize=ax_lab_fsize2) 216 twin.tick_params(labelsize=tick_fsize) 217 twin.legend(loc="lower left", fontsize=ax_lab_fsize2) 218 else: 219 twin = ax.twinx() 220 twin.plot( 221 time, 222 abs_diff_j_hom, 223 label=r"$\Delta J_{\mathrm{hom}}$", 224 color="gray", 225 linestyle="dashdot", 226 ) 227 twin.set_ylim(-10, 10) 228 twin.set_ylabel("relative error", fontsize=ax_lab_fsize2) 229 twin.tick_params(labelsize=tick_fsize) 230 twin.legend(loc="lower left", fontsize=ax_lab_fsize2) 231 else: 232 ax.scatter(radius * 1e6, multiplicity) 233 ax.set_yscale("log") 234 ax.set_xscale("log") 235 ax.set_xlim(1e-3, 5e-0) 236 ax.set_xlabel("initial radius [µm]", fontsize=ax_lab_fsize2) 237 ax.set_ylabel("multiplicity", fontsize=ax_lab_fsize2) 238 ax.set_title(title_add[iax] + r"CCN size distribution", fontsize=ax_lab_fsize) 239 240 ax.tick_params(labelsize=tick_fsize) 241 ax.grid(visible=True)
def
plot_freezing_temperatures_histogram(ax, simulation, plot_rhi=False):
244def plot_freezing_temperatures_histogram(ax, simulation, plot_rhi=False): 245 number_of_ensemble_runs = simulation["settings"]["number_of_ensemble_runs"] 246 247 for i in range(number_of_ensemble_runs): 248 output = simulation["ensemble_member_outputs"][i] 249 if plot_rhi: 250 var = np.asarray(output["RHi_frz"]) 251 bins = np.linspace(1, 1.6, num=60, endpoint=True) 252 cumulative = 1 253 else: 254 var = np.asarray(output["T_frz"]) 255 bins = T_frz_bins_kelvin 256 cumulative = -1 257 title = "Nucleation rate=" + simulation["settings"]["hom_freezing"] 258 259 ax.hist( 260 var, 261 bins=bins, 262 density=True, 263 cumulative=cumulative, 264 alpha=1.0, 265 histtype="step", 266 linewidth=1.5, 267 ) 268 269 if plot_rhi: 270 ax.set_xlim(left=1.0, right=1.6) 271 ax.set_xlabel("freezing supersaturation wrt ice", fontsize=ax_lab_fsize) 272 else: 273 ax.set_xlim(left=234, right=239) 274 ax.axvline(x=235, color="k", linestyle="--") 275 ax.set_xlabel("freezing temperature [K]", fontsize=ax_lab_fsize) 276 ax.set_title(title, fontsize=ax_lab_fsize) 277 ax.set_ylabel("frozen fraction", fontsize=ax_lab_fsize) 278 ax.tick_params(labelsize=tick_fsize) 279 return ax
def
plot_freezing_temperatures_histogram_allinone(ax, simulations, title, lloc='upper right'):
282def plot_freezing_temperatures_histogram_allinone( 283 ax, simulations, title, lloc="upper right" 284): 285 286 colors = ["black", "blue", "red"] 287 linestyles = ["-", "--", ":"] 288 289 for k, simulation in enumerate(simulations): 290 291 number_of_ensemble_runs = simulation["settings"]["number_of_ensemble_runs"] 292 n_sd = simulation["settings"]["n_sd"] 293 histogram_list = np.zeros((number_of_ensemble_runs, len(T_frz_bins_kelvin) - 1)) 294 for i in range(number_of_ensemble_runs): 295 output = simulation["ensemble_member_outputs"][i] 296 T_frz = np.asarray(output["T_frz"]) 297 298 hist, T_frz_bins_center = cumulative_histogram( 299 T_frz, T_frz_bins_kelvin, reverse=True 300 ) 301 histogram_list[i, :] = hist 302 303 max_line = np.max(histogram_list, axis=0) 304 mean_line = np.mean(histogram_list, axis=0) 305 min_line = np.min(histogram_list, axis=0) 306 307 ax.plot( 308 T_frz_bins_center, 309 mean_line, 310 color=colors[k], 311 linestyle=linestyles[k], 312 label=r"$n_\text{sd}$: " + f"{int(n_sd):5.0f}", 313 ) 314 ax.fill_between( 315 T_frz_bins_center, min_line, max_line, color=colors[k], alpha=0.2 316 ) 317 ax.set_xlim(left=234.5, right=239) 318 ax.set_title(title, fontsize=ax_lab_fsize) 319 ax.set_xlabel("freezing temperature [K]", fontsize=ax_lab_fsize) 320 ax.set_ylabel("frozen fraction", fontsize=ax_lab_fsize) 321 ax.tick_params(labelsize=tick_fsize) 322 ax.xaxis.set_major_locator(ticker.MultipleLocator(1)) 323 ax.xaxis.set_minor_locator(ticker.MultipleLocator(0.5)) 324 ax.grid(visible=True) 325 ax.legend(loc=lloc, fontsize=ax_lab_fsize) 326 return ax
def
plot_freezing_temperatures_2d_histogram_seaborn(ensemble_simulations, hom_freezing_type, title='', height=4, width=5):
329def plot_freezing_temperatures_2d_histogram_seaborn( 330 ensemble_simulations, 331 hom_freezing_type, 332 title="", 333 height=4, 334 width=5, 335): 336 sns.set_theme(style="ticks") 337 338 second_axis = True 339 y_log = True 340 341 ens_variable_name = ensemble_simulations["ens_variable_name"] 342 simulations = ensemble_simulations[hom_freezing_type] 343 344 ens_variable = [] 345 ens_variable_sec = [] 346 T_frz_hist = [] 347 multiplicity_hist = [] 348 349 for simulation in simulations: 350 if ens_variable_name == "sig": 351 ens_variable_value = simulation["settings"]["sigma_droplet_distribution"] 352 else: 353 ens_variable_value = simulation["settings"][ens_variable_name] 354 355 output = simulation["ensemble_member_outputs"][0] 356 qi = np.asarray(output["IWC"]) 357 T_frz = np.asarray(output["T_frz"]) 358 multiplicity = np.asarray(output["multiplicity"]) 359 first_ice_idx = np.where(qi > 1e-7)[0][0] 360 361 if ens_variable_name == "n_ccn": 362 rs = np.asarray(output["rs"]) 363 ens_variable_sec_value = rs[first_ice_idx - 1] 364 else: 365 time = np.asarray(output["t"]) 366 T = np.asarray(output["T"]) 367 ens_variable_sec_value = abs(np.mean(np.diff(T) / np.diff(time))) 368 369 T_frz_hist.extend(T_frz) 370 multiplicity_hist.extend(multiplicity) 371 ens_variable.extend(np.full_like(T_frz, ens_variable_value)) 372 ens_variable_sec.append(ens_variable_sec_value) 373 374 ens_variable = np.array(ens_variable) 375 ens_variable_sec = np.array(ens_variable_sec) 376 y_label, y_label_sec = "", "" 377 ylim = None 378 if ens_variable_name == "w_updraft": 379 y_label = r"vertical updraft [$\mathrm{m \, s^{-1}}$]" 380 y_label_sec = r"cooling rate [$\mathrm{mK \, s^{-1}}$]" 381 ens_variable_sec = ens_variable_sec * 1e3 382 binwidth = 0.25 383 elif ens_variable_name == "n_ccn": 384 y_label = r"ccn concentration at STP [$\mathrm{cm^{-3}}$]" 385 y_label_sec = r"radius [$\mathrm{\mu m}$]" 386 ens_variable = ens_variable / 1.0e6 387 ens_variable_sec = ens_variable_sec * 1.0e6 388 binwidth = 0.25 389 elif ens_variable_name == "sig": 390 y_label = r"standard deviation DSD" 391 second_axis = False 392 y_log = False 393 binwidth = 0.1 394 ylim = (1.5, 2.0) 395 ens_variable_label = np.unique(np.sort(ens_variable)) 396 397 xlim = (232.5, 240) 398 h = sns.JointGrid( 399 x=T_frz_hist, 400 y=ens_variable, 401 xlim=xlim, 402 ylim=ylim, 403 ) 404 if y_log: 405 h.ax_joint.set(yscale="log") 406 407 h.plot_joint( 408 sns.histplot, 409 stat="probability", 410 binwidth=binwidth, 411 discrete=(False, False), 412 weights=multiplicity_hist, 413 pmax=0.8, 414 ) 415 416 h.plot_marginals( 417 sns.histplot, 418 element="step", 419 ) 420 h.set_axis_labels("freezing temperature [K]", y_label, fontsize=ax_lab_fsize) 421 h.ax_joint.set_title( 422 title, 423 pad=50, 424 fontsize=ax_title_size, 425 ) 426 h.ax_joint.tick_params(labelsize=tick_fsize) 427 h.ax_joint.grid(True, axis="x", which="both", linestyle="-", linewidth=0.6) 428 h.ax_marg_x.tick_params(labelsize=tick_fsize) 429 h.ax_marg_y.tick_params(labelsize=tick_fsize) 430 431 h.ax_joint.xaxis.set_major_locator(ticker.MultipleLocator(2, offset=1)) 432 h.ax_joint.xaxis.set_minor_locator(ticker.MultipleLocator(1, offset=1)) 433 h.ax_marg_y.remove() 434 435 if second_axis: 436 ax2 = h.ax_joint.secondary_yaxis("right", functions=(lambda y: y, lambda y: y)) 437 ax2.set_yticks(ens_variable_label) 438 ax2.set_yticklabels([f"{v:.1f}" for v in ens_variable_sec]) 439 ax2.minorticks_off() 440 ax2.tick_params( 441 axis="y", 442 which="both", 443 direction="out", 444 length=h.ax_joint.yaxis.get_ticklines()[0].get_markersize(), 445 width=h.ax_joint.yaxis.get_ticklines()[0].get_markeredgewidth(), 446 labelsize=tick_fsize, 447 ) 448 ax2.set_ylabel(y_label_sec, fontsize=ax_lab_fsize) 449 450 h.fig.set_size_inches(width, height)
def
plot_ensemble_bulk(ax, ensemble_simulations, var_name, title_add=''):
453def plot_ensemble_bulk( 454 ax, ensemble_simulations, var_name, title_add="" 455): # pylint: disable=too-many-nested-blocks 456 colors = ["blue", "red", "cyan"] 457 linestyles = ["-", "--", ":"] 458 pyplot.rcParams["axes.prop_cycle"] = cycler(color=colors) + cycler( 459 linestyle=linestyles 460 ) 461 462 for ensemble_simulation in ensemble_simulations: 463 ens_var = np.asarray(ensemble_simulation["ens_variable"]) 464 ens_var_name = ensemble_simulation["ens_variable_name"] 465 hom_freezing_types = ensemble_simulation["hom_freezing_types"] 466 hom_freezing_labels = ensemble_simulation["hom_freezing_labels"] 467 len_ens_var = len(ens_var) 468 469 if ens_var_name == "n_ccn": 470 ens_var_scale = 1.0 / 1e6 471 else: 472 ens_var_scale = 1.0 473 if ens_var_name == "sig": 474 ens_var_name = "sigma_droplet_distribution" 475 476 for j, hom_freezing_type in enumerate(hom_freezing_types): 477 simulations = ensemble_simulation[hom_freezing_type] 478 number_of_ensemble_runs = simulations[0]["settings"][ 479 "number_of_ensemble_runs" 480 ] 481 var = np.zeros((len_ens_var, number_of_ensemble_runs)) 482 for i in range(len_ens_var): 483 for simulation in simulations: 484 if simulation["settings"][ens_var_name] == ens_var[i]: 485 for h in range(number_of_ensemble_runs): 486 output = simulation["ensemble_member_outputs"][h] 487 if var_name == "freezing_fraction": 488 ni = np.asarray(output["ni"])[-1] 489 nc = np.asarray(output["ns"])[0] 490 var[i, h] = (1 - (nc - ni) / nc) * 100 491 else: 492 var[i, h] = np.asarray(output[var_name])[-1] 493 494 if number_of_ensemble_runs > 1: 495 ax.fill_betweenx( 496 ens_var * ens_var_scale, 497 np.min(var, axis=1), 498 np.max(var, axis=1), 499 alpha=0.2, 500 ) 501 else: 502 ax.scatter( 503 var[:, 0], 504 ens_var * ens_var_scale, 505 ) 506 ax.plot( 507 var[:, 0], 508 ens_var * ens_var_scale, 509 label=hom_freezing_labels[j], 510 ) 511 512 title, x_label, y_label, ens_label = "", "", "", "" 513 if var_name == "ni": 514 ax.set_xscale("log") 515 x_label = r"ice number concentration [$\mathrm{kg^{-1}}$]" 516 title = r"$n_\text{i}$" 517 ax.set_xlim(1e6, 1e10) 518 elif var_name == "IWC": 519 ax.set_xscale("linear") 520 x_label = r"mixing ratio [$\mathrm{kg \, kg^{-1}}$]" 521 title = "ice mixing ratio" 522 ax.set_xlim(3e-4, 1e-3) 523 elif var_name == "freezing_fraction": 524 title = r"$n_\text{frz}$" 525 x_label = r"frozen fraction [$\mathrm{\%}$]" 526 ax.set_xlim(0, 20) 527 528 if ens_var_name == "n_ccn": 529 ax.set_yscale("log") 530 y_label = r"ccn concentration [$\mathrm{cm^{-3}}$]" 531 ens_label = r"$n_\text{ccn}$ ensemble" 532 elif ens_var_name == "w_updraft": 533 ax.set_yscale("linear") 534 y_label = r"vertical updraft [$\mathrm{m \, s^{-1}}$]" 535 ens_label = "w ensemble" 536 elif ens_var_name in ("sigma_droplet_distribution", "sig"): 537 y_label = r"standard deviation DSD" 538 ens_label = r"$\sigma$ ensemble" 539 elif ens_var_name == "n_sd": 540 ax.set_yscale("log") 541 y_label = "number of super-particles" 542 ens_label = r"$n_\text{sd}$ ensemble" 543 544 ax.set_title(title_add + " " + title + " for " + ens_label, fontsize=ax_lab_fsize) 545 ax.set_xlabel(x_label, fontsize=ax_lab_fsize2) 546 ax.set_ylabel(y_label, fontsize=ax_lab_fsize2) 547 ax.grid(visible=True) 548 ax.legend(fontsize=ax_lab_fsize) 549 ax.tick_params(labelsize=tick_fsize + 2) 550 return ax