Create interpolator objects#

In this example, we are going to build the basic objects allowing to carry out interpolations.

Before starting, we will examine the properties of a Cartesian grid and the different classes associated with these objects.

Step-by-step creation of grids#

import timeit

import numpy
import pandas

import pyinterp
import pyinterp.backends.xarray
import pyinterp.tests

ds = pyinterp.tests.load_grid3d()
lon, lat, time, tcw = (
    ds['longitude'].values,
    ds['latitude'].values,
    ds['time'].values,
    ds['tcw'].values,
)

This regular 3-dimensional grid is associated with three axes:

  • longitudes,

  • latitudes and

  • time.

To perform the calculations quickly, we will build three objects that will be used by the interpolator to search for the data to be used. Let’s start with the y-axis representing the latitude axis.

y_axis = pyinterp.Axis(lat)
y_axis
<pyinterp.core.Axis>
  min_value: -90
  max_value: 90
  step     : -2.5
  is_circle: false

For example, you can search for the closest point to 0.12 degrees north latitude.

y_axis.find_index([0.12])
array([36])

Then, the x-axis representing the longitudinal axis. In this case, the axis is an axis representing a 360 degree circle.

x_axis = pyinterp.Axis(lon, is_circle=True)
x_axis
<pyinterp.core.Axis>
  min_value: 0
  max_value: 357.5
  step     : 2.5
  is_circle: true

The values -180 and 180 degrees represent the same point on the axis.

x_axis.find_index([-180]) == x_axis.find_index([180])
array([ True])

Finally, we create the time axis

t_axis = pyinterp.TemporalAxis(time)
t_axis
<pyinterp.core.TemporalAxis>
  min_value: 2002-07-01T12:00:00.000000000
  max_value: 2002-07-31T12:00:00.000000000
  step     : 86400000000000 nanoseconds

As these objects must communicate in C++ memory space, we use objects specific to the library much faster than other data models and manage the axes representing a circle. For example if we compare these objects to Pandas indexes:

values = lon[10:20] + 1 / 3
index = pandas.Index(lon)
print('pandas.Index: %f' % timeit.timeit(
    'index.searchsorted(values)', globals=dict(index=index, values=values)))
print('pyinterp.Axis %f' % timeit.timeit(
    'x_axis.find_index(values)', globals=dict(x_axis=x_axis, values=values)))
pandas.Index: 2.324706
pyinterp.Axis 1.039706

This time axis is also very efficient compared to the pandas index.

index = pandas.Index(time)
values = time + numpy.timedelta64(1, 'ns')
print('pandas.Index: %f' % timeit.timeit(
    'index.searchsorted(values)', globals=dict(index=index, values=values)))
print('pyinterp.Axis %f' % timeit.timeit(
    't_axis.find_index(values)', globals=dict(t_axis=t_axis, values=values)))
pandas.Index: 41.536177
pyinterp.Axis 7.281678

Before constructing the tensor for pyinterp, we must begin to organize the tensor data so that it is properly stored in memory for pyinterp.

  • The shape of the tensor must be (len(x_axis), len(y_axis), len(t_axis))

Warning

If the array handled is a masked array, the masked values must be set to nan.

Now we can build the object handling the regular 3-dimensional grid.

Note

Grid data are not copied, the Grid3D class just keeps a reference on the handled array. Axis data are copied for non-uniform axes, and only examined for regular axes.

