Python双线性插值和双三次插值

python二维数组的插值

通过scipy.interpolate中的griddata可以进行针对坐标网格的二维插值，其调用方法为

griddata(points, values, xi, method='linear', fill_value=nan, rescale=False)
1

points, values构成了用于插值的原始数据，xi为插值的坐标格点，method为插值方法，fill_value为超出插值范围时的填充值。

griddata有三种插值方法，分别是线性插值、三次插值和临近插值。接下来借用官网文档的例程，来对比一下这三种不同的插值方案

import numpy as np
import matplotlib.pyplot as plt
from scipy.interpolate import griddata
# 插值函数
def func(x, y):
    return x*(1-x)*np.cos(4*np.pi*x) * np.sin(4*np.pi*y**2)**2

rng = np.random.default_rng()
pts = rng.random((1000, 2))
values = func(pts[:,0], pts[:,1])

xs, ys = np.indices([100,100])/100
zs = func(xs, ys)

grids = {}
keys = ["nearest", "linear", "cubic"]
for key in keys:
    grids[key] = griddata(pts, values, (xs,ys), method=key)

fig = plt.figure()
ax = fig.add_subplot(2,2,1)
ax.imshow(zs.T, extent=(0,1,0,1), origin='lower')
plt.scatter(pts[:,0], pts[:,1], marker='.', c='black')
ax.set_title('Original')
plt.axis('off')

for i,key in zip([2,3,4], keys):
    ax = fig.add_subplot(2,2,i,projection='3d')
    ax.plot_surface(xs, ys, grids[key])
    ax.set_title(key)

plt.show()
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32

效果如下

可以看到，三次插值最为光滑，线性插值次之，邻近插值则没有光滑性可言。

基本原理

之所以出现这种差别，当然是插值的原理不同所致。

所谓临近插值，就是选择距离最近的点作为插值点，以一维情况为例，若待插值序列为

那么在$[10,15)$这个区间中，所有插值点的值都将是10，相当于是0次插值。下图中，红、绿、蓝、灰分别代表0到3次插值，可见，尽管只有10个点，但分段的二次函数已经描绘出了三角函数的形状，其插值效果还是不错的。

0次插值在二维情况下是临近点插值，这很好理解，那么线性插值变成双线性插值就有些复杂，毕竟在一维数组中，插入一个点，会在左右各有一个最相近的点，由这两个点可以连成一条直线，从而完成某种近似的线性拟合。

但在二维的待插值网格中，插入一点，那么这个点的周围会有4个临近点，而一个平面只需要三个点，所以不能用平面来类比直线。

一般在双线性插值中，比较常用的方案类似于权重法，即通过距离来判定某个点的权重来进行计算。
在矩阵中，索引坐标是规整的整数，现有一浮点型的点$(x,y)$，则距离其左下角最近的坐标格点为

$$x_L=\lfloor x\rfloor y_L=\lfloor y\rfloor$$

记矩阵为$M$，则矩阵中对应的值为

$$M(x_L,y_L)$$

记$x_R=x_L+1,y_R=y_L+1$，则可得到$(x,y)$附近的四个点的坐标

[ M ( x L , y R ) M ( x R , y R ) M ( x L , y L ) M ( x R , y L ) ] [M(xL,yR)M(xR,yR)M(xL,yL)M(xR,yL)]

[M(xL​,yR​)M(xL​,yL​)​M(xR​,yR​)M(xR​,yL​)​]

记作

[ M L R M R R M L L M R L ] [MLRMRRMLLMRL]

[MLR​MLL​​MRR​MRL​​]

点$(x,y)$在这四个点中间，相当于分别受到这四个点的影响，即

M ( x , y ) = ∑ M i j ∣ ( x − x i ˉ ) ( y − y j ˉ ) ∣ i , j ∈ { L , R } M(x,y)=\sum M_{ij}\vert(x-x_{\bar i})(y-y_{\bar j})\vert\quad i,j\in\{L,R\} M(x,y)=∑Mij​∣(x−xiˉ​)(y−yjˉ​​)∣i,j∈{L,R}

其中，$\bar{L}=R,\bar{R}=L$，即对于左下角$x_L,y_L$而言，其对$(x,y)$的影响正比于$(x,y)$距离其右上方的点的距离。当$(x,y)$和$(x_L,y_L)$重合时，$|(x-x_R)(y-y_R)|=1$。

由于矩阵是规整的，所以$|(x-x_{\bar{i}})|=1-|(x-x_i)|$，则上式可以写为

$$M(x,y)=\sum M_{ij}(1-|x-x_i|)(1-|y-y_j|)$$