rustmatrix.spectra.beam#

Radar beam patterns and scene-integrated spectra.

When a radar beam has non-zero width the observed Doppler spectrum is a solid-angle weighted sum of per-pixel spectra. For a homogeneous scene across the beam, turbulence and fall speeds independent of position, and horizontal wind in the beam this reduces to the closed-form

σ_beam = |u_h| · θ_b / (2 √(2 ln 2))

that SpectralIntegrator already builds in through its u_h / beamwidth arguments. For a heterogeneous scene — isolated convective cells, updraft/downdraft couplets, horizontal reflectivity gradients — the analytic formula fails: the beam samples regions with different reflectivities, vertical motions, and wind speeds, and the observed spectrum is a genuine spatial integral.

This module provides the pieces needed for that integral:

  • BeamPattern base and three concrete subclasses:

    • GaussianBeam — canonical well-tapered main-lobe; no sidelobes.

    • AiryBeam — uniform circular aperture pattern [2 J₁(x)/x]² with realistic sidelobes (first at −17.6 dB).

    • TabulatedBeam — user-supplied (θ, gain) samples, for measured patterns or Taylor tapers.

  • Scene — spatial fields Z_dBZ(x,y,z), w(x,y,z), u_h(x,y,z) that the integrator evaluates at each beam sample.

  • BeamIntegrator — drives a scatterer + beam + scene + PSD factory through a Doppler-velocity grid and returns a SpectralResult.

Coordinate conventions#

Radar at the origin, boresight along -ẑ (down-looking) is the primary use case. A beam sample offset by (θ, φ) from boresight points along (sin θ cos φ, sin θ sin φ, -cos θ). The line-of-sight Doppler of a particle with velocity (u_x, u_y, -v_fall) — where v_fall = v_t(D) + w is the downward fall rate — is

v_obs = (v_t + w) cos θ + u_h sin θ cos(φ − φ_wind)

which collapses to v_t + w at boresight (θ = 0). Positive v_obs means the particle is moving away from a down-looking radar, i.e. in the fall direction — same convention as the rest of rustmatrix.spectra.

References

Doviak, R. J. & Zrnić, D. S. (1993). Doppler Radar and Weather Observations, 2nd ed., Academic Press. §3.2, §5.3.

Balanis, C. A. (2016). Antenna Theory: Analysis and Design, 4th ed., Wiley. §15.4 (circular aperture Airy pattern).

class rustmatrix.spectra.beam.BeamPattern[source]#

Bases: object

Base class for circularly symmetric one-way power patterns.

Subclasses define gain(), the normalized one-way power pattern G(θ) with G(0) = 1. sample() returns a set of (θ, φ, weight) triples suitable for driving a BeamIntegrator; default weights are proportional to G²(θ) sin θ (two-way pattern × solid-angle element) and sum to 1 over the sampled cone.

hpbw#

One-way half-power full-width [rad].

Type:

float

sample(n_theta=32, n_phi=24, max_theta=None)[source]#

Return (theta, phi, weight) flat arrays on the beam cone.

Parameters:

  • n_theta (int) – Number of polar / azimuthal samples. The default 32 × 24 is generally enough to recover the analytic σ_beam to ~1 %.

  • n_phi (int) – Number of polar / azimuthal samples. The default 32 × 24 is generally enough to recover the analytic σ_beam to ~1 %.

  • max_theta (float, optional) – Outer cone half-angle [rad]. Default 3 × hpbw — captures >99.9 % of a Gaussian two-way pattern and a few Airy sidelobes.

Return type:

Tuple[ndarray, ndarray, ndarray]

Notes

Weights are two-way, G²(θ) sin θ, normalized so they sum to 1. The line-of-sight velocity projection is the only thing that knows about φ, so for a circularly symmetric beam the azimuth integration is exactly Simpson’s rule on a periodic grid.

class rustmatrix.spectra.beam.GaussianBeam(hpbw)[source]#

Bases: BeamPattern

Gaussian one-way power pattern.

G(θ) = exp(-θ² / (2 σ²)) with σ = hpbw / (2 √(2 ln 2)). No sidelobes; well-tapered aperture illumination. This is the pattern implicit in SpectralIntegrator when beamwidth > 0 is supplied.

Parameters:

hpbw (float)

class rustmatrix.spectra.beam.AiryBeam(hpbw)[source]#

Bases: BeamPattern

Uniformly illuminated circular aperture (Airy) pattern.

G(θ) = [2 J₁(x) / x]² with x = α sin θ and α chosen so the pattern crosses half-power at θ = hpbw / 2. Gives the canonical parabolic-dish response — a main lobe slightly narrower than a matched-HPBW Gaussian, followed by the first sidelobe at −17.6 dB and nulls between lobes.

Parameters:

hpbw (float) – One-way half-power beamwidth [rad].

X_HALFPOWER = 1.6163399#

Argument x where [2 J₁(x)/x]² = 0.5.

class rustmatrix.spectra.beam.TabulatedBeam(theta, gain, hpbw=None)[source]#

Bases: BeamPattern

Beam pattern interpolated from a user-supplied table.

