我有以下三个列表:
x = [100.0, 75.0, 50.0, 0.0, 100.0, 75.0, 50.0, 0.0, 100.0, 75.0, 50.0, 0.0, 100.0, 75.0, 50.0, 0.0, 100.0, 75.0, 50.0, 0.0, 100.0, 75.0, 50.0, 0.0, 100.0, 75.0, 50.0, 0.0, 100.0, 75.0, 50.0, 0.0, 100.0, 75.0, 50.0, 0.0, 100.0, 75.0, 50.0, 0.0]
y = [300.0, 300.0, 300.0, 300.0, 500.0, 500.0, 500.0, 500.0, 700.0, 700.0, 700.0, 700.0, 1000.0, 1000.0, 1000.0, 1000.0, 1500.0, 1500.0, 1500.0, 1500.0, 2000.0, 2000.0, 2000.0, 2000.0, 3000.0, 3000.0, 3000.0, 3000.0, 5000.0, 5000.0, 5000.0, 5000.0, 7500.0, 7500.0, 7500.0, 75000.0, 10000.0, 10000.0, 10000.0, 10000.0]
z = [100.0, 95.0, 87.5, 77.5, 60.0, 57.0, 52.5, 46.5, 40.0, 38.0, 35.0, 31.0, 30.0, 28.5, 26.25, 23.25, 23.0, 21.85, 20.125, 17.825, 17.0, 16.15, 14.875, 13.175, 13.0, 12.35, 11.375, 10.075, 10.0, 9.5, 8.75, 7.75, 7.0, 6.65, 6.125, 5.425, 5.0, 4.75, 4.375, 3.875]
每个列表的每个条目都被读取为一个点,因此点 0 是 (100,300,100),点 1 是 (75,300,95),依此类推。
我正在尝试进行 2d 插值,以便可以计算任何给定输入 (x0, y0)
点的 z 值。
我读到,使用 meshgrid 我可以使用 scipy 中的 RegularGridInterpolator
进行插值,但我不知道如何设置它:
x_,y_,z_ = np.meshgrid(x,y,z) # both indexing ij or xy
我没有得到有意义的 x_,y_,z_
值,并且我不知道如何从那里开始。
我正在尝试使用上面的数据点来查找中间值,因此类似于 scipy 的 interp1d
其中
f = interp1d(x, y, kind='cubic')
稍后我可以在其中调用f(范围内的任何(x,y)点)
并获取相应的z
值。
最佳答案
您需要对分散数据进行二维插值。在这种情况下,我默认使用 scipy.interpolate.griddata ,但您似乎想要一个可调用的插值器,而 griddata 需要给定的它将插入的点集。
不用担心:具有 2d 三次插值的 griddata
使用 CloughTocher2DInterpolator
。所以我们完全可以这样做:
import numpy as np
import scipy.interpolate as interp
x = [100.0, 75.0, 50.0, 0.0, 100.0, 75.0, 50.0, 0.0, 100.0, 75.0, 50.0, 0.0, 100.0, 75.0, 50.0, 0.0, 100.0, 75.0, 50.0, 0.0, 100.0, 75.0, 50.0, 0.0, 100.0, 75.0, 50.0, 0.0, 100.0, 75.0, 50.0, 0.0, 100.0, 75.0, 50.0, 0.0, 100.0, 75.0, 50.0, 0.0]
y = [300.0, 300.0, 300.0, 300.0, 500.0, 500.0, 500.0, 500.0, 700.0, 700.0, 700.0, 700.0, 1000.0, 1000.0, 1000.0, 1000.0, 1500.0, 1500.0, 1500.0, 1500.0, 2000.0, 2000.0, 2000.0, 2000.0, 3000.0, 3000.0, 3000.0, 3000.0, 5000.0, 5000.0, 5000.0, 5000.0, 7500.0, 7500.0, 7500.0, 75000.0, 10000.0, 10000.0, 10000.0, 10000.0]
z = [100.0, 95.0, 87.5, 77.5, 60.0, 57.0, 52.5, 46.5, 40.0, 38.0, 35.0, 31.0, 30.0, 28.5, 26.25, 23.25, 23.0, 21.85, 20.125, 17.825, 17.0, 16.15, 14.875, 13.175, 13.0, 12.35, 11.375, 10.075, 10.0, 9.5, 8.75, 7.75, 7.0, 6.65, 6.125, 5.425, 5.0, 4.75, 4.375, 3.875]
interpolator = interp.CloughTocher2DInterpolator(np.array([x,y]).T, z)
现在您可以使用 2 个坐标调用此插值器来为您提供相应的插值数据点:
>>> interpolator(x[10], y[10]) == z[10]
True
>>> interpolator(2, 300)
array(77.81343)
请注意,您必须留在输入点的凸包内,否则您将得到 nan
(或作为 fill_value
关键字传递给插值器):
>>> interpolator(2, 30)
array(nan)
无论如何,外推通常毫无意义,并且您的输入点以有点不稳定的方式分散:
所以即使推断是可能的,我也不会相信。
<小时/>只是为了演示生成的插值器如何约束到输入点的凸包,这是我们仅为绘图而创建的网格上的数据曲面图:
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
# go linearly in the x grid
xline = np.linspace(min(x), max(x), 30)
# go logarithmically in the y grid (considering y distribution)
yline = np.logspace(np.log10(min(y)), np.log10(max(y)), 30)
# construct 2d grid from these
xgrid,ygrid = np.meshgrid(xline, yline)
# interpolate z data; same shape as xgrid and ygrid
z_interp = interpolator(xgrid, ygrid)
# create 3d Axes and plot surface and base points
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.plot_surface(xgrid, ygrid, z_interp, cmap='viridis',
vmin=min(z), vmax=max(z))
ax.plot(x, y, z, 'ro')
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('z')
plt.show()
这是两个角度的输出(最好以交互方式旋转;这样的静态图像无法正确呈现 3D 表示):
有两个主要特点需要注意:
- 表面很好地拟合了红点,这是插值所期望的。幸运的是,输入点很好且平滑,因此插值一切顺利。 (红点通常被表面隐藏的事实只是由于 pyplot 的渲染器错误地处理了复杂 3D 对象的相对位置)
- 表面是沿着输入点的凸包切割的(由于
nan
值),因此即使我们的网格数组定义了一个矩形网格,我们也只能得到插值产生的表面切割感觉。
关于python - 给定分散的输入数据构建二维插值器,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/53036365/