Note
Go to the end to download the full example code
2D ndsplines vs. scipy.interpolateΒΆ
7 import numpy as np
8 import gc
9 import time
10
11 from scipy import interpolate
12 import ndsplines
13
14 import matplotlib.pyplot as plt
15
16 # number of time measurements per input/query size
17 n_iter = 10
18
19
20 def timeit(func, n_iter=1, return_samps=True, **func_kwargs):
21 results = np.empty(n_iter, dtype=np.double)
22 for i in range(n_iter):
23 # gc.collect()
24
25 tstart = time.time()
26 func(**func_kwargs)
27 delta = time.time() - tstart
28
29 results[i] = delta
30
31 if return_samps:
32 return results
33 else:
34 return np.mean(results)
35
36
37 def gen_xyz(sizex, sizey):
38 x = np.pi * np.linspace(-1, 1, sizex)
39 y = np.pi * np.linspace(-1, 1, sizey)
40 xx, yy = np.meshgrid(x, y, indexing='ij')
41
42 zz = np.sin(xx) * np.cos(yy)
43 return x, y, zz
44
45
46 def gen_xxyy(sizex, sizey):
47 x = np.pi * np.linspace(-1, 1, sizex)
48 y = np.pi * np.linspace(-1, 1, sizey)
49 xx, yy = np.meshgrid(x, y, indexing='ij')
50 return xx, yy
51
52
53 # make_interp_spline timing
54 x_sizes = np.logspace(1, 3, 10, dtype=int)
55 t_scipy_build = np.empty((2, x_sizes.size))
56 t_ndspl_build = np.empty((2, x_sizes.size))
57 for i, size in enumerate(x_sizes):
58 x, y, z = gen_xyz(size, size)
59 t_scipy = 10e3 * timeit(interpolate.RectBivariateSpline, x=x.copy(), y=y.copy(), z=z.copy(),
60 n_iter=n_iter)
61 t_ndspl = 10e3 * timeit(ndsplines.make_interp_spline, x=[x,y], y=z,
62 n_iter=n_iter)
63 t_scipy_build[:, i] = np.mean(t_scipy), np.std(t_scipy)
64 t_ndspl_build[:, i] = np.mean(t_ndspl), np.std(t_ndspl)
65
66 # spline query timing
67 x, y, z = gen_xyz(7, 5)
68 xx_sizes = np.logspace(0, 2, 10, dtype=int)
69 t_scipy_call = np.empty((2, xx_sizes.size))
70 t_ndspl_npy_call = np.empty((2, xx_sizes.size))
71 t_ndspl_pyx_call = np.empty((2, xx_sizes.size))
72 for i, size in enumerate(xx_sizes):
73 xx, yy = gen_xxyy(size, size)
74 xxyy = np.stack((xx, yy), axis=-1)
75 spl_scipy = interpolate.RectBivariateSpline(x.copy(), y.copy(), z)
76 spl_ndspl = ndsplines.make_interp_spline((x,y), z)
77 spl_ndspl.allocate_workspace_arrays(size)
78 t_scipy = 10e3 * timeit(spl_scipy, x=xx.copy(), y=yy.copy(), grid=False, n_iter=n_iter)
79 ndsplines.set_impl('cython')
80 t_ndspl_pyx = 10e3 * timeit(spl_ndspl, x=xxyy,
81 n_iter=n_iter)
82 ndsplines.set_impl('numpy')
83 t_ndspl_npy = 10e3 * timeit(spl_ndspl, x=xxyy,
84 n_iter=n_iter)
85 t_scipy_call[:, i] = np.mean(t_scipy), np.std(t_scipy)
86 t_ndspl_pyx_call[:, i] = np.mean(t_ndspl_pyx), np.std(t_ndspl_pyx)
87 t_ndspl_npy_call[:, i] = np.mean(t_ndspl_npy), np.std(t_ndspl_npy)
88
89 # plot results
90 fig, axes = plt.subplots(nrows=2)
91
92 axes[0].errorbar(x_sizes, t_scipy_build[0], capsize=3, yerr=t_scipy_build[1],
93 label='scipy')
94 axes[0].errorbar(x_sizes, t_ndspl_build[0], capsize=3, yerr=t_ndspl_build[1],
95 label='ndsplines')
96 axes[0].set_title('make_interp_spline')
97
98 axes[1].errorbar(xx_sizes, t_scipy_call[0], capsize=3, yerr=t_scipy_call[1],
99 label='scipy.interpolate')
100 axes[1].errorbar(xx_sizes, t_ndspl_npy_call[0], capsize=3, yerr=t_ndspl_npy_call[1],
101 label='ndsplines npy')
102 axes[1].errorbar(xx_sizes, t_ndspl_pyx_call[0], capsize=3, yerr=t_ndspl_pyx_call[1],
103 label='ndsplines pyx')
104 axes[1].set_title('spline.__call__')
105
106 for ax in axes:
107 ax.set_xlabel('input array size')
108 ax.set_ylabel('time [ms]')
109 ax.set_xscale('log')
110 ax.grid()
111
112 axes[-1].legend()
113 fig.tight_layout()
114
115 plt.show()
Total running time of the script: (0 minutes 7.679 seconds)