Dual-Pol and Differential Phase#
Overview#
This module provides algorithms to process polarimetric radar moments, namely the differential phase, \(Phi_{DP}\), and, based on successful \(Phi_{DP}\) retrieval, also the specific differential phase, \(K_{DP}\). Please note that the actual application of polarimetric moments is implemented in the corresponding wradlib modules, e.g.:
fuzzy echo classification from polarimetric moments (
wradlib.classify.classify_echo_fuzzy
)attenuation correction (
wradlib.atten.pia_from_kdp
)direct precipitation retrieval from Kdp (
wradlib.trafo.kdp_to_r
)
Establishing a valid \(Phi_{DP}\) profile for \(K_{DP}\) retrieval
involves despeckling (wradlib.util.despeckle
), phase unfolding, and iterative
retrieval of \(Phi_{DP}\) form \(K_{DP}\).
The main workflow and its single steps is based on a publication by
[Vulpiani et al., 2012]. For convenience, the entire workflow has been
put together in the function wradlib.dp._phidp_vulpiani
.
Once a valid \(Phi_{DP}\) profile has been established, the kdp_from_phidp functions can be used to retrieve \(K_{DP}\).
Please note that so far, the functions in this module were designed to increase performance. This was mainly achieved by allowing the simultaneous application of functions over multiple array dimensions. The only requirement to apply these function is that the range dimension must be the last dimension of all input arrays.
Compute the depolarization ration. |
|
Retrieves \(K_{DP}\) from \(Phi_{DP}\). |
|
Establish consistent \(Phi_{DP}\) profiles from raw data. |
|
Compute the texture of data. |
|
Unfolds differential phase by adjusting values that exceeded maximum ambiguous range. |
|
Alternative phase unfolding which completely relies on \(K_{DP}\). |
|
wradlib xarray SubAccessor methods for DualPol. |
- wradlib.dp.depolarization(zdr, rho)[source]#
- wradlib.dp.depolarization(obj: Dataset, **kwargs)
Compute the depolarization ration.
Compute the depolarization ration using differential reflectivity \(Z_{DR}\) and crosscorrelation coefficient \(Rho_{HV}\) of a radar sweep ([Kilambi et al., 2018], [Melnikov et al., 2013], [Ryzhkov et al., 2017]).
- Parameters
zdr (
float
ornumpy.ndarray
) – differential reflectivityrho (
float
ornumpy.ndarray
) – crosscorrelation coefficient
- Returns
depolarization (
numpy.ndarray
) – array of depolarization ratios with the same shape as input data, numpy broadcasting rules apply
- wradlib.dp.kdp_from_phidp(phidp, *, winlen=7, dr=1.0, method='lanczos_conv', **kwargs)[source]#
- wradlib.dp.kdp_from_phidp(obj: DataArray, *, winlen=7, **kwargs)
Retrieves \(K_{DP}\) from \(Phi_{DP}\).
In normal operation the method uses convolution to estimate \(K_{DP}\) (the derivative of \(Phi_{DP}\)) with Low-noise Lanczos differentiators (method=’lanczos_conv’). The results are very similar to the moving window linear regression (method=’lstsq’), but the former is much faster.
The \(K_{DP}\) retrieval will return NaNs in case at least one value in the moving window is NaN. By default, the remaining gates are treated by using local linear regression where possible.
Please note that the moving window size
winlen
is specified as the number of range gates. Thus, this argument might need adjustment in case the range resolution changes. In the original publication ([Vulpiani et al., 2012]), the valuewinlen=7
was chosen for a range resolution of 1km.Uses
derivate
to calculate the derivation. See for additional kwargs.Warning
The function is designed for speed by allowing to process multiple dimensions in one step. For this purpose, the RANGE dimension needs to be the LAST dimension of the input array.
- Parameters
phidp (
numpy.ndarray
) – multidimensional array, note that the range dimension must be the last dimension of the input array.winlen (
int
, optional) – Width of the window (as number of range gates), defaults to 7dr (
float
, optional) – gate length in km, defaults to 1method (
str
, optional) – Defaults to ‘lanczos_conv’. Can also take one of ‘lanczos_dot’, ‘lstsq’, ‘cov’, ‘cov_nan’, ‘matrix_inv’.
- Keyword Arguments
- Returns
out (
numpy.ndarray
) – array of \(K_{DP}\) with the same shape as phidp
Examples
>>> import wradlib >>> import numpy as np >>> import matplotlib.pyplot as plt >>> plt.interactive(True) >>> kdp_true = np.sin(3 * np.arange(0, 10, 0.1)) >>> phidp_true = np.cumsum(kdp_true) >>> phidp_raw = phidp_true + np.random.uniform(-1, 1, len(phidp_true)) >>> gaps = np.concatenate([range(10, 20), range(30, 40), range(60, 80)]) >>> phidp_raw[gaps] = np.nan >>> kdp_re = wradlib.dp.kdp_from_phidp(phidp_raw) >>> line1 = plt.plot(np.ma.masked_invalid(phidp_true), "b--", label="phidp_true") # noqa >>> line2 = plt.plot(np.ma.masked_invalid(phidp_raw), "b-", label="phidp_raw") # noqa >>> line3 = plt.plot(kdp_true, "g-", label="kdp_true") >>> line4 = plt.plot(np.ma.masked_invalid(kdp_re), "r-", label="kdp_reconstructed") # noqa >>> lgnd = plt.legend(("phidp_true", "phidp_raw", "kdp_true", "kdp_reconstructed")) # noqa >>> plt.show()
- wradlib.dp.phidp_kdp_vulpiani(obj, dr, *, ndespeckle=5, winlen=7, niter=2, copy=False, **kwargs)[source]#
- wradlib.dp.phidp_kdp_vulpiani(obj: DataArray, *, winlen=7, **kwargs)
Establish consistent \(Phi_{DP}\) profiles from raw data.
This approach is based on [Vulpiani et al., 2012] and involves a two-step procedure of \(Phi_{DP}\) reconstruction.
Processing of raw \(Phi_{DP}\) data contains the following steps:
Despeckle
Initial \(K_{DP}\) estimation
Removal of artifacts
Phase unfolding
\(Phi_{DP}\) reconstruction using iterative estimation of \(K_{DP}\)
- Parameters
obj (
numpy.ndarray
) – array of shape (n azimuth angles, n range gates)dr (
float
) – gate length in kmndespeckle (
int
, optional) –ndespeckle
parameter ofdespeckle
, defaults to 5winlen (
int
, optional) –winlen
parameter ofkdp_from_phidp
, defaults to 7niter (
int
, optional) – Number of iterations in which \(Phi_{DP}\) is retrieved from \(K_{DP}\) and vice versa, defaults to 2.copy (
bool
, optional) – if True, the original \(Phi_{DP}\) array will remain unchanged, defaults to False
- Keyword Arguments
- Returns
phidp (
numpy.ndarray
) – array of shape (…, n azimuth angles, n range gates) reconstructed \(Phi_{DP}\)kdp (
numpy.ndarray
) – array of shape (…, n azimuth angles, n range gates)kdp
estimate corresponding tophidp
output
Examples
See Routine verification measures for radar-based precipitation estimates.
- wradlib.dp.texture(data)[source]#
- wradlib.dp.texture(obj: DataArray)
- wradlib.dp.texture(obj: Dataset)
Compute the texture of data.
Compute the texture of the data by comparing values with a 3x3 neighborhood (based on [Gourley et al., 2007]). NaN values in the original array have NaN textures.
- Parameters
data (
numpy.ndarray
) – multidimensional array with shape (…, number of beams, number of range bins)- Returns
texture (
numpy.ndarray
) – array of textures with the same shape as data
- wradlib.dp.unfold_phi(phidp, rho, *, width=5, copy=False)[source]#
- wradlib.dp.unfold_phi(obj: Dataset, **kwargs)
Unfolds differential phase by adjusting values that exceeded maximum ambiguous range.
Accepts arbitrarily dimensioned arrays, but THE LAST DIMENSION MUST BE THE RANGE.
Uses the fast Fortran-based implementation if the speedup module is compiled.
The algorithm is based on the paper of [Wang et al., 2009].
- Parameters
phidp (
numpy.ndarray
) – array of shape (…,nr) with nr being the number of range binsrho (
numpy.ndarray
) – array of same shape asphidp
width (
int
, optional) – Width of the analysis window, defaults to 5.copy (
bool
, optional) – Leaves original phidp array unchanged if set to True (default: False)
- Returns
phidp (
numpy.ndarray
) – array of shape (…, n azimuth angles, n range gates) reconstructed \(Phi_{DP}\)
- wradlib.dp.unfold_phi_vulpiani(phidp, kdp, *, th=-20, winlen=7)[source]#
- wradlib.dp.unfold_phi_vulpiani(obj: Dataset, **kwargs)
Alternative phase unfolding which completely relies on \(K_{DP}\).
This unfolding should be used in oder to iteratively reconstruct \(Phi_{DP}\) and \(K_{DP}\) (see [Vulpiani et al., 2012]).
Note
\(Phi_{DP}\) is assumed to be in the interval [-180, 180] degree. From experience the window for calculation of \(K_{DP}\) should not be too large to catch possible phase wraps.
- Parameters
phidp (
numpy.ndarray
) – array of floatskdp (
numpy.ndarray
) – array of floatsth (
float
, optional) – Threshold th3 in the above citation, defaults to -20.winlen (
int
, optional) – Length of window to fix possible phase over-correction. Normally should take the value of the length of the processing window in the above citation, defaults to 7.
- Returns
phidp (
numpy.ndarray
) – array of floats
- class wradlib.dp.DpMethods(obj)[source]#
Bases:
XarrayMethods
wradlib xarray SubAccessor methods for DualPol.
- depolarization(**kwargs)[source]#
Compute the depolarization ration.
Compute the depolarization ration using differential reflectivity \(Z_{DR}\) and crosscorrelation coefficient \(Rho_{HV}\) of a radar sweep ([Kilambi et al., 2018], [Melnikov et al., 2013], [Ryzhkov et al., 2017]).
Parameter#
obj :
xarray.Dataset
- keyword zdr
name of differential reflectivity
- kwtype zdr
- keyword rho
name crosscorrelation coefficient
- kwtype rho
- returns
depolarization (
xarray.DataArray
) – array of depolarization ratios with the same shape as input data, numpy broadcasting rules apply
- kdp_from_phidp(*, winlen=7, **kwargs)[source]#
Retrieves \(K_{DP}\) from \(Phi_{DP}\).
Parameter#
- obj
xarray.DataArray
DataArray containing differential phase
- keyword winlen
window length
- kwtype winlen
- keyword method
Defaults to ‘lanczos_conv’. Can also take one of ‘lanczos_dot’, ‘lstsq’, ‘cov’, ‘cov_nan’, ‘matrix_inv’.
- kwtype method
- keyword skipna
Defaults to True. Local Linear regression removing NaN values using valid neighbors > min_periods
- kwtype skipna
- keyword min_periods
Minimum number of valid values in moving window for linear regression. Defaults to winlen // 2 + 1.
- kwtype min_periods
- returns
out (
xarray.DataArray
) – DataArray
- obj
- phidp_kdp_vulpiani(*, winlen=7, **kwargs)[source]#
Retrieves \(K_{DP}\) from \(Phi_{DP}\).
Parameter#
- obj
xarray.DataArray
DataArray containing differential phase
- winlenint
window length
- keyword method
Defaults to ‘lanczos_conv’. Can also take one of ‘lanczos_dot’, ‘lstsq’, ‘cov’, ‘cov_nan’, ‘matrix_inv’.
- kwtype method
- keyword skipna
Defaults to True. Local Linear regression removing NaN values using valid neighbors > min_periods
- kwtype skipna
- keyword min_periods
Minimum number of valid values in moving window for linear regression. Defaults to winlen // 2 + 1.
- kwtype min_periods
- returns
phidp (
xarray.DataArray
) – DataArraykdp (
xarray.DataArray
) – DataArray
- obj
- texture()[source]#
Compute the texture of data.
Compute the texture of the data by comparing values with a 3x3 neighborhood (based on [Gourley et al., 2007]). NaN values in the original array have NaN textures.
- Parameters
obj (
xarray.DataArray
) – DataArray- Returns
texture (
xarray.DataArray
) – DataArray
- unfold_phi(**kwargs)[source]#
Unfolds differential phase by adjusting values that exceeded maximum ambiguous range.
Accepts arbitrarily dimensioned arrays, but THE LAST DIMENSION MUST BE THE RANGE.
This is the fast Fortran-based implementation (RECOMMENDED).
The algorithm is based on the paper of [Wang et al., 2009].
- Parameters
obj (
xarray.Dataset
)- Keyword Arguments
- Returns
out (
xarray.DataArray
) – DataArray
- unfold_phi_vulpiani(**kwargs)[source]#
Alternative phase unfolding which completely relies on \(K_{DP}\).
This unfolding should be used in oder to iteratively reconstruct \(Phi_{DP}\) and \(K_{DP}\) (see [Vulpiani et al., 2012]).
Note
\(Phi_{DP}\) is assumed to be in the interval [-180, 180] degree. From experience the window for calculation of \(K_{DP}\) should not be too large to catch possible phase wraps.
- Parameters
obj (
xarray.Dataset
) – Dataset- Keyword Arguments
- Returns
out (
xarray.DataArray
) – DataArray