4种插值算法

如果有一些稀疏的轨迹，如何将这些轨迹平滑插值处理呢？

方法1：线性插值

方法2：三次样条插值

方法3：贝塞尔曲线插值

方法4：拉格朗日插值

线性插值：在两两相邻的点之间差值，两点间全部插值点都在一条直线上。

贝塞尔曲线：贝塞尔曲线不会经过全部的坐标点，会根据坐标点的排列趋势去拟合出一条相对平滑的从第1个点到末了一个点之间的曲线。

三次样条插值：插值函数会经过全部的坐标点，拟合函数平滑。

拉格朗日插值：点太多，会出现不稳定的结果。见下图：

前三种插值算法都有特定的使用场景，按需使用就好了。

1. import time
3. import numpy as np
4. from scipy.interpolate import interp1d
5. from scipy.special import comb
8. def linear_interpolation(route, num_points):
9.     # 1. 线性插值
10.     # 将经纬度分开
11.     lons = np.array([point[0] for point in route])
12.     lats = np.array([point[1] for point in route])
14.     # 创建插值函数
15.     distance = np.cumsum(np.sqrt(np.ediff1d(lons, to_begin=0) ** 2 + np.ediff1d(lats, to_begin=0) ** 2))
16.     distance /= distance[-1]
18.     # 创建插值函数
19.     lon_interp = interp1d(distance, lons, kind='linear')
20.     lat_interp = interp1d(distance, lats, kind='linear')
22.     # 生成新的距离点
23.     new_distance = np.linspace(0, 1, num_points)
25.     # 插值
26.     new_lons = lon_interp(new_distance)
27.     new_lats = lat_interp(new_distance)
29.     return list(zip(new_lons, new_lats))
32. def cubic_spline_interpolation(route, num_points):
33.     # 2. 三次样条插值
34.     lons = np.array([point[0] for point in route])
35.     lats = np.array([point[1] for point in route])
37.     distance = np.cumsum(np.sqrt(np.ediff1d(lons, to_begin=0) ** 2 + np.ediff1d(lats, to_begin=0) ** 2))
38.     distance /= distance[-1]
40.     lon_interp = interp1d(distance, lons, kind='cubic')
41.     lat_interp = interp1d(distance, lats, kind='cubic')
43.     new_distance = np.linspace(0, 1, num_points)
45.     new_lons = lon_interp(new_distance)
46.     new_lats = lat_interp(new_distance)
48.     return list(zip(new_lons, new_lats))
51. def bernstein_poly(i, n, t):
52.     return comb(n, i) * (t ** i) * ((1 - t) ** (n - i))
55. def bezier_curve(route, num_points=100):
56.     # 3. 贝塞尔曲线插值
57.     n = len(route) - 1
58.     t = np.linspace(0, 1, num_points)
60.     lons = np.array([point[0] for point in route])
61.     lats = np.array([point[1] for point in route])
63.     new_lons = np.zeros(num_points)
64.     new_lats = np.zeros(num_points)
66.     for i in range(n + 1):
67.         new_lons += bernstein_poly(i, n, t) * lons[i]
68.         new_lats += bernstein_poly(i, n, t) * lats[i]
70.     return list(zip(new_lons, new_lats))
73. import matplotlib.pyplot as plt
76. def plot_routes(original, linear, cubic, bezier):
77.     plt.figure(figsize=(12, 8))
79.     # 原始轨迹
80.     orig_lons, orig_lats = zip(*original)
81.     plt.plot(orig_lons, orig_lats, 'ro-', label='Original', alpha=0.5)
83.     # 线性插值
84.     lin_lons, lin_lats = zip(*linear)
85.     plt.plot(lin_lons, lin_lats, 'b-', label='Linear', alpha=0.7)
87.     # 三次样条插值
88.     cub_lons, cub_lats = zip(*cubic)
89.     plt.plot(cub_lons, cub_lats, 'g-', label='Cubic Spline', alpha=0.7)
91.     # 贝塞尔曲线
92.     bez_lons, bez_lats = zip(*bezier)
93.     plt.plot(bez_lons, bez_lats, 'm-', label='Bezier', alpha=0.7)
95.     plt.legend()
96.     plt.xlabel('Longitude')
97.     plt.ylabel('Latitude')
98.     plt.title('Trajectory Interpolation Comparison')
99.     plt.grid(True)
100.     plt.show()
103. if __name__ == '__main__':
104.     # 原始轨迹数据
105.     route = [
106.         [122.123456, 31.123456],
107.         [122.234567, 31.234567],
108.         [122.345678, 31.345678],
109.         [122.456789, 31.456789],
110.         [122.567890, 31.567890],
111.         [122.678901, 31.578901],
112.         [122.789012, 31.789012],
113.         [122.890123, 31.890123],
114.         [122.901234, 31.901234],
115.     ]
116.     start_time = time.time()
117.     # 线性插值
118.     linear_route = linear_interpolation(route, 1000)
119.     print("线性插值结果 (前5个点):", linear_route[:5])
120.     print("线性插值用时:", time.time() - start_time, "秒")
122.     start_time = time.time()
123.     # 三次样条插值
124.     cubic_route = cubic_spline_interpolation(route, 1000)
125.     print("三次样条插值结果 (前5个点):", cubic_route[:5])
126.     print("三次样条插值用时:", time.time() - start_time, "秒")
128.     start_time = time.time()
129.     # 贝塞尔曲线插值
130.     bezier_route = bezier_curve(route, 1000)
131.     print("贝塞尔曲线插值结果 (前5个点):", bezier_route[:5])
132.     print("贝塞尔曲线插值用时:", time.time() - start_time, "秒")
134.     # 绘制比较图
135.     plot_routes(route, linear_route, cubic_route, bezier_route)

复制代码

1. # 4、拉格朗日插值算法
2. import time
4. from scipy.interpolate import lagrange
5. import numpy as np
8. def lagrange_interp(x, y, x_new):
9.     """
10.     Lagrange interpolation
11.     :param x: x coordinates of data points
12.     :param y: y coordinates of data points
13.     :param x_new: x coordinates of new interpolated points
14.     :return: y coordinates of new interpolated points
15.     """
16.     f = lagrange(x, y)
17.     y_new = f(x_new)
18.     return y_new
21. if __name__ == '__main__':
22.     # 原始数据
23.     route = [
24.         [122.123456, 31.123456],
25.         [122.234567, 31.234567],
26.         [122.345678, 31.345678],
27.         [122.456789, 31.456789],
28.         [122.567890, 31.567890],
29.         [122.678901, 31.678901],
30.         [122.789012, 31.789012],
31.         [122.890123, 31.890123],
32.         [122.901234, 31.901234],
34.     ]
36.     x_list = [i[0] for i in route]
37.     y_list = [i[1] for i in route]
38.     # 新数据
39.     x_new = np.arange(122.123456, 122.990123, 0.01)
40.     y_new = lagrange_interp(x_list, y_list, x_new)
41.     # 绘图
42.     import matplotlib.pyplot as plt
44.     plt.plot(x_list, y_list, 'o', label='original data')
45.     plt.plot(x_new, y_new, label='interpolated data')
46.     plt.legend()
47.     plt.show()

复制代码

免责声明：如果侵犯了您的权益，请联系站长，我们会及时删除侵权内容，谢谢合作！更多信息从访问主页：qidao123.com:ToB企服之家，中国第一个企服评测及商务社交产业平台。