插值方法很多,本文只是总结下 scipy 库中插值用法
一维插值 - 拉格朗日插值
import numpy as np from scipy import interpolate import matplotlib.pylab as plt import pylab as mpl mpl.rcParams['font.sans-serif'] = ['FangSong'] # 指定默认字体 mpl.rcParams['axes.unicode_minus'] = False # 解决保存图像是负号'-'显示为方块的问题 ##################################### 一维数据插值 拉格朗日插值 x = np.linspace(0, 2*np.pi+np.pi/4, 10) y = np.sin(x) f_linear = interpolate.lagrange(x, y) # 拉格朗日插值 x_new = np.linspace(0, 2*np.pi+np.pi/4, 100) plt.plot(x, y, "o", label=u"原始数据") plt.plot(x_new, f_linear(x_new), label=u"拉格朗日插值") plt.xlabel(u'安培/A') plt.ylabel(u'伏特/V') plt.legend() plt.show()
输出
一维插值 - B样条插值
并且与 线性插值 做了对比
x = np.linspace(0, 2*np.pi+np.pi/4, 10) y = np.sin(x) f_linear = interpolate.interp1d(x, y) # 线性插值 tck = interpolate.splrep(x, y) # B-spline插值 x_new = np.linspace(0, 2*np.pi+np.pi/4, 100) y_bspline = interpolate.splev(x_new, tck) plt.plot(x, y, "o", label=u"原始数据") plt.plot(x_new, f_linear(x_new), label=u"线性插值") plt.plot(x_new, y_bspline, label=u"B-spline插值") plt.xlabel(u'安培/A') plt.ylabel(u'伏特/V') plt.legend() plt.show()
输出
一维插值 - 其他插值汇总
由于用法一样,这里统一记录
x = np.linspace(0, 10, 11) y = np.sin(x) plt.figure(figsize=(12,9)) plt.plot(x, y, 'ro') #建立插值数据点 xnew = np.linspace(0, 10, 101) for kind in ["nearest","zero","slinear","quadratic","cubic","linear"]:#插值方式 #"nearest","zero"为阶梯插值 #slinear 线性插值 #"quadratic","cubic" 为2阶、3阶B样条曲线插值 f = interpolate.interp1d(x, y, kind = kind) ynew = f(xnew)#计算插值结果 plt.plot(xnew, ynew, label = str(kind)) plt.xticks(fontsize=20) plt.yticks(fontsize=20) plt.legend(loc = 'lower right') plt.show()
输出
二维插值 - 样条插值
这里只是拿 样条 做个 demo,其他类型的插值 类似
def func(x, y): return (x+y)*np.exp(-5.0*(x**2 + y**2)) # X-Y轴分为15*15的网格 y, x = np.mgrid[-1:1:15j, -1:1:15j] fvals = func(x,y) # 计算每个网格点上的函数值 15*15的值 print(len(fvals[0])) # 三次样条二维插值 newfunc = interpolate.interp2d(x, y, fvals, kind='cubic') # x y z # 计算100*100的网格上的插值 xnew = np.linspace(-1,1,100)#x ynew = np.linspace(-1,1,100)#y fnew = newfunc(xnew, ynew)#仅仅是y值 100*100的值 ##### 绘图 # 为了更明显地比较插值前后的区别,使用关键字参数interpolation='nearest' # 关闭 imshow()内置的插值运算 plt.subplot(121) im1 = plt.imshow(fvals, extent=[-1,1,-1,1], cmap=mpl.cm.hot, interpolation='nearest', origin="lower") # pl.cm.jet # extent=[-1,1,-1,1]为x,y范围 plt.colorbar(im1) plt.subplot(122) im2 = plt.imshow(fnew, extent=[-1,1,-1,1], cmap=mpl.cm.hot, interpolation='nearest', origin="lower") plt.colorbar(im2) plt.show()
输出
二维插值 - 三维展示
from mpl_toolkits.mplot3d import Axes3D import matplotlib.cm as cm def func(x, y): return (x+y)*np.exp(-5.0*(x**2 + y**2)) # X-Y轴分为20*20的网格 x = np.linspace(-1, 1, 20) y = np.linspace(-1,1,20) x, y = np.meshgrid(x, y) # 20*20的网格数据 fvals = func(x,y) # 计算每个网格点上的函数值 15*15的值 fig = plt.figure(figsize=(9, 6)) ax = plt.subplot(1, 2, 1,projection = '3d') surf = ax.plot_surface(x, y, fvals, rstride=2, cstride=2, cmap=cm.coolwarm,linewidth=0.5, antialiased=True) ax.set_xlabel('x') ax.set_ylabel('y') ax.set_zlabel('f(x, y)') plt.colorbar(surf, shrink=0.5, aspect=5)# 标注 ##### 二维插值 newfunc = interpolate.interp2d(x, y, fvals, kind='cubic') # newfunc为一个函数 # 计算100*100的网格上的插值 xnew = np.linspace(-1,1,100)#x ynew = np.linspace(-1,1,100)#y fnew = newfunc(xnew, ynew)#仅仅是y值 100*100的值 np.shape(fnew) is 100*100 xnew, ynew = np.meshgrid(xnew, ynew) ax2 = plt.subplot(1, 2, 2,projection = '3d') surf2 = ax2.plot_surface(xnew, ynew, fnew, rstride=2, cstride=2, cmap=cm.coolwarm,linewidth=0.5, antialiased=True) ax2.set_xlabel('xnew') ax2.set_ylabel('ynew') ax2.set_zlabel('fnew(x, y)') plt.colorbar(surf2, shrink=0.5, aspect=5)#标注 plt.show()
输出
参考资料: