2

我有一个 3D 数据集,我想对其进行线性插值和外推。插值可以很容易地完成scipy.interpolate.LinearNDInterpolator。该模块只能为参数范围之外的值填充常量/nan,但我不明白为什么它不提供打开外推的选项。

查看代码,我看到模块是用 cython 编写的。由于没有 cython 经验,因此很难使用代码来实现外推。我可以用纯python代码编写它,但也许这里的其他人有更好的主意?我的特殊情况涉及一个恒定的 xy 网格,但 z 值不断变化很大(-100,000),因此插值必须很快,因为每次 z 值变化时都会运行插值。

根据要求给出一个基本示例,假设我有一个像

xyPairs = [[-1.0, 0.0], [-1.0, 4.0],
           [-0.5, 0.0], [-0.5, 4.0],
           [-0.3, 0.0], [-0.3, 4.0],
           [+0.0, 0.0], [+0.0, 4.0],
           [+0.2, 0.0], [+0.2, 4.0]]

假设我想计算x = -1.5, -0.8, +0.5和的值y = -0.2, +0.2, +0.5。目前,我对每个 y 值沿 x 轴执行 1d 插值/外推,然后对每个 x 值沿 y 轴执行一维插值/外插。外推由 中的第二个函数完成ryggyr's answer

4

3 回答 3

4

我提出了一个方法,代码很糟糕,但我希望它会对你有所帮助。这个想法是,如果您提前知道必须外推的界限,您可以在数组边缘添加额外的列/行,并使用线性外推值,然后在新数组上进行内插。这是一个示例,其中一些数据将被外推到 x=+-50 和 y=+-40:

import numpy as np
x,y=np.meshgrid(np.linspace(0,6,7),np.linspace(0,8,9)) # create x,y grid
z=x**2*y # and z values
# create larger versions with two more columns/rows
xlarge=np.zeros((x.shape[0]+2,x.shape[1]+2))
ylarge=np.zeros((x.shape[0]+2,x.shape[1]+2))
zlarge=np.zeros((x.shape[0]+2,x.shape[1]+2))
xlarge[1:-1,1:-1]=x # copy data on centre
ylarge[1:-1,1:-1]=y
zlarge[1:-1,1:-1]=z
# fill extra columns/rows
xmin,xmax=-50,50
ymin,ymax=-40,40
xlarge[:,0]=xmin;xlarge[:,-1]=xmax # fill first/last column
xlarge[0,:]=xlarge[1,:];xlarge[-1,:]=xlarge[-2,:] # copy first/last row
ylarge[0,:]=ymin;ylarge[-1,:]=ymax
ylarge[:,0]=ylarge[:,1];ylarge[:,-1]=ylarge[:,-2]
# for speed gain: store factor of first/last column/row
first_column_factor=(xlarge[:,0]-xlarge[:,1])/(xlarge[:,1]-xlarge[:,2]) 
last_column_factor=(xlarge[:,-1]-xlarge[:,-2])/(xlarge[:,-2]-xlarge[:,-3])
first_row_factor=(ylarge[0,:]-ylarge[1,:])/(ylarge[1,:]-ylarge[2,:])
last_row_factor=(ylarge[-1,:]-ylarge[-2,:])/(ylarge[-2,:]-ylarge[-3,:])
# extrapolate z; this operation only needs to be repeated when zlarge[1:-1,1:-1] is updated
zlarge[:,0]=zlarge[:,1]+first_column_factor*(zlarge[:,1]-zlarge[:,2]) # extrapolate first column
zlarge[:,-1]=zlarge[:,-2]+last_column_factor*(zlarge[:,-2]-zlarge[:,-3]) # extrapolate last column
zlarge[0,:]=zlarge[1,:]+first_row_factor*(zlarge[1,:]-zlarge[2,:]) # extrapolate first row
zlarge[-1,:]=zlarge[-2,:]+last_row_factor*(zlarge[-2,:]-zlarge[-3,:]) #extrapolate last row

然后您可以在 (xlarge,ylarge,zlarge) 上进行插值。由于所有操作都是 numpy 切片操作,我希望它对您来说足够快。当 z 数据更新时,将它们复制进去zlarge[1:-1,1:-1]并重新执行最后 4 行。

于 2014-01-31T15:04:15.660 回答
4

使用最近插值和线性插值的组合。如果LinearNDInterpolator 插值失败,则返回np.nan,否则返回数组大小(1) NearestNDInterpolator 返回浮点数

import scipy.interpolate
import numpy
class LinearNDInterpolatorExt(object):
  def __init__(self, points,values):
    self.funcinterp=scipy.interpolate.LinearNDInterpolator(points,values)
    self.funcnearest=scipy.interpolate.NearestNDInterpolator(points,values)
  def __call__(self,*args):
    t=self.funcinterp(*args)
    if not numpy.isnan(t):
      return t.item(0)
    else:
      return self.funcnearest(*args)
于 2018-02-12T23:54:17.737 回答
0

我稍微修改了@Keith Williams 的回答,这对我来说效果很好(注意它不能线性推断——它只使用最近的邻居):

import numpy as np
from scipy.interpolate import LinearNDInterpolator as linterp
from scipy.interpolate import NearestNDInterpolator as nearest

class LinearNDInterpolatorExt(object):
    def __init__(self, points, values):
        self.funcinterp = linterp(points, values)
        self.funcnearest = nearest(points, values)
    
    def __call__(self, *args):
        z = self.funcinterp(*args)
        chk = np.isnan(z)
        if chk.any():
            return np.where(chk, self.funcnearest(*args), z)
        else:
            return z
于 2020-12-07T05:45:15.327 回答