wradlib.dp.kdp_from_phidp#
- wradlib.dp.kdp_from_phidp(phidp, *, winlen=7, dr=1.0, method='lanczos_conv', skipna=True, **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
) – multi-dimensional array, note that the range dimension must be the last dimension of the input array.- Keyword Arguments
winlen (
int
) – Width of the window (as number of range gates)dr (
float
) – gate length in kmmethod (
str
) – Defaults to ‘lanczos_conv’. Can also take one of ‘lanczos_dot’, ‘lstsq’, ‘cov’, ‘cov_nan’, ‘matrix_inv’.skipna (
bool
) – Defaults to True. Local Linear regression removing NaN values using valid neighbors > min_periodsmin_periods (
int
) – Minimum number of valid values in moving window for linear regression. Defaults to winlen // 2 + 1.
- 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 pl >>> pl.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 = pl.plot(np.ma.masked_invalid(phidp_true), "b--", label="phidp_true") # noqa >>> line2 = pl.plot(np.ma.masked_invalid(phidp_raw), "b-", label="phidp_raw") # noqa >>> line3 = pl.plot(kdp_true, "g-", label="kdp_true") >>> line4 = pl.plot(np.ma.masked_invalid(kdp_re), "r-", label="kdp_reconstructed") # noqa >>> lgnd = pl.legend(("phidp_true", "phidp_raw", "kdp_true", "kdp_reconstructed")) # noqa >>> pl.show()