Parameters:

  • theta (array_like) – Polar-angle samples [rad], strictly increasing and starting at 0.

  • gain (array_like) – One-way power pattern at those angles, gain[0] = 1 assumed.

  • hpbw (float, optional) – Advertised one-way half-power beamwidth [rad]. If omitted it is inferred from the table by linear interpolation at G = 0.5.

class rustmatrix.spectra.beam.Scene(Z_dBZ, w, u_h, u_h_azimuth=0.0)[source]#

Bases: object

Spatial fields evaluated at beam-sample positions.

Each field is a callable f(x, y, z) -> array that accepts numpy arrays of shape (N,) (positions of a batch of beam samples) and returns an (N,) array.

Parameters:
Z_dBZ#

Equivalent reflectivity field [dBZ].

Type:

callable

w#

Mean vertical air motion [m/s], positive-downward.

Type:

callable

u_h#

Horizontal wind speed magnitude [m/s].

Type:

callable

u_h_azimuth#

Horizontal wind azimuth [rad], measured from +x axis. If a scalar, treated as constant; if callable, evaluated like the other fields. Default 0 (wind blows in +x direction).

Type:

float or callable

evaluate(x, y, z)[source]#

Return (Z_dBZ, w, u_h, u_h_azimuth) at each (x, y, z).

All inputs and outputs are 1-D arrays of the same length.

class rustmatrix.spectra.beam.BeamIntegrator(scatterer, beam, scene, psd_factory, fall_speed, turbulence=None, *, radar_position=(0.0, 0.0, 0.0), boresight=(0.0, 0.0, -1.0), range_m=1.0, v_bins=None, v_min=None, v_max=None, n_bins=None, n_theta=32, n_phi=24, max_theta_over_hpbw=3.0, geometry_backscatter=(0.0, 180.0, 0.0, 0.0, 0.0, 0.0), geometry_forward=None)[source]#

Bases: object

Integrate a spectral S/Z response over a beam pattern + scene.

Parameters:

  • scatterer (Scatterer) – Single-species scatterer with a PSDIntegrator already populated for the backscatter (and, optionally, forward) geometry. The per-diameter _S_table / _Z_table are reused at every beam sample — only the PSD weights change.

  • beam (BeamPattern) – Beam pattern whose sample() method produces the angular quadrature grid.

  • scene (Scene) – Scene fields evaluated at each beam-sample pixel.

  • psd_factory (callable) – Z_dBZ -> PSD mapping used to build the per-pixel PSD.

  • fall_speed (callable) – D_mm -> v_t_m_per_s, same contract as SpectralIntegrator.

  • turbulence (callable or _TurbulenceModel, optional) – Per-pixel turbulence. Default zero turbulence.

  • radar_position (3-tuple) – (x, y, z) of the radar [m]. Default (0, 0, 0).

  • boresight (3-tuple) – Unit vector along the beam axis. Default (0, 0, -1) (down-looking).

  • range_m (float) – Slant range to the range-gate centre [m].

  • v_bins (array_like, optional) – Explicit 1-D velocity grid [m/s]. Mutually exclusive with (v_min, v_max, n_bins).

  • v_min (float, float, int, optional) – Convenience triple.

  • v_max (float, float, int, optional) – Convenience triple.

  • n_bins (float, float, int, optional) – Convenience triple.

  • n_theta (int) – Beam-sample counts. Defaults 32 × 24.

  • n_phi (int) – Beam-sample counts. Defaults 32 × 24.

  • max_theta_over_hpbw (float) – Outer sampling cone in multiples of hpbw. Default 3.

  • geometry_backscatter (6-tuples) – Same meaning as SpectralIntegrator.

  • geometry_forward (6-tuples) – Same meaning as SpectralIntegrator.

Notes

This integrator is single-species. For a mixture, run one integrator per species and add the resulting S_spec / Z_spec on the shared velocity grid (the same linearity argument as HydroMix).

run()[source]#

Evaluate the beam-integrated spectrum.

Returns:

Same container used by SpectralIntegrator. All bulk-sum identities still hold (the spectrum is linear in the beam weights and per-pixel PSD).

Return type:

SpectralResult

rustmatrix.spectra.beam.marshall_palmer_psd_factory(N0=8000.0, D_max=6.0, Z_floor_dBZ=-30.0)[source]#

Return a PSD factory that maps dBZ to a Marshall–Palmer exponential PSD.

Uses the Rayleigh-equivalent Z-Λ relation for an exponential PSD, Z ≈ 720 · N₀ / Λ⁷ (mm⁶ m⁻³ with D in mm). This is a self-consistent analytic inversion — accurate at rain-band wavelengths where Rayleigh scattering dominates — and avoids a root-finding step at each pixel. Below Z_floor_dBZ the factory returns a PSD with effectively zero concentration (Λ very large), so empty pixels contribute nothing.

Parameters:

  • N0 (float) – Intercept parameter [mm⁻¹ m⁻³]. Default 8000 (Marshall–Palmer).

  • D_max (float) – Upper integration cutoff [mm].

  • Z_floor_dBZ (float) – Reflectivities below this are treated as “empty”.

Returns:

factory – factory(Z_dBZ) -> ExponentialPSD suitable as the BeamIntegrator PSD mapper.

Return type:

callable