Main Content

# weatherTimeSeries

Simulate I/Q signals for weather returns using a Monte Carlo approach

Since R2024a

## Description

This function generates I/Q signals for monostatic polarimetric weather radar systems. Each weatherTimeSeries simulation produces an independent radar return from a particular resolution volume, or cell, that is modeled as a complex stationary Gaussian random process using Monte Carlo method (see Algorithms).

[Vh,Vv] = weatherTimeSeries(freq,Pr,Vr,SW,ZDR,Rho,Phi) returns I/Q signals as complex voltages for horizontal and vertical polarizations for a monostatic polarimetric weather radar system. Resolution cells are characterized by polarimetric variables specified by Pr, Vr, SW, ZDR, Rho, and Phi.

example

[Vh,Vv] = weatherTimeSeries(freq,Pr,Vr,SW,ZDR,Rho,Phi,Name=Value) returns I/Q signals for a monostatic polarimetric weather radar with additional options specified using one or more name-value arguments. For example, NumPulses sets the number of pulses and NumTrials controls the number of trials.

## Examples

collapse all

This example generates weatherTimeSeries simulations and displays the radar estimates alongside the truth values.

Simulate Weather Time Series I/Q Data

Simulate 100 trials of time series I/Q signals for a resolution volume collected using a polarimetric weather radar operating at a frequency of 3 GHz with a PRF of 1000 Hz and 64 pulses in the dwell time. Assume the resolution volume is comprised of pure rain which has a signal power of -80 dBW, a radial velocity of 10 m/s, a spectrum width of 4 m/s, a differential reflectivity of 2 dB, a correlation coefficient of 0.99, and a differential phase of 60 degrees.

% Define radar and weather parameters
freq = 3e9;                     % Operating frequency (Hz)
Pr = -80;                       % Signal power (dBW)
Vr = 10;                        % Radial velocity (m/s)
SW = 4;                         % Spectrum width (m/s)
ZDR = 2;                        % Differential reflectivity (dB)
Rho = 0.99;                     % Correlation coefficient
Phi = 60;                       % Differential phase (degree)
prf = 1000;                     % PRF (Hz)
numPulses = 64;                 % Number of pulses per trial
numTrials = 100;                % Number of trials

% Generate time series
[Vh,Vv] = weatherTimeSeries(freq,Pr,Vr,SW,ZDR,Rho,Phi,...
PRF=prf,NumPulses=numPulses,NumTrials=numTrials);

Verify Simulation Results

Verify that the radar estimates derived from the simulated I/Q are consistent with the truth.

% Calculate covariance matrix
R0_hh = mean(abs(Vh.^2));
R0_vv = mean(abs(Vv.^2));
C0_hv = mean(conj(Vh).*Vv);
R1_hh = mean(conj(Vh(1:numPulses-1,:)).*(Vh(2:numPulses,:)));

% Pulse Pair Processing
lambda = freq2wavelen(freq);
Pr1 = pow2db(R0_hh);
Vr1 = -lambda*prf/(4*pi)*angle(R1_hh);
SW1 = abs(lambda*prf/(4*pi).*(2*log(abs(R0_hh)./abs(R1_hh))).^0.5);
ZDR1 = pow2db(R0_hh./R0_vv);
Rho1 = abs(C0_hv)./(R0_hh).^0.5./(R0_vv).^0.5;
Phi1 = mod(angle(C0_hv)*180/pi,360);

% Display estimation result and compare with truth
VariableName = {'Pr';'Vr';'SW';'ZDR';'Rho';'Phi'};
Truth = [Pr;Vr;SW;ZDR;Rho;Phi];
Estimate = [round(mean(Pr1),2);round(mean(Vr1),2);round(mean(SW1),2);...
round(mean(ZDR1),2);round(mean(Rho1),3);round(mean(Phi1),2)];
Unit = {'dBW';'m/s';'m/s';'dB';'dimensionless';'degree'};
T = table(VariableName,Truth,Estimate,Unit)
T=6×4 table
VariableName    Truth    Estimate          Unit
____________    _____    ________    _________________

