{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Tutorial 05 \u2014 One PSD, six radar bands\n",
"\n",
"Retrieval algorithms that use more than one frequency need a\n",
"forward model that works across the whole instrument suite. This\n",
"notebook runs the same moderate-convective gamma PSD through\n",
"rustmatrix at S, C, X, Ku, Ka, and W band and plots the frequency\n",
"response of the dual-pol observables.\n"
]
},
{
"cell_type": "code",
"metadata": {},
"execution_count": null,
"outputs": [],
"source": [
"import numpy as np\n",
"import matplotlib.pyplot as plt\n",
"\n",
"from rustmatrix import Scatterer, radar, psd as rs_psd\n",
"from rustmatrix.tmatrix_aux import (dsr_thurai_2007, geom_horiz_back,\n",
" geom_horiz_forw, K_w_sqr,\n",
" wl_S, wl_C, wl_X, wl_Ku,\n",
" wl_Ka, wl_W)\n",
"from rustmatrix.refractive import m_w_20C\n",
"\n",
"BANDS = [('S', wl_S), ('C', wl_C), ('X', wl_X),\n",
" ('Ku', wl_Ku), ('Ka', wl_Ka), ('W', wl_W)]\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Run the six bands\n"
]
},
{
"cell_type": "code",
"metadata": {},
"execution_count": null,
"outputs": [],
"source": [
"def run(wl):\n",
" s = Scatterer(wavelength=wl, m=m_w_20C[wl],\n",
" Kw_sqr=K_w_sqr[wl], ddelt=1e-4, ndgs=2)\n",
" integ = rs_psd.PSDIntegrator()\n",
" integ.D_max = 8.0\n",
" integ.num_points = 64\n",
" integ.axis_ratio_func = lambda D: 1.0 / dsr_thurai_2007(D)\n",
" integ.geometries = (geom_horiz_back, geom_horiz_forw)\n",
" s.psd_integrator = integ\n",
" s.psd_integrator.init_scatter_table(s)\n",
" s.psd = rs_psd.GammaPSD(D0=1.5, Nw=8e3, mu=4)\n",
" s.set_geometry(geom_horiz_back)\n",
" Zh = 10 * np.log10(radar.refl(s))\n",
" Zdr = 10 * np.log10(radar.Zdr(s))\n",
" s.set_geometry(geom_horiz_forw)\n",
" return Zh, Zdr, radar.Kdp(s), radar.Ai(s)\n",
"\n",
"Zh = {}; Zdr = {}; Kdp = {}; Ai = {}\n",
"for name, wl in BANDS:\n",
" Zh[name], Zdr[name], Kdp[name], Ai[name] = run(wl)\n",
" print(f'{name:<3} Z_h={Zh[name]:6.2f} Z_dr={Zdr[name]:+.3f} '\n",
" f'K_dp={Kdp[name]:+.3f} A_i={Ai[name]:.4f}')\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Plot the frequency response\n"
]
},
{
"cell_type": "code",
"metadata": {},
"execution_count": null,
"outputs": [],
"source": [
"names = [b[0] for b in BANDS]\n",
"wls = np.array([b[1] for b in BANDS])\n",
"fig, axes = plt.subplots(2, 2, figsize=(9, 6))\n",
"axes[0, 0].semilogx(wls, [Zh[n] for n in names], 'o-'); axes[0, 0].set_ylabel('Z_h [dBZ]')\n",
"axes[0, 1].semilogx(wls, [Zdr[n] for n in names], 'o-'); axes[0, 1].set_ylabel('Z_dr [dB]')\n",
"axes[1, 0].semilogx(wls, [Kdp[n] for n in names], 'o-'); axes[1, 0].set_ylabel('K_dp [\u00b0/km]')\n",
"axes[1, 1].semilogx(wls, [Ai[n] for n in names], 'o-'); axes[1, 1].set_ylabel('A_i [dB/km]')\n",
"for ax, label in zip(axes.flat, names*4):\n",
" ax.set_xlabel('wavelength [mm]')\n",
" ax.grid(True, alpha=0.3, which='both')\n",
" for name, wl in BANDS:\n",
" ax.axvline(wl, color='k', alpha=0.1)\n",
"fig.suptitle('Moderate convective rain (D0=1.5 mm): frequency response')\n",
"fig.tight_layout();\n"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"name": "python",
"version": "3.11"
}
},
"nbformat": 4,
"nbformat_minor": 5
}