Note
Go to the end to download the full example code or to run this example in your browser via Binder
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#
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.129623
pyinterp.Axis 1.138777
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: 43.244655
pyinterp.Axis 6.945435
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 53.662 seconds)