{'Pr' }        -80      -80.1      {'dBW'          }
{'Vr' }         10       9.85      {'m/s'          }
{'SW' }          4       3.96      {'m/s'          }
{'ZDR'}          2       1.99      {'dB'           }
{'Rho'}       0.99      0.991      {'dimensionless'}
{'Phi'}         60      59.97      {'degree'       }

## Input Arguments

collapse all

Radar operating frequency, specified as a positive scalar. Units are in hertz (Hz).

Data Types: double

Signal power, specified as a scalar or length-N vector where N is the number of resolution cells. Units are in decibel-watt (dBW). Pr is related to the radar constant C, target range R, and reflectivity factor Z by Pr = Z -20log10(R) + C, where C is given in dB, R is in km, and Z is in dBZ [1].

Data Types: double

Radial velocity of the weather returns relative to the radar system, specified as a scalar or length-N vector where N is the number of resolution cells. Units are in meters per second (m/s). The absolute value of Vr should not exceed the maximum unambiguous velocity of the radar. The maximum unambiguous velocity is equal to the radar wavelength (in m) times the pulse repetition frequency (in Hz) divided by 4.

Data Types: double

Spectrum width, a measure of the velocity dispersion, specified as a scalar or length-N vector, where N is the number of resolution cells. Units are in meters per second (m/s).

Data Types: double

Differential reflectivity between horizontally and vertically polarized signals, specified as a scalar or length-N vector, where N is the number of resolution cells. Units are in decibel (dB).

Data Types: double

Correlation coefficient between horizontally and vertically polarized signals specified as a scalar or length-N vector, where N is the number of resolution cells. The value of Rho, which is commonly denoted as ρhv, must be greater than 0 and less than or equal to 1. Units are dimensionless.

Data Types: double

Differential phase between horizontally and vertically polarized signals, specified as a scalar or length-N vector, where N is the number of resolution cells. The value of Phi, which is commonly denoted as ${\varphi }_{DP}$, must be greater than or equal to 0 and less than 360. Units are in degrees (deg).

Data Types: double

### Name-Value Arguments

Specify optional pairs of arguments as Name1=Value1,...,NameN=ValueN, where Name is the argument name and Value is the corresponding value. Name-value arguments must appear after other arguments, but the order of the pairs does not matter.

Example: PRF=1000,NumPulses=64,NumTrials=100

Pulse repetition frequency (PRF), specified as a positive scalar. Units are in hertz (Hz).

Data Types: double

Number of pulses per trial, specified as a positive integer scalar. Units are dimensionless.

Data Types: double

Number of trials, specified as a positive integer scalar. Units are dimensionless.

Data Types: double

## Output Arguments

collapse all

Complex voltage for horizontal polarization returned as a complex-valued N-by-M-by-L array where N is the number of resolution cells, M is the number of pulses (NumPulses), and L is the number of trials (NumTrials). By default, the number of pulses is set to 64, and the number of trials is set to 1. The resolution cells are characterized by polarimetric variables as input arguments Pr, Vr, SW, ZDR, Rho, and Phi. Each trial produces an independent radar return from a particular resolution volume, or cell, that is modeled as a complex stationary Gaussian random process using Monte Carlo method.

Data Types: double

Complex voltage for vertical polarization returned as a complex-valued N-by-M-by-L array where N is the number of resolution cells, M is the number of pulses (NumPulses), and L is the number of trials (NumTrials). By default, the number of pulses is set to 64, and the number of trials is set to 1. The resolution cells are characterized by polarimetric variables as input arguments Pr, Vr, SW, ZDR, Rho, and Phi. Each trial produces an independent radar return from a particular resolution volume, or cell, that is modeled as a complex stationary Gaussian random process using Monte Carlo method.

Data Types: double

## Algorithms

weatherTimeSeries simulates I/Q signals using a numerical Monte Carlo approach combined with a statistical model that is based on the expected behavior of radar returns from weather phenomena [2] [3]. weatherTimeSeries calculates scattering amplitudes for every time step using relevant scattering parameters specified as input arguments to generate a time series.

