忽略NaN值的二维插值

5

如何告诉interp2d忽略NaN值?

我有一个带有任意值z的表面x和y。

x =   np.array([[9.19632, 9.62141, 10.0829, np.isnan, np.isnan],
        [9.21164, 9.64347, 10.1392, 10.5698,  np.isnan],
        [9.22175, 9.65439, 10.1423, 10.6301, 11.0323],
        [9.21632, 9.67060, 10.1474, 10.6230, 11.0818]])

y =  np.array([[11.5466,11.6485,11.7619, np.isnan, np.isnan],
        [12.4771, 12.5460, 12.5453, 12.7142, np.isnan],
        [13.5578, 13.5581, 13.5505, 13.5309, 13.6081],
        [14.5653, 14.5504, 14.5036, 14.5145, 14.5060]])

z = np.array([[0.466113, 0.0484404, -0.385355, np.isnan, np.isnan],
        [0.366125, -0.160165, -0.548668, -0.888301,np.isnan],
        [-0.0970777, -0.346734, -0.826576, -1.08412, -1.33129],
        [-0.259981, -0.586938, -1.03477, -1.32384, -1.61500]])

这里输入图像描述

我已经使用掩膜数组生成了上述颜色网格,但是当我尝试使用二维插值创建更精细的网格时失败了。下面是我的代码,注意我将nan值设置为零以获得这个结果,因此它显然会破坏“正确”的插值。相反,我想忽略它们并留空参数空间。

f = interp.interp2d(x,y,z, kind='linear')
xnew = np.arange(9,11.5, 0.01)
ynew = np.arange(9,15, 0.01)
znew = f(xnew, ynew)


levels = np.linspace(zmin, zmax, 15)
plt.ylabel('Y', size=15)
plt.xlabel('X', size=15)
cmap = plt.cm.jet_r
cmap.set_bad('white',0.1) # set nan to white
cs = plt.contourf(xnew, ynew, znew, levels=levels, cmap=cmap)
cbar = plt.colorbar(cs)
cbar.set_label('Z', rotation=90, fontsize=15) # gas fraction
plt.show()

我希望能够简单地创建一个平滑的颜色图,其中由x和y界定的区域根据z进行着色。

contour with nan set to zero

1个回答

7

我不知道为什么interp2d在处理不规则间隔数据时会出现问题,我建议使用griddata,您可以使用ravel将输入数据平铺成向量,然后消除NaN,并将其作为griddata的输入,您会得到类似这样的结果。

enter image description here

代码与您所拥有的并没有太大的区别

import numpy as np
import matplotlib.pyplot as plt
from scipy.interpolate import griddata

x =   np.array([[9.19632, 9.62141, 10.0829,np.isnan,np.isnan],
    [9.21164, 9.64347, 10.1392, 10.5698,np.isnan],
    [9.22175, 9.65439, 10.1423, 10.6301, 11.0323],
    [9.21632, 9.67060, 10.1474, 10.6230, 11.0818]])

y =  np.array([[11.5466,11.6485,11.7619,np.isnan,np.isnan],
    [12.4771, 12.5460, 12.5453, 12.7142,np.isnan],
    [13.5578, 13.5581, 13.5505, 13.5309, 13.6081],
    [14.5653, 14.5504, 14.5036, 14.5145, 14.5060]])

z = np.array([[0.466113, 0.0484404, -0.385355,np.isnan,np.isnan],
    [0.366125, -0.160165, -0.548668, -0.888301,np.isnan],
    [-0.0970777, -0.346734, -0.826576, -1.08412, -1.33129],
    [-0.259981, -0.586938, -1.03477, -1.32384, -1.61500]])

x=x.ravel()              #Flat input into 1d vector
x=list(x[x!=np.isnan])   #eliminate any NaN
y=y.ravel()
y=list(y[y!=np.isnan])
z=z.ravel()
z=list(z[z!=np.isnan])


xnew = np.arange(9,11.5, 0.01)
ynew = np.arange(9,15, 0.01)
znew = griddata((x, y), z, (xnew[None,:], ynew[:,None]), method='linear')



levels = np.linspace(min(z), max(z), 15)
plt.ylabel('Y', size=15)
plt.xlabel('X', size=15)
cmap = plt.cm.jet_r
cs = plt.contourf(xnew, ynew, znew, levels=levels, cmap=cmap)
cbar = plt.colorbar(cs)
cbar.set_label('Z', rotation=90, fontsize=15) # gas fraction
plt.show()

如果您必须外推数据(请查看我的下面的评论),您可以使用 SmoothBivariateSpline 并尝试调整样条的阶数,但我不建议这样做,我会向您展示原因。
代码更接近于您最初的代码。
from scipy.interpolate import SmoothBivariateSpline

x=x.ravel()
x=(x[x!=np.isnan])
y=y.ravel()
y=(y[y!=np.isnan])
z=z.ravel()
z=(z[z!=np.isnan])

xnew = np.arange(9,11.5, 0.01)
ynew = np.arange(10.5,15, 0.01)

f = SmoothBivariateSpline(x,y,z,kx=1,ky=1)

znew=np.transpose(f(xnew, ynew))

当kx=1,ky=1时,使用f = SmoothBivariateSpline(x,y,z,kx=1,ky=1),您将得到:

enter image description here

使用kx=2和ky=2,您会得到:

enter image description here

当kx=3且ky=3时,你会得到:

enter image description here

我调整了这三张图片的级别,以便更容易查看。但是请检查比例尺,数值在采样区域之外可能会非常快速地变得疯狂,因此如果您必须进行外推,请务必谨慎处理。


在二维外推方面,这是一个有点棘手的问题,具体取决于你正在做什么。如果你需要知道某个值是否在采样区域之外,那么使用NaN可能更好。但如果必须要外推,线性外推可能是一个比较安全的选择,但你很快就会超出范围。我将编辑答案,包括一些使用SmoothBivariateSpline进行外推的图像。 - Noel Segura Meraz
我同意,你必须谨慎处理二维外推问题,这是一个很好的例子。你能否编辑一下,展示一下如何使用SmoothBivariateSpline类呢?感谢你的帮助。 - smashbro
啊,转置(!)。现在有意义了 :) - smashbro
你为什么要从列表中删除 nan?我认为griddata可以处理它们。 - FooBar
@FooBar griddata 可以处理函数插值中的 nan 值,但不能作为坐标点。由于原始问题中也出现了 nan 作为坐标,因此将它们全部删除更容易一些。 - Noel Segura Meraz
显示剩余3条评论

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接