Interpolation#
Interpolation allows to transfer data from one set of locations to another. This includes for example:
interpolating the data from a polar grid to a cartesian grid or irregular points
interpolating point observations to a grid or a set of irregular points
filling missing values, e.g. filling clutters
The base class for interpolation in N dimensions. |
|
Nearest-neighbour interpolation in N dimensions. |
|
Inverse distance weighting interpolation in N dimensions. |
|
Interface to the |
|
Interpolate using Ordinary Kriging |
|
ExternalDriftKriging(src, trg, cov='1.0 Exp(10000.)', nnearest=12, |
|
Interpolation on a 2d grid in arbitrary dimensions. |
|
Bin points values to regular grid cells |
|
Map values representing quadrilateral grid cells to another quadrilateral grid. |
|
Convenience function to use the interpolation classes in an efficient way |
|
Convenience function to interpolate polar data |
|
Interpolate array |
|
Map array |
|
wradlib xarray SubAccessor methods for Ipol Methods. |
- class wradlib.ipol.IpolBase(src, trg)[source]#
Bases:
object
The base class for interpolation in N dimensions. Provides the basic interface for all other classes.
- Parameters
src (
numpy.ndarray
) – ndarray of floats, shape (npoints, ndims) Data point coordinates of the source points.trg (
numpy.ndarray
) – ndarray of floats, shape (npoints, ndims) Data point coordinates of the target points.
- class wradlib.ipol.Nearest(src, trg)[source]#
Bases:
IpolBase
Nearest-neighbour interpolation in N dimensions.
- Parameters
src (
numpy.ndarray
) – ndarray of floats, shape (npoints, ndims) or cKDTree object Data point coordinates of the source points.trg (
numpy.ndarray
) – ndarray of floats, shape (npoints, ndims) Data point coordinates of the target points.remove_missing (
int
) – Number of neighbours to consider in the presence of NaN, defaults to 0.
- Keyword Arguments
**kwargs (
dict
) – keyword arguments of ipclass (see class documentation)
Examples
See How to use wradlib's ipol module for interpolation tasks?.
Note
- class wradlib.ipol.Idw(src, trg, nnearest=4, p=2.)[source]#
Bases:
IpolBase
Inverse distance weighting interpolation in N dimensions.
- Parameters
src (
numpy.ndarray
) – ndarray of floats, shape (npoints, ndims) of cKDTree object Data point coordinates of the source points.trg (
numpy.ndarray
) – ndarray of floats, shape (npoints, ndims) Data point coordinates of the target points.nnearest (
int
) – max. number of neighbours to be consideredp (
float
) – inverse distance power used in 1/dist**premove_missing (
bool
) – If True masks NaN values in the data values, defaults to False
- Keyword Arguments
**kwargs (
dict
) – keyword arguments of ipclass (see class documentation)
Examples
See How to use wradlib's ipol module for interpolation tasks?.
Note
- class wradlib.ipol.Linear(src, trg, *, remove_missing=False)[source]#
Bases:
IpolBase
Interface to the
scipy.interpolate.LinearNDInterpolator
class.We provide this class in order to achieve a uniform interface for all Interpolator classes
- Parameters
src (
numpy.ndarray
) – ndarray of floats, shape (npoints, ndims) Data point coordinates of the source points.trg (
numpy.ndarray
) – ndarray of floats, shape (npoints, ndims) Data point coordinates of the target points.
Examples
See How to use wradlib's ipol module for interpolation tasks?.
- class wradlib.ipol.OrdinaryKriging(src, trg, cov='1.0 Exp(10000.)', nnearest=12)[source]#
Bases:
IpolBase
Interpolate using Ordinary Kriging
(Co-)Variogram definitions are given in the syntax that
gstat
uses. It allows nesting of different basic variogram types using linear combinations. Each basic covariogram is usually defined by Note that, strictly speaking, this implementation doesn’t allow Kriging of fields for which the covariance does not exist. While this is mathematically possible, it is rather rare for fields encountered in reality. Therefore, this should not be a severe limitation.Most (co-)variograms are characterized by a sill parameter (which is the (co-)variance at separation distance 0) a range parameter (which indicates a separation distance after which the covariance drops close to zero) a sometimes additional parameters governing the shape of the function. In the following range is given by the variable r and the sill by the variable s. Currently implemented are:
Pure Nugget
Exponential
Spherical
Gaussian
Linear
Matern
Power
Cauchy
- Parameters
src (
numpy.ndarray
) – ndarray of floats, shape (npoints, ndims) Data point coordinates of the source points.trg (
numpy.ndarray
) – ndarray of floats, shape (npoints, ndims) Data point coordinates of the target points.cov (
str
) – covariance (variogram) model string in the syntaxgstat
uses.nnearest (
int
) – max. number of neighbours to be consideredremove_missing (
bool
) – If True masks NaN values in the data values, defaults to False
Note
The class calculates the Kriging weights during initialization, because these only depend on the configuration of the points.
The call method is then only used to calculate estimated values at the target points based on those at the source points. Therefore, the main computational load is experienced during initialization. This behavior is different from that of the Idw or Nearest Interpolators.
After initialization the estimation variance at each interpolation target may be retrieved from the attribute estimation_variance.
Examples
See How to use wradlib's ipol module for interpolation tasks?.
- class wradlib.ipol.ExternalDriftKriging(src, trg, cov='1.0 Exp(10000.)', *, nnearest=12, src_drift=None, trg_drift=None, remove_missing=False, **kwargs)[source]#
Bases:
IpolBase
- ExternalDriftKriging(src, trg, cov=’1.0 Exp(10000.)’, nnearest=12,
drift_src=None, drift_trg=None)
Kriging with external drift
- Parameters
src (
numpy.ndarray
) – ndarray of floats, shape (nsrcpoints, ndims) Data point coordinates of the source points.trg (
numpy.ndarray
) – ndarray of floats, shape (ntrgpoints, ndims) Data point coordinates of the target points.cov (
str
) – covariance (variogram) model string in the syntaxgstat
uses.nnearest (
int
) – max. number of neighbours to be consideredsrc_drift (
numpy.ndarray
) – ndarray of floats, shape (nsrcpoints, ) values of the external drift at each source pointtrg_drift (
numpy.ndarray
) – ndarray of floats, shape (ntrgpoints, ) values of the external drift at each target point
See also
Note
After calling the object in order to get the interpolated values, the estimation variance of the system may be retrieved from the attribute estimation_variance. Accordingly, the interpolation weights can be retrieved from the attribute weights
If drift_src or drift_trg are not given on initialization, they must be provided when using the __call__ method. If any of them is given on initialization its values may be overridden by passing new values to the __call__ method.
Examples
See How to use wradlib's ipol module for interpolation tasks?.
- class wradlib.ipol.RectGrid(src, trg, *, method='linear')[source]#
Bases:
RectGridBase
Interpolation on a 2d grid in arbitrary dimensions.
The source data must be defined on a regular grid, the grid spacing however may be uneven. Linear, nearest-neighbour and spline interpolation are supported.
Based on
scipy.interpolate.interpn
, uses: - nearestscipy.interpolate.RegularGridInterpolator
- linearscipy.interpolate.RegularGridInterpolator
- splinef2dscipy.interpolate.RectBivariateSpline
- Parameters
src (
numpy.ndarray
) – 3d array of shape (…, 2) The points defining the regular grid in n dimensions.trg (
numpy.ndarray
) – Array of shape (…, ndim) The coordinates to sample the gridded data at
- Keyword Arguments
method (
str
) – Method of interpolation used, defaults to ‘linear’.
- class wradlib.ipol.RectBin(src, trg)[source]#
Bases:
RectGridBase
Bin points values to regular grid cells
Based on
scipy.stats.binned_statistic_dd
- Parameters
src (
numpy.ndarray
) – data point coordinates of the source points with shape (…, ndims).trg (
numpy.ndarray
) – rectangular grid coordinates (center) with shape (…, 2)
- property binned_stats#
- class wradlib.ipol.QuadriArea(src, trg, **kwargs)[source]#
Bases:
PolyArea
Map values representing quadrilateral grid cells to another quadrilateral grid.
Based on
wradlib.zonalstats
- Parameters
src (
numpy.ndarray
) – Source grid edge coordinates with shape (n+1, m+1, 2).trg (
numpy.ndarray
) – Target grid edge coordinates with shape (o+1, p+1, 2).
- wradlib.ipol.interpolate(src, trg, vals, ipclass, *args, **kwargs)[source]#
Convenience function to use the interpolation classes in an efficient way
The interpolation classes in
wradlib.ipol
are computationally very efficient if they are applied on large multidimensional arrays of which the first dimension must be the locations’ dimension (1d or 2d coordinates) and the following dimensions can be anything (e.g. time or ensembles). This way, the weights need to be computed only once. However, this can only be done with success if all source values for the interpolation are valid numbers. If the source values contain let’s say np.nan types, the result of the interpolation will be np.nan in the vicinity of the corresponding points, too. Imagine that you have a time series of observations at points and in each time step one observation is missing. You would still like to efficiently apply the interpolation classes, but you will need to account for the resulting np.nan values in the interpolation output.In order to still allow for the efficient application, you have to take care of the remaining np.nan in your interpolation result. This is done by this convenience function.
Alternatively, you have to make sure that your vals argument does not contain any np.nan values OR you have to post-process missing values in your interpolation result in another way.
Warning
Works only for one- and two-dimensional vals arrays, yet.
- Parameters
src (
numpy.ndarray
) – ndarray of floats, shape (npoints, ndims) Data point coordinates of the source points.trg (
numpy.ndarray
) – ndarray of floats, shape (npoints, ndims) Data point coordinates of the target points.vals (
numpy.ndarray
) – ndarray of float, shape (numsourcepoints, …) Values at the source points which to interpolateipclass (
wradlib.ipol.IpolBase
) – A class which inherits from IpolBase.
- Other Parameters
*args (
list
) – arguments of ipclass (see class documentation)- Keyword Arguments
**kwargs (
dict
) – keyword arguments of ipclass (see class documentation)
Examples
>>> # test for 1 dimension in space and two value dimensions >>> src = np.arange(10)[:,None] >>> trg = np.linspace(0,20,40)[:,None] >>> vals = np.hstack((np.sin(src), 10.+np.sin(src))) >>> # here we introduce missing values only in the second dimension >>> vals[3:5,1] = np.nan >>> ipol_result = interpolate(src, trg, vals, Idw, nnearest=2) >>> import matplotlib.pyplot as plt >>> plt.interactive(True) >>> line1 = plt.plot(trg, ipol_result, 'b+') >>> line2 = plt.plot(src, vals, 'ro')
- wradlib.ipol.interpolate_polar(data, *, mask=None, ipclass=<class 'wradlib.ipol.Nearest'>)[source]#
- wradlib.ipol.interpolate_polar(obj: DataArray, mask, **kwargs)
Convenience function to interpolate polar data
- Parameters
data (
numpy.ndarray
) – 2-dimensional array (azimuth, ranges) of floats;if no mask is assigned explicitly polar data should be a masked array
mask (
numpy.ndarray
, optional) – boolean array with pixels to be interpolated set to True, must have the same shape as data, defaults to None.ipclass (
wradlib.ipol.IpolBase
) – A class which inherits from IpolBase, defaults towradlib.ipol.Nearest
.
- Returns
filled_data (
numpy.ndarray
) – 2D array with interpolated values for the values set to True in the mask
Examples
>>> import numpy as np # noqa >>> import wradlib as wrl >>> # creating a data array and mask some values >>> data = np.arange(12.).reshape(4,3) >>> masked_values = (data==2) | (data==9) >>> # interpolate the masked data based on ''masked_values'' >>> filled_a = wrl.ipol.interpolate_polar(data, mask = masked_values, ipclass = wrl.ipol.Linear) # noqa >>> da = wrl.georef.create_xarray_dataarray(filled_a) >>> da = da.wrl.georef.georeference() >>> pm = wrl.vis.plot(da) >>> # the same result can be achieved by using a masked array instead of an explicit mask # noqa >>> mdata = np.ma.array(data, mask = masked_values) >>> filled_b = wrl.ipol.interpolate_polar(mdata, ipclass = wrl.ipol.Linear) # noqa >>> da = wrl.georef.create_xarray_dataarray(filled_b) >>> da = da.wrl.georef.georeference() >>> pm = wrl.vis.plot(da)
- wradlib.ipol.cart_to_irregular_interp(cartgrid, values, newgrid, **kwargs)[source]#
Interpolate array
values
defined by cartesian coordinate arraycartgrid
to new coordinates defined bynewgrid
using nearest neighbour, linear or cubic interpolationSlow for large arrays
Keyword arguments are fed to
scipy.interpolate.griddata
- Parameters
cartgrid (
numpy.ndarray
) – 3-dimensional array (nx, ny, lon/lat) of floats;values (
numpy.ndarray
) – 2-dimensional array (nx, ny) of data valuesnewgrid (
numpy.ndarray
) – Nx2-dimensional array (…, lon/lat) of floatskwargs (
scipy.interpolate.griddata
)
- Returns
interp (
numpy.ndarray
) – array with interpolated values of size N
- wradlib.ipol.cart_to_irregular_spline(cartgrid, values, newgrid, **kwargs)[source]#
Map array
values
defined by cartesian coordinate arraycartgrid
to new coordinates defined bynewgrid
using spline interpolation.Keyword arguments are fed through to
scipy.ndimage.map_coordinates
- Parameters
cartgrid (
numpy.ndarray
) – 3-dimensional array (nx, ny, lon/lat) of floatsvalues (
numpy.ndarray
) – 2-dimensional array (nx, ny) of data valuesnewgrid (
numpy.ndarray
) – Nx2-dimensional array (…, lon/lat) of floatskwargs (
scipy.ndimage.map_coordinates
)
- Returns
interp (
numpy.ndarray
) – array with interpolated values of size N
Examples
- class wradlib.ipol.IpolMethods(obj)[source]#
Bases:
XarrayMethods
wradlib xarray SubAccessor methods for Ipol Methods.
- interpolate_polar(*, mask=None, ipclass=<class 'wradlib.ipol.Nearest'>)[source]#
Convenience function to interpolate polar data
- Parameters
data (
numpy.ndarray
) – 2-dimensional array (azimuth, ranges) of floats;if no mask is assigned explicitly polar data should be a masked array
mask (
numpy.ndarray
, optional) – boolean array with pixels to be interpolated set to True, must have the same shape as data, defaults to None.ipclass (
wradlib.ipol.IpolBase
) – A class which inherits from IpolBase, defaults towradlib.ipol.Nearest
.
- Returns
filled_data (
numpy.ndarray
) – 2D array with interpolated values for the values set to True in the mask
Examples
>>> import numpy as np # noqa >>> import wradlib as wrl >>> # creating a data array and mask some values >>> data = np.arange(12.).reshape(4,3) >>> masked_values = (data==2) | (data==9) >>> # interpolate the masked data based on ''masked_values'' >>> filled_a = wrl.ipol.interpolate_polar(data, mask = masked_values, ipclass = wrl.ipol.Linear) # noqa >>> da = wrl.georef.create_xarray_dataarray(filled_a) >>> da = da.wrl.georef.georeference() >>> pm = wrl.vis.plot(da) >>> # the same result can be achieved by using a masked array instead of an explicit mask # noqa >>> mdata = np.ma.array(data, mask = masked_values) >>> filled_b = wrl.ipol.interpolate_polar(mdata, ipclass = wrl.ipol.Linear) # noqa >>> da = wrl.georef.create_xarray_dataarray(filled_b) >>> da = da.wrl.georef.georeference() >>> pm = wrl.vis.plot(da)