Tutorial 03 — Gamma-PSD rain at C-band

Radar volumes sample thousands of drops. The observed Z_h, Z_dr, K_dp, and specific attenuation A_i are PSD-weighted integrals of the single-drop quantities. This notebook tabulates S(D) and Z(D) once, then sweeps the integrated observables across a range of normalised gamma PSDs parameterised by the median volume diameter D0.

import numpy as np
import matplotlib.pyplot as plt

from rustmatrix import Scatterer, radar, psd as rs_psd
from rustmatrix.tmatrix_aux import (dsr_thurai_2007, geom_horiz_back,
                                      geom_horiz_forw, K_w_sqr, wl_C)
from rustmatrix.refractive import m_w_10C

Tabulate once, evaluate many

s = Scatterer(wavelength=wl_C, m=m_w_10C[wl_C],
              Kw_sqr=K_w_sqr[wl_C], ddelt=1e-4, ndgs=2)
integ = rs_psd.PSDIntegrator()
integ.D_max = 8.0
integ.num_points = 64
integ.axis_ratio_func = lambda D: 1.0 / dsr_thurai_2007(D)
integ.geometries = (geom_horiz_back, geom_horiz_forw)
s.psd_integrator = integ
s.psd_integrator.init_scatter_table(s)

Sweep median diameter D0

D0s = np.linspace(0.5, 3.0, 12)
Zh = np.empty_like(D0s)
Zdr = np.empty_like(D0s)
Kdp = np.empty_like(D0s)
Ai = np.empty_like(D0s)

for i, D0 in enumerate(D0s):
    s.psd = rs_psd.GammaPSD(D0=D0, Nw=8e3, mu=4)
    s.set_geometry(geom_horiz_back)
    Zh[i] = 10 * np.log10(radar.refl(s))
    Zdr[i] = 10 * np.log10(radar.Zdr(s))
    s.set_geometry(geom_horiz_forw)
    Kdp[i] = radar.Kdp(s)
    Ai[i] = radar.Ai(s)

fig, axes = plt.subplots(2, 2, figsize=(9, 6), sharex=True)
axes[0, 0].plot(D0s, Zh, 'C0-o'); axes[0, 0].set_ylabel('Z_h [dBZ]')
axes[0, 1].plot(D0s, Zdr, 'C1-o'); axes[0, 1].set_ylabel('Z_dr [dB]')
axes[1, 0].plot(D0s, Kdp, 'C2-o'); axes[1, 0].set_ylabel('K_dp [°/km]')
axes[1, 1].plot(D0s, Ai, 'C3-o'); axes[1, 1].set_ylabel('A_i [dB/km]')
for ax in axes.flat:
    ax.set_xlabel('D0 [mm]')
    ax.grid(True, alpha=0.3)
fig.suptitle('C-band gamma-PSD rain observables (Nw=8e3, mu=4)')
fig.tight_layout();