SciPy - 插值


插值方法很多,本文只是總結下 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()

輸出

 

 

參考資料:

https://www.cnblogs.com/xiuercui/p/12292563.html


免責聲明!

本站轉載的文章為個人學習借鑒使用,本站對版權不負任何法律責任。如果侵犯了您的隱私權益,請聯系本站郵箱yoyou2525@163.com刪除。



 
粵ICP備18138465號   © 2018-2025 CODEPRJ.COM