grid_3d = pyinterp.Grid3D(x_axis, y_axis, t_axis, tcw)
grid_3d
<pyinterp.grid.Grid3D>
array([[[10.15271338, 10.66710078, 11.29219624, ..., 10.60634637,
         15.88388013, 19.57369844],
        [12.08470387,  9.64372637, 12.96766807, ..., 10.47133655,
         12.06850269, 24.25583891],
        [15.77317208, 10.12436132, 11.59731843, ...,  9.53976881,
         11.00327523, 21.59209521],
        ...,
        [ 0.38745328,  0.44010711,  0.65747292, ...,  0.19438924,
          0.57376683,  0.36045132],
        [ 0.39420377,  0.35910122,  0.57106663, ...,  0.18763875,
          0.3739523 ,  0.27269494],
        [ 0.33344936,  0.28214562,  0.35640102, ...,  0.2227413 ,
          0.23759238,  0.20113973]],

       [[10.15271338, 10.66710078, 11.29219624, ..., 10.60634637,
         15.88388013, 19.57369844],
        [12.19541192,  9.65722735, 12.9460665 , ..., 10.34847761,
         12.54643744, 24.19508449],
        [15.72996894, 10.18511573, 11.41775537, ...,  9.63562578,
         11.89839032, 21.96202211],
        ...,
        [ 0.38340299,  0.43065642,  0.64937233, ...,  0.19573934,
          0.56026585,  0.34695034],
        [ 0.38880338,  0.35505093,  0.56566624, ...,  0.18763875,
          0.36855191,  0.26729454],
        [ 0.33344936,  0.28214562,  0.35640102, ...,  0.2227413 ,
          0.23759238,  0.20113973]],

       [[10.15271338, 10.66710078, 11.29219624, ..., 10.60634637,
         15.88388013, 19.57369844],
        [12.30611997,  9.67072833, 12.92311483, ..., 10.22561868,
         13.02707239, 24.13433007],
        [15.68946599, 10.24316996, 11.23954241, ...,  9.73283285,
         12.79485551, 22.33194901],
        ...,
        [ 0.37935269,  0.41985564,  0.63992164, ...,  0.19978964,
          0.54406467,  0.33344936],
        [ 0.38475309,  0.34830043,  0.56026585, ...,  0.18763875,
          0.36315151,  0.26324425],
        [ 0.33344936,  0.28214562,  0.35640102, ...,  0.2227413 ,
          0.23759238,  0.20113973]],

       ...,

       [[10.15271338, 10.66710078, 11.29219624, ..., 10.60634637,
         15.88388013, 19.57369844],
        [11.72422765,  9.73148275, 12.9447164 , ..., 10.88716679,
         11.01542611, 24.1167788 ],
        [15.4585992 , 10.30257428, 11.62972078, ...,  9.36560614,
          9.72743245, 19.62230197],
        ...,
        [ 0.40500456,  0.47250947,  0.71417704, ...,  0.21599081,
          0.62642066,  0.40770475],
        [ 0.40635466,  0.3739523 ,  0.57376683, ...,  0.21329062,
          0.38880338,  0.28889612],
        [ 0.33344936,  0.28214562,  0.35640102, ...,  0.2227413 ,
          0.23759238,  0.20113973]],

       [[10.15271338, 10.66710078, 11.29219624, ..., 10.60634637,
         15.88388013, 19.57369844],
        [11.84438639,  9.70178059, 12.95281699, ..., 10.74810667,
         11.36645164, 24.16403223],
        [15.56390686, 10.24316996, 11.61892   , ...,  9.42501046,
         10.15406348, 20.27844969],
        ...,
        [ 0.39825407,  0.46170868,  0.69527567, ...,  0.20789022,
          0.60886938,  0.39285367],
        [ 0.40230436,  0.36855191,  0.57376683, ...,  0.20383993,
          0.38475309,  0.28484582],
        [ 0.33344936,  0.28214562,  0.35640102, ...,  0.2227413 ,
          0.23759238,  0.20113973]],

       [[10.15271338, 10.66710078, 11.29219624, ..., 10.60634637,
         15.88388013, 19.57369844],
        [11.96454513,  9.67342853, 12.96091758, ..., 10.61039666,
         11.71612706, 24.20858547],
        [15.66786442, 10.18241554, 11.60676911, ...,  9.48306468,
         10.5779943 , 20.9345974 ],
        ...,
        [ 0.39285367,  0.45225799,  0.67502419, ...,  0.20248983,
          0.59131811,  0.3766525 ],
        [ 0.39825407,  0.36450161,  0.57376683, ...,  0.19573934,
          0.37800259,  0.27944543],
        [ 0.33344936,  0.28214562,  0.35640102, ...,  0.2227413 ,
          0.23759238,  0.20113973]]])
Axis:
 * x: <pyinterp.core.Axis>
       min_value: 0
       max_value: 357.5
       step     : 2.5
       is_circle: true
 * y: <pyinterp.core.Axis>
       min_value: -90
       max_value: 90
       step     : -2.5
       is_circle: false
 * z: <pyinterp.core.TemporalAxis>
       min_value: 2002-07-01T12:00:00.000000000
       max_value: 2002-07-31T12:00:00.000000000
       step     : 86400000000000 nanoseconds

xarray backend#

The construction of these objects manipulating the regular grids can be done more easily using the xarray library and CF convention usually found in NetCDF files.

interpolator = pyinterp.backends.xarray.RegularGridInterpolator(
    pyinterp.tests.load_grid3d().tcw)
