Python样条插值和坐标映射

样条插值基本原理

由于数组本身是格点化了的，所以对数组的旋转、平移和缩放，并不像实数空间中那么简单。以一维的平移为例，现有三个点，坐标为0,1,2，值对应为a,b,c，现在将这个数轴向右平移0.5，则这a,b,c三个点的坐标就变成了0.5, 1.5, 2.5，但数组中不允许出现非整数的坐标，所以显示的仍然是0,1,2位置处的值，所以就需要通过0.5, 1.5, 2.5位置处的值，来插值得到0,1,2位置处的值。

几何变换中所用到的插值方法，是样条插值。

所谓样条插值，简单来说，就是根据插值次数，每$N$个点确定一个分段函数，每个分段函数都是一个样条。将所有的样条组合起来，就是样条插值。最简单的样条插值函数是线性样条插值，可表示为

Si(x)=yi+yi+1−yixi+1−xi(x−xi)S_i(x)=y_i+\\frac{y_{i+1}-y_i}{x_{i+1}-x_i}(x-x_i) Si​(x)=yi​+xi+1​−xi​yi+1​−yi​​(x−xi​)

scipy实现

在scipy.ndimage中，提供了一维样条插值和多维样条插值函数，分别是spline_filter1d和spline_filter，二者参数如下，其中order, output, mode均与几何变换中的定义相同

spline_filter1d(input, order=3, axis=-1, output=np.float64, mode='mirror')
spline_filter(input, order=3, output=np.float64, mode='mirror')

由于插值本身

import numpy as np
import matplotlib.pyplot as plt
from scipy.ndimage import *
ori = np.eye(20)  # create an image
ori[10, :] = 1.0fDct = {"ori" : ori,"axis_0" : spline_filter1d(ori, axis=0),"axis_1" : spline_filter1d(ori, axis=1),"multi axis, order=3" : spline_filter(ori, order=3),"multi axis, order=4" : spline_filter(ori, order=4),"multi axis, order=5" : spline_filter(ori, order=5)
}fig = plt.figure()
for i, key in enumerate(fDct):fig.add_subplot(2, 3, i+1)plt.imshow(fDct[key], cmap='gray_r')plt.title(key)plt.show()

效果如下

坐标映射

坐标映射的逻辑是，将输入的数组，通过插值的方式，映射到新的坐标中。例如现在有一个矩阵数组

[01234567891011]\\begin{bmatrix} 0 & 1 & 2 \\\\ 3 & 4 & 5 \\\\ 6 & 7 & 8 \\\\ 9 &10 &11 \\end{bmatrix} ​0369​14710​25811​​

在python中坐标原点从0开始，所以(0,0)处的值是0；(1,1)处的值是4，那么(0.5,2)这个点并不在矩阵的格点上，需要通过插值得到

import numpy as np
from scipy.ndimage import *
a = np.arange(12.).reshape((4, 3))
map_coordinates(a, [[0.5, 2], [0.5, 1]], order=1)
# 返回array([ 2.,  7.])

其中函数map_coordinates极为坐标映射函数，除了输入的数组a，以及要提取的坐标位置之外，还支持output, order, mode, cval以及prefilter这5个参数，参数作用与旋转平移等函数相同。