一个1D的numpy/scipy插值并不完全是1D

Question

Numpy数组data有。。。。数据。Numpy数组z有距离。数据和z具有相同的形状，z的每个点都是测量数据对应点的距离。更复杂的是，用户将提供3维、4维或5维的data/z数组。在

我想将数据插值到1D numpy数组dists中的一组距离。由于数据结构的原因，插值轴始终是距末端的两个轴，即，如果数组有3维，则插值轴为0；如果数组有4维，则插值轴为1，依此类推。由于AFAICT，所有numpy/scipy插值例程都希望在1D数组中给出原始距离，插值数据和z超距离似乎是一个有点复杂的任务。这就是我所拥有的：

def dist_interp(data, z, dists):
    # construct array to hold interpolation results
    num_dims = len(data.shape)
    interp_axis = num_dims-3
    interp_shape = list(data.shape)
    interp_shape[interp_axis] = dists.shape[0]
    interp_res = np.zeros(shape=interp_shape)
    # depending on usage, data could have between 3 and five dimensions.
    # add dims to generalize.  I hate doing it this way.  Must be
    # some other way.
    for n in range(num_dims, 5) :
        data = np.expand_dims(data, axis=0)
        z = np.expand_dims(z, axis=0)
        interp_res = np.expand_dims(interp_res, axis=0)
    for m in range(data.shape[0]):
        for l in range(data.shape[1]):
            for j in range(data.shape[3]):
                for i in range(data.shape[4]):
                    interp_res[m,l,:,j,i]=(
                                  np.interp(dists,z[m,l,:,j,i],
                                            data[m,l,:,j,i]))
    # now remove extra "wrapping" dimensions
    for n in range(0,5-num_dims):
        interp_res = interp_res[0]
    return(interp_res)

但我认为去除这些额外的维度并不是很好。有更好的主意吗？谢谢。在

Answer 1

对于这种情况，一般的经验法则是：

1. 将要执行操作的轴移动到形状元组的末尾
2. 将生成的阵列重塑为二维
3. 从这个二维数组创建新数组
4. 撤消新阵列上的步骤2和1

对于你的情况，它可能看起来像：

# Create some random test data
axis = -2
shape = np.random.randint(10, size=(5,))
data = np.random.rand(*shape)
data = np.sort(data, axis=axis)
z = np.random.rand(*shape)
dists = np.linspace(0,1, num=100)

data = np.rollaxis(data, axis, data.ndim)
new_shape = data.shape
data = data.reshape(-1, data.shape[-1])
z = np.rollaxis(z, axis, z.ndim)
z = z.reshape(-1, z.shape[-1])

out = np.empty(z.shape[:1]+dists.shape)
for o, x, f in zip(out, data, z):
    o[:] = np.interp(dists, x, f)

out = out.reshape(new_shape[:-1]+dists.shape)
out = np.rollaxis(out, -1, axis)