薄板样条插值法

Python实现薄板样条插值

最近在做有关卫星高度计的课题，其中涉及到同潮图的绘制。需要利用某海域测高卫星星下点调和常数，进行插值得到目标海域的网格数据，然后绘制同潮图。
某天翻阅论文发现某篇文章对薄板样条插值方法描述比较详细，就试着用Python实现了一下。（水平不高，仅供参考）

以下是代码:
代码中主要采用矩阵计算（因为矩阵计算真的比常规的循环方法要快很多）

    def TPS_interpolation(lon,lat,height,loc,num):
        """薄板插值(采用lon, lat = np.mgrid[99:122:0.5, 2:25:0.5]构造网格)
           lon:格网经度；lat:格网纬度；height:已有高程；
           loc:已有高程点的坐标数据[[lon1,lat1],[lon2,lat2].....]；
           num:loc数据集的列数，比如我的数据就是2列"""

        #构造矩阵
          #T矩阵
        T_a = np.ones(len(loc))
        T_b1 = np.hsplit(loc,num)
        T_b2 = np.mat(np.hstack((T_b1[0],T_b1[1]))).T
        T = np.vstack((T_a,T_b2))
          #E矩阵
        t1 = np.vsplit(loc,len(loc))
        t2 = list()
        for i in t1:
            t3 = list()
            
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17