xarray.DataArray.interp#
- DataArray.interp(coords=None, method='linear', assume_sorted=False, kwargs=None, **coords_kwargs)[source]#
Interpolate a DataArray onto new coordinates.
Performs univariate or multivariate interpolation of a Dataset onto new coordinates, utilizing either NumPy or SciPy interpolation routines.
Out-of-range values are filled with NaN, unless specified otherwise via kwargs to the numpy/scipy interpolant.
- Parameters:
coords (
dict
, optional) – Mapping from dimension names to the new coordinates. New coordinate can be a scalar, array-like or DataArray. If DataArrays are passed as new coordinates, their dimensions are used for the broadcasting. Missing values are skipped.method (
{"linear", "nearest", "zero", "slinear", "quadratic", "cubic", "quintic", "polynomial", "pchip", "barycentric", "krogh", "akima", "makima"}
) – Interpolation method to use (see descriptions above).assume_sorted (
bool
, default:False
) – If False, values of x can be in any order and they are sorted first. If True, x has to be an array of monotonically increasing values.kwargs (dict-like or
None
, default:None
) – Additional keyword arguments passed to scipy’s interpolator. Valid options and their behavior depend whetherinterp1d
orinterpn
is used.**coords_kwargs (
{dim: coordinate, ...}
, optional) – The keyword arguments form ofcoords
. One of coords or coords_kwargs must be provided.
- Returns:
interpolated (
DataArray
) – New dataarray on the new coordinates.
Notes
SciPy is required for certain interpolation methods.
- When interpolating along multiple dimensions with methods linear and nearest,
the process attempts to decompose the interpolation into independent interpolations along one dimension at a time.
- The specific interpolation method and dimensionality determine which
interpolant is used:
- Interpolation along one dimension of 1D data (`method=’linear’`)
Uses
numpy.interp()
, unless fill_value=’extrapolate’ is provided via kwargs.
- Interpolation along one dimension of N-dimensional data (N ≥ 1)
- Methods {“linear”, “nearest”, “zero”, “slinear”, “quadratic”, “cubic”, “quintic”, “polynomial”}
use
scipy.interpolate.interp1d()
, unless conditions permit the use ofnumpy.interp()
(as in the case of method=’linear’ for 1D data).
If method=’polynomial’, the order keyword argument must also be provided.
- Special interpolants for interpolation along one dimension of N-dimensional data (N ≥ 1)
- Depending on the method, the following interpolants from
scipy.interpolate
are used: “pchip”:
scipy.interpolate.PchipInterpolator
“barycentric”:
scipy.interpolate.BarycentricInterpolator
“krogh”:
scipy.interpolate.KroghInterpolator
- “akima” or “makima”:
scipy.interpolate.Akima1dInterpolator
(makima is handled by passing the makima flag).
- “akima” or “makima”:
- Depending on the method, the following interpolants from
- Interpolation along multiple dimensions of multi-dimensional data
- Uses
scipy.interpolate.interpn()
for methods {“linear”, “nearest”, “slinear”, “cubic”, “quintic”, “pchip”}.
- Uses
See also
- Manipulating Dimensions (Data Resolution)
Tutorial material on manipulating data resolution using
interp()
Examples
>>> da = xr.DataArray( ... data=[[1, 4, 2, 9], [2, 7, 6, np.nan], [6, np.nan, 5, 8]], ... dims=("x", "y"), ... coords={"x": [0, 1, 2], "y": [10, 12, 14, 16]}, ... ) >>> da <xarray.DataArray (x: 3, y: 4)> Size: 96B array([[ 1., 4., 2., 9.], [ 2., 7., 6., nan], [ 6., nan, 5., 8.]]) Coordinates: * x (x) int64 24B 0 1 2 * y (y) int64 32B 10 12 14 16
1D linear interpolation (the default):
>>> da.interp(x=[0, 0.75, 1.25, 1.75]) <xarray.DataArray (x: 4, y: 4)> Size: 128B array([[1. , 4. , 2. , nan], [1.75, 6.25, 5. , nan], [3. , nan, 5.75, nan], [5. , nan, 5.25, nan]]) Coordinates: * y (y) int64 32B 10 12 14 16 * x (x) float64 32B 0.0 0.75 1.25 1.75
1D nearest interpolation:
>>> da.interp(x=[0, 0.75, 1.25, 1.75], method="nearest") <xarray.DataArray (x: 4, y: 4)> Size: 128B array([[ 1., 4., 2., 9.], [ 2., 7., 6., nan], [ 2., 7., 6., nan], [ 6., nan, 5., 8.]]) Coordinates: * y (y) int64 32B 10 12 14 16 * x (x) float64 32B 0.0 0.75 1.25 1.75
1D linear extrapolation:
>>> da.interp( ... x=[1, 1.5, 2.5, 3.5], ... method="linear", ... kwargs={"fill_value": "extrapolate"}, ... ) <xarray.DataArray (x: 4, y: 4)> Size: 128B array([[ 2. , 7. , 6. , nan], [ 4. , nan, 5.5, nan], [ 8. , nan, 4.5, nan], [12. , nan, 3.5, nan]]) Coordinates: * y (y) int64 32B 10 12 14 16 * x (x) float64 32B 1.0 1.5 2.5 3.5
2D linear interpolation:
>>> da.interp(x=[0, 0.75, 1.25, 1.75], y=[11, 13, 15], method="linear") <xarray.DataArray (x: 4, y: 3)> Size: 96B array([[2.5 , 3. , nan], [4. , 5.625, nan], [ nan, nan, nan], [ nan, nan, nan]]) Coordinates: * x (x) float64 32B 0.0 0.75 1.25 1.75 * y (y) int64 24B 11 13 15