Statistical model assumptions for radar returns from weather phenomena:

• I/Q signals follow Gaussian distribution

• Signal amplitude follows Rayleigh distribution

• Signal phase follows uniform distribution

The first step in the Monte Carlo approach is to specify scattering amplitudes and generate random numbers for particle position and motion. Scattering amplitudes for a given resolution volume are modeled as combinations from random multiple scattering centers. Assume that scatterers within the resolution volume have the same size and a canting angle of zero (which implies that off-diagonal elements of the backscattering matrices have a value of zero). Then we have

${P}_{r}=N{|{V}_{hh}|}^{2},$

where Pr is the received power for horizontal polarization, N is the number of scatterers in the resolution volume, and Vhh is the complex voltage for an individual scatterer for horizontal polarization. At least 20 scatterers (N20) are necessary to adequately model Gaussian random signals.

The amplitude of Vhh can be obtained as

$|{V}_{hh}|=\sqrt{{P}_{r}/N}\text{\hspace{0.17em}}.$

Accordingly, the amplitude of the complex voltage for an individual scatterer for vertical polarization can be calculated as

$|{V}_{vv}|=|{V}_{hh}|\text{\hspace{0.17em}}/\sqrt{{Z}_{dr}}\text{\hspace{0.17em}},$

where Zdr is the differential reflectivity on a linear scale.

In addition, the correlation coefficient ρhv is assumed to be reduced by a factor of ${\text{e}}^{{}^{-{\sigma }_{\delta }^{2}/2}}$ due to the random scattering phase difference, where σδ is the standard deviation of random scattering phase difference, that is, ${\rho }_{hv}={\text{e}}^{{}^{-{\sigma }_{\delta }^{2}/2}}$. Then we have

${\sigma }_{\delta }=\sqrt{-2\mathrm{ln}\left({\rho }_{hv}\right)}=\sqrt{{\sigma }^{2}{}_{\delta h}+{\sigma }^{2}{}_{\delta v}}\text{\hspace{0.17em}}.$

For simplicity, it is assumed that ${\sigma }_{\delta h}={\sigma }_{\delta v}={\sigma }_{\delta }\text{\hspace{0.17em}}/\sqrt{2}$, where σδh and σδv are standard deviations of the backscattering phase δh for horizontal polarization and the backscattering phase δv for vertical polarization, respectively.

Next, the scattered wave field for each particle is calculated. The amplitudes of complex voltage for an individual scatterer can be expressed as

${V}_{hh}=|{V}_{hh}|\cdot {\text{e}}^{-j{\delta }_{h}}\cdot {\text{e}}^{{}^{-j{\varphi }_{DP}}}$

${V}_{vv}=|{V}_{vv}|\cdot {\text{e}}^{-j{\delta }_{v}},$

where ${\varphi }_{DP}$ is the differential phase.

The final step is to calculate the total scattered wave fields. The total complex voltage of the resolution volume is

${V}_{h}={\sum }_{l=1}^{N}{V}_{hh}\cdot {\text{e}}^{-j2{k}_{i}\cdot {r}_{l}}$

${V}_{v}={\sum }_{l=1}^{N}{V}_{vv}\cdot {\text{e}}^{-j2{k}_{i}\cdot {r}_{l}},$

where ${k}_{i}$ is the incident wave vector and ${r}_{l}$ is the random position of each scatterer.

## References

[1] Doviak, R. and D. S. Zrnic. Doppler Radar and Weather Observations, 2nd Ed. New York: Dover Publications, 2006.

[2] Zhang, G. Weather Radar Polarimetry. Boca Raton: CRC Press, 2016.

[3] Li, Z, et al. Polarimetric phased array weather radar data quality evaluation through combined analysis, simulation, and measurements. IEEE Geosci. Remote Sens. Lett., vol. 18, no. 6, pp. 1029–1033, Jun. 2021.

## Version History

Introduced in R2024a