cupyx.scipy.interpolate.RegularGridInterpolator#
- class cupyx.scipy.interpolate.RegularGridInterpolator(points, values, method='linear', bounds_error=True, fill_value=nan)[source]#
Interpolator on a regular or rectilinear grid in arbitrary dimensions.
The data must be defined on a rectilinear grid; that is, a rectangular grid with even or uneven spacing. Linear, nearest-neighbor, spline interpolations are supported. After setting up the interpolator object, the interpolation method may be chosen at each evaluation.
- Parameters:
points (tuple of ndarray of float, with shapes (m1, ), ..., (mn, )) – The points defining the regular grid in n dimensions. The points in each dimension (i.e. every elements of the points tuple) must be strictly ascending or descending.
values (ndarray, shape (m1, ..., mn, ...)) – The data on the regular grid in n dimensions.
method (str, optional) – The method of interpolation to perform. Supported are “linear”, “nearest”, “slinear”, “cubic”, “quintic” and “pchip”. This parameter will become the default for the object’s
__call__
method. Default is “linear”.bounds_error (bool, optional) – If True, when interpolated values are requested outside of the domain of the input data, a ValueError is raised. If False, then fill_value is used. Default is True.
fill_value (float or None, optional) – The value to use for points outside of the interpolation domain. If None, values outside the domain are extrapolated. Default is
cp.nan
.
Notes
Contrary to scipy’s LinearNDInterpolator and NearestNDInterpolator, this class avoids expensive triangulation of the input data by taking advantage of the regular grid structure.
In other words, this class assumes that the data is defined on a rectilinear grid.
The ‘slinear’(k=1), ‘cubic’(k=3), and ‘quintic’(k=5) methods are tensor-product spline interpolators, where k is the spline degree, If any dimension has fewer points than k + 1, an error will be raised.
If the input data is such that dimensions have incommensurate units and differ by many orders of magnitude, the interpolant may have numerical artifacts. Consider rescaling the data before interpolating.
** Choosing a spline method **
Spline methods, “slinear”, “cubic” and “quintic” involve solving a large sparse linear system at instantiation time. Alternatively, you may instead use the legacy methods, “slinear_legacy”, “cubic_legacy” and “quintic_legacy”. These methods allow faster construction but evaluations will be much slower.
Examples
Evaluate a function on the points of a 3-D grid
As a first example, we evaluate a simple example function on the points of a 3-D grid:
>>> from cupyx.scipy.interpolate import RegularGridInterpolator >>> import cupy as cp >>> def f(x, y, z): ... return 2 * x**3 + 3 * y**2 - z >>> x = cp.linspace(1, 4, 11) >>> y = cp.linspace(4, 7, 22) >>> z = cp.linspace(7, 9, 33) >>> xg, yg ,zg = cp.meshgrid(x, y, z, indexing='ij', sparse=True) >>> data = f(xg, yg, zg)
data
is now a 3-D array withdata[i, j, k] = f(x[i], y[j], z[k])
. Next, define an interpolating function from this data:>>> interp = RegularGridInterpolator((x, y, z), data)
Evaluate the interpolating function at the two points
(x,y,z) = (2.1, 6.2, 8.3)
and(3.3, 5.2, 7.1)
:>>> pts = cp.array([[2.1, 6.2, 8.3], ... [3.3, 5.2, 7.1]]) >>> interp(pts) array([ 125.80469388, 146.30069388])
which is indeed a close approximation to
>>> f(2.1, 6.2, 8.3), f(3.3, 5.2, 7.1) (125.54200000000002, 145.894)
Interpolate and extrapolate a 2D dataset
As a second example, we interpolate and extrapolate a 2D data set:
>>> x, y = cp.array([-2, 0, 4]), cp.array([-2, 0, 2, 5]) >>> def ff(x, y): ... return x**2 + y**2
>>> xg, yg = cp.meshgrid(x, y, indexing='ij') >>> data = ff(xg, yg) >>> interp = RegularGridInterpolator((x, y), data, ... bounds_error=False, fill_value=None)
>>> import matplotlib.pyplot as plt >>> fig = plt.figure() >>> ax = fig.add_subplot(projection='3d') >>> ax.scatter(xg.ravel().get(), yg.ravel().get(), data.ravel().get(), ... s=60, c='k', label='data')
Evaluate and plot the interpolator on a finer grid
>>> xx = cp.linspace(-4, 9, 31) >>> yy = cp.linspace(-4, 9, 31) >>> X, Y = cp.meshgrid(xx, yy, indexing='ij')
>>> # interpolator >>> ax.plot_wireframe(X.get(), Y.get(), interp((X, Y)).get(), rstride=3, cstride=3, alpha=0.4, color='m', label='linear interp')
>>> # ground truth >>> ax.plot_wireframe(X.get(), Y.get(), ff(X, Y).get(), rstride=3, cstride=3, ... alpha=0.4, label='ground truth') >>> plt.legend() >>> plt.show()
See also
scipy.interpolate.RegularGridInterpolator
interpn
a convenience function which wraps RegularGridInterpolator
scipy.ndimage.map_coordinates
interpolation on grids with equal spacing (suitable for e.g., N-D image resampling)
References
- [1] Python package regulargrid by Johannes Buchner, see
- [2] Wikipedia, “Trilinear interpolation”,
- [3] Weiser, Alan, and Sergio E. Zarantonello. “A note on piecewise
linear and multilinear table interpolation in many dimensions.” MATH. COMPUT. 50.181 (1988): 189-196. https://www.ams.org/journals/mcom/1988-50-181/S0025-5718-1988-0917826-0/S0025-5718-1988-0917826-0.pdf
Methods
- __call__(xi, method=None, *, nu=None)[source]#
Interpolation at coordinates.
- Parameters:
xi (cupy.ndarray of shape (..., ndim)) – The coordinates to evaluate the interpolator at.
method (str, optional) – The method of interpolation to perform. Supported are “linear”, “nearest”, “slinear”, “cubic”, “quintic” and “pchip”. Default is the method chosen when the interpolator was created.
nu (sequence of ints, length ndim, optional) – If not None, the orders of the derivatives to evaluate. Each entry must be non-negative. Only allowed for methods “slinear”, “cubic” and “quintic”.
- Returns:
values_x – Interpolated values at xi. See notes for behaviour when
xi.ndim == 1
.- Return type:
cupy.ndarray, shape xi.shape[:-1] + values.shape[ndim:]
Notes
In the case that
xi.ndim == 1
a new axis is inserted into the 0 position of the returned array, values_x, so its shape is instead(1,) + values.shape[ndim:]
.Examples
Here we define a nearest-neighbor interpolator of a simple function
>>> import cupy as cp >>> x, y = cp.array([0, 1, 2]), cp.array([1, 3, 7]) >>> def f(x, y): ... return x**2 + y**2 >>> data = f(*cp.meshgrid(x, y, indexing='ij', sparse=True)) >>> from cupyx.scipy.interpolate import RegularGridInterpolator >>> interp = RegularGridInterpolator((x, y), data, method='nearest')
By construction, the interpolator uses the nearest-neighbor interpolation
>>> interp([[1.5, 1.3], [0.3, 4.5]]) array([2., 9.])
We can however evaluate the linear interpolant by overriding the method parameter
>>> interp([[1.5, 1.3], [0.3, 4.5]], method='linear') array([ 4.7, 24.3])
- __eq__(value, /)#
Return self==value.
- __ne__(value, /)#
Return self!=value.
- __lt__(value, /)#
Return self<value.
- __le__(value, /)#
Return self<=value.
- __gt__(value, /)#
Return self>value.
- __ge__(value, /)#
Return self>=value.