interpolator.grid
<pyinterp.backends.xarray.Grid3D>
array([[[ 0.33344936,  0.28214562,  0.35640102, ...,  0.2227413 ,
          0.23759238,  0.20113973],
        [ 0.39420377,  0.35910122,  0.57106663, ...,  0.18763875,
          0.3739523 ,  0.27269494],
        [ 0.38745328,  0.44010711,  0.65747292, ...,  0.19438924,
          0.57376683,  0.36045132],
        ...,
        [15.77317208, 10.12436132, 11.59731843, ...,  9.53976881,
         11.00327523, 21.59209521],
        [12.08470387,  9.64372637, 12.96766807, ..., 10.47133655,
         12.06850269, 24.25583891],
        [10.15271338, 10.66710078, 11.29219624, ..., 10.60634637,
         15.88388013, 19.57369844]],

       [[ 0.33344936,  0.28214562,  0.35640102, ...,  0.2227413 ,
          0.23759238,  0.20113973],
        [ 0.38880338,  0.35505093,  0.56566624, ...,  0.18763875,
          0.36855191,  0.26729454],
        [ 0.38340299,  0.43065642,  0.64937233, ...,  0.19573934,
          0.56026585,  0.34695034],
        ...,
        [15.72996894, 10.18511573, 11.41775537, ...,  9.63562578,
         11.89839032, 21.96202211],
        [12.19541192,  9.65722735, 12.9460665 , ..., 10.34847761,
         12.54643744, 24.19508449],
        [10.15271338, 10.66710078, 11.29219624, ..., 10.60634637,
         15.88388013, 19.57369844]],

       [[ 0.33344936,  0.28214562,  0.35640102, ...,  0.2227413 ,
          0.23759238,  0.20113973],
        [ 0.38475309,  0.34830043,  0.56026585, ...,  0.18763875,
          0.36315151,  0.26324425],
        [ 0.37935269,  0.41985564,  0.63992164, ...,  0.19978964,
          0.54406467,  0.33344936],
        ...,
        [15.68946599, 10.24316996, 11.23954241, ...,  9.73283285,
         12.79485551, 22.33194901],
        [12.30611997,  9.67072833, 12.92311483, ..., 10.22561868,
         13.02707239, 24.13433007],
        [10.15271338, 10.66710078, 11.29219624, ..., 10.60634637,
         15.88388013, 19.57369844]],

       ...,

       [[ 0.33344936,  0.28214562,  0.35640102, ...,  0.2227413 ,
          0.23759238,  0.20113973],
        [ 0.40635466,  0.3739523 ,  0.57376683, ...,  0.21329062,
          0.38880338,  0.28889612],
        [ 0.40500456,  0.47250947,  0.71417704, ...,  0.21599081,
          0.62642066,  0.40770475],
        ...,
        [15.4585992 , 10.30257428, 11.62972078, ...,  9.36560614,
          9.72743245, 19.62230197],
        [11.72422765,  9.73148275, 12.9447164 , ..., 10.88716679,
         11.01542611, 24.1167788 ],
        [10.15271338, 10.66710078, 11.29219624, ..., 10.60634637,
         15.88388013, 19.57369844]],

       [[ 0.33344936,  0.28214562,  0.35640102, ...,  0.2227413 ,
          0.23759238,  0.20113973],
        [ 0.40230436,  0.36855191,  0.57376683, ...,  0.20383993,
          0.38475309,  0.28484582],
        [ 0.39825407,  0.46170868,  0.69527567, ...,  0.20789022,
          0.60886938,  0.39285367],
        ...,
        [15.56390686, 10.24316996, 11.61892   , ...,  9.42501046,
         10.15406348, 20.27844969],
        [11.84438639,  9.70178059, 12.95281699, ..., 10.74810667,
         11.36645164, 24.16403223],
        [10.15271338, 10.66710078, 11.29219624, ..., 10.60634637,
         15.88388013, 19.57369844]],

       [[ 0.33344936,  0.28214562,  0.35640102, ...,  0.2227413 ,
          0.23759238,  0.20113973],
        [ 0.39825407,  0.36450161,  0.57376683, ...,  0.19573934,
          0.37800259,  0.27944543],
        [ 0.39285367,  0.45225799,  0.67502419, ...,  0.20248983,
          0.59131811,  0.3766525 ],
        ...,
        [15.66786442, 10.18241554, 11.60676911, ...,  9.48306468,
         10.5779943 , 20.9345974 ],
        [11.96454513,  9.67342853, 12.96091758, ..., 10.61039666,
         11.71612706, 24.20858547],
        [10.15271338, 10.66710078, 11.29219624, ..., 10.60634637,
         15.88388013, 19.57369844]]])
Axis:
 * x: <pyinterp.core.Axis>
       min_value: 0
       max_value: 357.5
       step     : 2.5
       is_circle: true
 * y: <pyinterp.core.Axis>
       min_value: -90
       max_value: 90
       step     : 2.5
       is_circle: false
 * z: <pyinterp.core.TemporalAxis>
       min_value: 2002-07-01T12:00:00.000000000
       max_value: 2002-07-31T12:00:00.000000000
       step     : 86400000000000 nanoseconds

Total running time of the script: ( 0 minutes 52.391 seconds)

Gallery generated by Sphinx-Gallery