Python/Scipy 插值（map_coordinates）

Question

Python/Scipy 插值（map_coordinates）

12

我正在尝试使用scipy进行插值。我查看了许多例子，但没有找到完全符合我的要求的。

假设我有一些数据，其中行和列变量可以从0到1变化。每行和每列之间的间隔并不总是相同的（如下所示）。

      | 0.00   0.25  0.80  1.00
------|----------------------------
0.00  | 1.40   6.50  1.50  1.80
0.60  | 8.90   7.30  1.10  1.09
1.00  | 4.50   9.20  1.80  1.20

现在我想要能够通过一组x，y坐标确定插值的值。我知道可以使用map_coordinates来实现这一点。我想知道是否有任何简单/聪明的方法将x，y值转换为数据数组中相应的索引。例如，如果我输入x，y=0.60，0.25，则应该返回正确的插值索引。在这种情况下，这将是1.0，1.0，因为0.60，0.25会精确地映射到第二行和第二列。x=0.3会映射到0.5，因为它恰好在0.00和0.60之间的中间位置。我知道如何获得我想要的结果，但我确信有一个非常快速/清晰的一两行代码（或已经存在的函数）可以执行此操作，以使我的代码更清晰。基本上，它需要在某个数组中分段插值。这里是一个例子（严重基于Scipy interpolation on a numpy array的代码）-我在新函数将放置TODO：

from scipy.ndimage import map_coordinates
from numpy import arange
import numpy as np
#            0.000,  0.175,  0.817,  1.000
z = array([ [ 3.6,    6.5,    9.1,    11.5],    # 0.0000
            [ 3.9,   -7.3,    10.0,   13.1],    # 0.2620
            [ 1.9,    8.3,   -15.0,  -12.1],    # 0.6121
            [-4.5,    9.2,    12.2,   14.8] ])  # 1.0000
ny, nx = z.shape
xmin, xmax = 0., 1.
ymin, ymax = 0., 1.

xrange = array([0.000,  0.175,  0.817,  1.000 ])
yrange = array([0.0000, 0.2620, 0.6121, 1.0000])

# Points we want to interpolate at
x1, y1 = 0.20, 0.45
x2, y2 = 0.30, 0.85
x3, y3 = 0.95, 1.00

# To make our lives easier down the road, let's
# turn these into arrays of x & y coords
xi = np.array([x1, x2, x3], dtype=np.float)
yi = np.array([y1, y2, y3], dtype=np.float)

# Now, we'll set points outside the boundaries to lie along an edge
xi[xi > xmax] = xmax
xi[xi < xmin] = xmin

yi[yi > ymax] = ymax
yi[yi < ymin] = ymin

# We need to convert these to (float) indicies
#   (xi should range from 0 to (nx - 1), etc)
xi = (nx - 1) * (xi - xmin) / (xmax - xmin)
yi = (ny - 1) * (yi - ymin) / (ymax - ymin)
# TODO: Instead, xi and yi need to be mapped as described.  This can only work with
# even spacing...something like:
#xi = SomeInterpFunction(xi, xrange)
#yi = SomeInterpFunction(yi, yrange)

# Now we actually interpolate
# map_coordinates does cubic interpolation by default,
# use "order=1" to preform bilinear interpolation instead...
print xi
print yi
z1, z2, z3 = map_coordinates(z, [yi, xi], order=1)

# Display the results
for X, Y, Z in zip((x1, x2, x3), (y1, y2, y3), (z1, z2, z3)):
    print X, ',', Y, '-->', Z

- Scott B

1个回答

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Paul · Accepted Answer

我认为你想要在一个矩形结构网格上进行双变量样条插值：

import numpy
from scipy import interpolate
x = numpy.array([0.0, 0.60, 1.0])
y = numpy.array([0.0, 0.25, 0.80, 1.0])
z = numpy.array([ 
   [ 1.4 ,  6.5 ,  1.5 ,  1.8 ],
   [ 8.9 ,  7.3 ,  1.1 ,  1.09],
   [ 4.5 ,  9.2 ,  1.8 ,  1.2 ]])
# you have to set kx and ky small for this small example dataset
# 3 is more usual and is the default
# s=0 will ensure this interpolates.  s>0 will smooth the data
# you can also specify a bounding box outside the data limits
# if you want to extrapolate
sp = interpolate.RectBivariateSpline(x, y, z, kx=2, ky=2, s=0)

sp([0.60], [0.25])  # array([[ 7.3]])
sp([0.25], [0.60])  # array([[ 2.66427408]])