interpolate

  • Interpolators for spatial data

Calling Sequence

import pyTMD.interpolate
data = pyTMD.interpolate.bilinear(ilon, ilat, idata, lon, lat)

Source code

pyTMD.interpolate.bilinear(ilon: ~numpy.ndarray, ilat: ~numpy.ndarray, idata: ~numpy.ndarray, lon: ~numpy.ndarray, lat: ~numpy.ndarray, fill_value: float = nan, dtype: str | numpy.dtype = <class 'numpy.float64'>)[source]

Bilinear interpolation of input data to output coordinates

Parameters
ilon: np.ndarray

longitude of tidal model

ilat: np.ndarray

latitude of tidal model

idata: np.ndarray

tide model data

lat: np.ndarray

output latitude

lon: np.ndarray

output longitude

fill_value: float, default np.nan

invalid value

dtype: np.dtype, default np.float64

output data type

Returns
data: np.ndarray

interpolated data

pyTMD.interpolate.spline(ilon: ~numpy.ndarray, ilat: ~numpy.ndarray, idata: ~numpy.ndarray, lon: ~numpy.ndarray, lat: ~numpy.ndarray, fill_value: float = None, dtype: str | numpy.dtype = <class 'numpy.float64'>, reducer=<ufunc 'ceil'>, **kwargs)[source]

Bivariate spline interpolation of input data to output coordinates

Parameters
ilon: np.ndarray

longitude of tidal model

ilat: np.ndarray

latitude of tidal model

idata: np.ndarray

tide model data

lat: np.ndarray

output latitude

lon: np.ndarray

output longitude

fill_value: float or NoneType, default None

invalid value

dtype: np.dtype, default np.float64

output data type

reducer: obj, default np.ceil

operation for converting mask to boolean

kx: int, default 1

degree of the bivariate spline in the x-dimension

ky: int, default 1

degree of the bivariate spline in the y-dimension

kwargs: dict

additional arguments for scipy.interpolate.RectBivariateSpline

Returns
data: np.ndarray

interpolated data

pyTMD.interpolate.regulargrid(ilon: ~numpy.ndarray, ilat: ~numpy.ndarray, idata: ~numpy.ndarray, lon: ~numpy.ndarray, lat: ~numpy.ndarray, fill_value: float = None, dtype: str | numpy.dtype = <class 'numpy.float64'>, reducer=<ufunc 'ceil'>, **kwargs)[source]

Regular grid interpolation of input data to output coordinates

Parameters
ilon: np.ndarray

longitude of tidal model

ilat: np.ndarray

latitude of tidal model

idata: np.ndarray

tide model data

lat: np.ndarray

output latitude

lon: np.ndarray

output longitude

fill_value: float or NoneType, default None

invalid value

dtype: np.dtype, default np.float64

output data type

reducer: obj, default np.ceil

operation for converting mask to boolean

bounds_error: bool, default False

raise Exception when values are requested outside domain

method: str, default ‘linear’

Method of interpolation

  • 'linear'

  • 'nearest'

  • 'slinear'

  • 'cubic'

  • 'quintic'

kwargs: dict

additional arguments for scipy.interpolate.RegularGridInterpolator

Returns
data: np.ndarray

interpolated data

pyTMD.interpolate.extrapolate(ilon: ~numpy.ndarray, ilat: ~numpy.ndarray, idata: ~numpy.ndarray, lon: ~numpy.ndarray, lat: ~numpy.ndarray, fill_value: float = None, dtype: str | numpy.dtype = <class 'numpy.float64'>, cutoff: int | float = inf, EPSG: str | int = '4326', **kwargs)[source]

Nearest-neighbor (NN) extrapolation of valid model data using kd-trees

Parameters
x: np.ndarray

x-coordinates of tidal model

y: np.ndarray

y-coordinates of tidal model

data: np.ndarray

Tide model data

XI: np.ndarray

Output x-coordinates

YI: np.ndarray

Output y-coordinates

fill_value: float, default np.nan

Invalid value

dtype: np.dtype, default np.float64

Output data type

cutoff: float, default np.inf

return only neighbors within distance [km]

Set to np.inf to extrapolate for all points

EPSG: str, default ‘4326’

projection of tide model data

Returns
DATA: np.ndarray

interpolated data