我想用scipy做一些插值。我看过很多例子,但我没有找到我想要的。
假设我有一些数据,其中的row和column变量可以在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点并确定插值值。我知道我可以用地图坐标做这个。我想知道是否有任何简单/聪明的方法可以为数据数组中的适当索引生成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
我想你想要一个bivariate spline on a rectangular structured mesh:
相关问题 更多 >
编程相关推荐