Python测量值插值生成效率图：插值谷值问题求助

问题概述 #

需基于(x,y,z)格式的测量数据生成符合实际物理规律的效率图，测试了RBF、LinearND、NearestND等插值方法后，均出现不符合预期的异常谷值，推测核心问题是平台两侧需插值至零的边界约束未被插值算法正确识别。

现有测试结果 #

• RBF插值：插值曲面出现明显不符合实际的谷值；
• LinearND插值：结果最接近预期，但仍存在局部异常谷值；
• NearestND插值：曲面过渡生硬，不符合效率图的平滑要求。

原始测量数据为(x,y,z)矩阵，现有测试代码如下：

import numpy as np  
import matplotlib.pyplot as plt  
import pandas as pd 
from scipy.interpolate import RBFInterpolator, LinearNDInterpolator, NearestNDInterpolator

df = pd.read_csv("Data.csv") 
mat = df.to_numpy()

x = mat[:,0]
y = mat[:,1]
z = mat[:,2]

# 绘制原始数据点
fig, ax = plt.subplots(subplot_kw={"projection": "3d"})
ax.set_title('Data points')
ax.scatter(x, y, np.transpose(z), marker='o', c='black', s=2)
ax.view_init(30, -115)
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('z')

num_points = 500
xi = np.linspace(min(x),max(x),num_points)
yi = np.linspace(min(y),max(y),num_points)
XI, YI = np.meshgrid(xi, yi) 

# RBF插值
measured_points = np.stack([x.ravel(), y.ravel()], -1)
interpolated_points = np.stack([XI.ravel(), YI.ravel()], -1)
interpolation_rbf = RBFInterpolator(measured_points, z.ravel(), smoothing=0, kernel='linear')
z_rbf = interpolation_rbf(interpolated_points).reshape(num_points, num_points)    
         
# 绘制RBF插值结果
fig, ax = plt.subplots(subplot_kw={"projection": "3d"})
ax.set_title('RBFInterpolator')
surf = ax.plot_surface(XI, YI, z_rbf, cmap='jet')   
ax.scatter(x, y, np.transpose(z), marker='o', c='black', s=2)
ax.view_init(30, -115)
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('z')

# LinearND插值
interp = LinearNDInterpolator(list(zip(x, y)), z)
ZI = interp(XI, YI)

# 绘制LinearND插值结果
fig, ax = plt.subplots(subplot_kw={"projection": "3d"})
ax.set_title('LinearNDInterpolator')
surf = ax.plot_surface(XI, YI, ZI, cmap='jet')   
ax.scatter(x, y, np.transpose(z), marker='o', c='black', s=2)
ax.view_init(30, -115)
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('z')

# NearestND插值
interp = NearestNDInterpolator(list(zip(x, y)), z)
ZI = interp(XI, YI)

# 绘制NearestND插值结果
fig, ax = plt.subplots(subplot_kw={"projection": "3d"})
ax.set_title('NearestNDInterpolator')
surf = ax.plot_surface(XI, YI, ZI, cmap='jet')   
ax.scatter(x, y, np.transpose(z), marker='o', c='black', s=2)
ax.view_init(30, -115)
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('z')

相关可视化结果：

• 测量值分布：3D散点图展示原始测量点的分布情况；
• RBF插值结果：曲面存在明显异常谷值；
• LinearND插值结果：整体趋势接近目标，但局部仍有异常谷值；
• 目标效率图样式：平滑连续的曲面，无不符合实际的波动，符合物理规律。

优化方案建议 #

1. 手动添加边界约束点 #

明确告知插值算法平台两侧的零值边界，避免算法无约束地生成异常值。示例代码如下：

# 获取原始数据的边界范围
x_min, x_max = np.min(x), np.max(x)
y_min, y_max = np.min(y), np.max(y)

# 生成边界点：在x的最小/最大值处，生成一系列y值对应的z=0点
boundary_y = np.linspace(y_min, y_max, 50)
# x_min处的边界点
boundary_x_min = np.full_like(boundary_y, x_min)
boundary_z_min = np.zeros_like(boundary_y)
# x_max处的边界点
boundary_x_max = np.full_like(boundary_y, x_max)
boundary_z_max = np.zeros_like(boundary_y)

# 合并原始数据与边界点
x_extended = np.concatenate([x, boundary_x_min, boundary_x_max])
y_extended = np.concatenate([y, boundary_y, boundary_y])
z_extended = np.concatenate([z, boundary_z_min, boundary_z_max])

# 使用扩展后的数据进行插值（以RBF为例）
measured_points_extended = np.stack([x_extended.ravel(), y_extended.ravel()], -1)
interpolation_rbf = RBFInterpolator(measured_points_extended, z_extended.ravel(), smoothing=0.2, kernel='cubic')
z_rbf = interpolation_rbf(interpolated_points).reshape(num_points, num_points)

2. 调整RBF插值的内核与平滑参数 #

更换非线性内核（如cubic或quintic）并设置合理的平滑值，减少过度拟合导致的异常波动：

# 示例：使用三次内核+0.5的平滑参数
interpolation_rbf = RBFInterpolator(measured_points_extended, z_extended.ravel(), smoothing=0.5, kernel='cubic')

3. 采用克里金插值（Kriging） #

克里金插值更适合空间数据的插值，能考虑数据的空间相关性，有效避免异常谷值。需先安装pykrige库：

pip install pykrige

示例代码：

from pykrige.ok import OrdinaryKriging

# 创建普通克里金模型，选择合适的变异函数模型（如gaussian）
OK = OrdinaryKriging(x_extended, y_extended, z_extended, variogram_model='gaussian')
# 执行网格插值
z_krige, _ = OK.execute('grid', xi, yi)

# 绘制克里金插值结果
fig, ax = plt.subplots(subplot_kw={"projection": "3d"})
ax.set_title('Kriging Interpolation')
surf = ax.plot_surface(XI, YI, z_krige, cmap='jet')   
ax.scatter(x, y, np.transpose(z), marker='o', c='black', s=2)
ax.view_init(30, -115)
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('z')

4. 对插值结果进行后处理平滑 #

使用高斯滤波对插值后的曲面进行平滑处理，消除局部异常波动：

from scipy.ndimage import gaussian_filter

# 对RBF插值结果进行平滑，sigma控制平滑程度
z_smoothed = gaussian_filter(z_rbf, sigma=2)

# 绘制平滑后的结果
fig, ax = plt.subplots(subplot_kw={"projection": "3d"})
ax.set_title('Smoothed RBF Interpolation')
surf = ax.plot_surface(XI, YI, z_smoothed, cmap='jet')   
ax.scatter(x, y, np.transpose(z), marker='o', c='black', s=2)
ax.view_init(30, -115)
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('z')

内容的提问来源于stack exchange，提问作者Patrick