我正在开发一个项目,该项目在很大程度上依赖于2D插值。经过一些初步搜索,我发现很多文章都说地图坐标或可能的RegularGridInterpolator可以提供最快的速度,只要我能接受结构化网格。我编写了以下代码试用scipy中的所有多变量插值函数,结果表明RectBivariateSpline是最佳选择。
from __future__ import division
import numpy as np
import time
import math
from scipy.interpolate import RectBivariateSpline, RegularGridInterpolator, interp2d, CloughTocher2DInterpolator, Rbf, griddata, interpn
from scipy.ndimage.interpolation import map_coordinates
import random
grid_size = 100
xrange = np.linspace(0, 1, grid_size)
yrange = np.linspace(0, 1, grid_size)
x2d, y2D = np.meshgrid(xrange, yrange, indexing='ij')
z2D = x2d ** 2.0 + y2D
i2d = interp2d(xrange, yrange, z2D, 'linear')
rbs = RectBivariateSpline(xrange, yrange, z2D)
RGI = RegularGridInterpolator((xrange, yrange), z2D)
CT2d = CloughTocher2DInterpolator(np.array(list(zip(x2d.flatten(), y2D.flatten()))), z2D.flatten())
rbfi = Rbf(x2d, y2D, z2D, function='linear')
grd_inputs1 = np.array(list(zip(x2d.flatten(), y2D.flatten())))
grd_inputs2 = z2D.flatten()
print 'Exact: ',0.5**2.0 + 0.5
print 'interp2d',i2d(0.5, 0.5)[0]
print 'RectBivariateSpline',rbs(0.5, 0.5)[0][0]
print 'RegularGridInterpolator',RGI([0.5, 0.5])[0]
coords = np.asarray(([0.5], [0.5]))
coords = [(c - lo) * (n - 1) / (hi - lo) for (lo, hi), c, n in zip([[0,1],[0,1]], coords, np.shape(z2D.flatten()))]
print 'Map Coordinates:', map_coordinates(z2D.flatten(), coords, order=1)[0]
print 'Rbf', rbfi([0.5, 0.5])[0]
print 'CloughTocher2DInterpolator', CT2d([0.5, 0.5])[0]
print 'griddata', griddata(grd_inputs1, grd_inputs2, [0.5, 0.5], 'linear')[0]
print 'interpn', interpn([xrange, yrange], z2D, [0.5, 0.5])[0]
print ''
samples = 25000
start = time.time()
for n in range(samples):
i2d(0.5, 0.5)
end = time.time()
print 'interp2d: %0.4f [us]' % (((end-start)/samples)*1e6)
start = time.time()
for n in range(samples):
rand = random.uniform(0,1)
rbs(0.5, 0.5)
end = time.time()
print 'RectBivariateSpline: %0.4f [us]' % (((end-start)/samples)*1e6)
start = time.time()
for n in range(samples):
rand = random.uniform(0,1)
RGI([0.5, 0.5])
end = time.time()
print 'RegularGridInterpolator: %0.4f [us]' % (((end-start)/samples)*1e6)
start = time.time()
for n in range(samples):
map_coordinates(z2D.flatten(), coords, order=1)
end = time.time()
print 'Map Coordiantes: %0.4f [us]' % (((end-start)/samples)*1e6)
start = time.time()
for n in range(samples):
rand = random.uniform(0,1)
rbfi([0.5, 0.5])[0]
end = time.time()
print 'Rbf: %0.4f [us]' % (((end-start)/samples)*1e6)
start = time.time()
for n in range(samples):
rand = random.uniform(0,1)
CT2d([0.5, 0.5])[0]
end = time.time()
print 'CloughTocher2DInterpolator: %0.4f [us]' % (((end-start)/samples)*1e6)
start = time.time()
for n in range(int(samples/100)):
rand = random.uniform(0,1)
griddata(grd_inputs1, grd_inputs2, [0.5, 0.5], 'linear')
end = time.time()
print 'griddata: %0.4f [us]' % (((end-start)/int(samples/100))*1e6)
start = time.time()
for n in range(int(samples)):
rand = random.uniform(0,1)
#grdd(250e3+rand, bar2Pa(125+rand))[0]
interpn([xrange, yrange], z2D, [0.5, 0.5])[0]
end = time.time()
print 'interpn: %0.4f [us]' % (((end-start)/int(samples))*1e6)
产生以下结果:
Exact: 0.75
interp2d 0.7500255076012651
RectBivariateSpline 0.7499999999999997
RegularGridInterpolator 0.750025507601265
Map Coordinates: 0.7500255076012652
Rbf 0.7500000434270742
CloughTocher2DInterpolator 0.750003857158422
griddata 0.7500255076012653
interpn 0.750025507601265
interp2d: 14.4000 [us]
RectBivariateSpline: 2.8400 [us]
RegularGridInterpolator: 90.9200 [us]
Map Coordiantes: 8.9600 [us]
Rbf: 105.8000 [us]
CloughTocher2DInterpolator: 18.0000 [us]
griddata: 138824.0004 [us]
interpn: 158.1600 [us]
我的实现有什么问题吗?这似乎难以置信,因为RectBivariateSpline是高阶拟合,并且相对于地图坐标允许非均匀网格间距。还有其他更快的二维插值方法吗?请原谅这段代码中可能存在的任何显著问题,我和Python刚刚开始接触。
谢谢,