1

我有一个由间距定义的立方网格xi,yi,zi：

xi,yi,zi = [linspace(ox,ox+s*d,s) for ox,s,d in zip(origin,size,delta)]

我还在W该网格上设置了一组标量值。W.shape() == size. 我想使用scipy 的线性插值，它需要作为输入：

类scipy.interpolate.LinearNDInterpolator(points, values)：

参数 ：

points: ndarray of floats, shape(npoints, ndims)数据点坐标。

values：浮点数或复数的ndarray，形状(npoints, ...)数据值。

我如何创建一组假的points（通过神奇的广播）xi,yi,zi？现在我正在创建一个中间数组来馈送到插值函数 - 有没有更好的方法？

相关问题：3D 中的 Numpy meshgrid。这篇文章中的答案实际上创建了网格 - 我只想将其模拟为另一个函数的输入（首选纯 numpy 解决方案）。

4

3 回答 3

3 
>>> xi, yi, zi = [np.arange(3) for i in range(3)]
>>> xx, yy, zz = np.broadcast_arrays(xi,yi[:,np.newaxis],zi[:,np.newaxis,np.newaxis])
>>> xx.shape
(3, 3, 3)
>>> xx.strides
(0, 0, 8)

您可以看到它没有创建新副本，因为前两个维度的步幅为 0。

我也写了一个维度版本：

def ndmesh(*args):
   args = map(np.asarray,args)
   return np.broadcast_arrays(*[x[(slice(None),)+(None,)*i] for i, x in enumerate(args)])
于 2012-08-16T20:55:24.200 回答
1

points您可以按照其他答案中解释的类似方式构造必要的数组：

xx, yy, zz = np.broadcast_arrays(xi[:,None,None], yi[None,:,None], zi[None,None,:])
points = (xx.ravel(), yy.ravel(), zz.ravel())
ip = LinearNDInterpolator(points, data.ravel())

但是，如果您有一个规则网格，那么使用LinearNDInterpolator很可能不是最佳选择，因为它是为分散数据插值而设计的。它构建了数据点的 Delaunay 三角剖分，但在这种情况下，原始数据已经具有非常规则的结构，可以更有效地利用。

由于您的网格是矩形的，您可以将插值构建为三个一维插值的张量积。Scipy 没有这个内置的（到目前为止），但它很容易做到，看到这个线程： http: //mail.scipy.org/pipermail/scipy-user/2012-June/032314.html （使用例如 interp1d 而不是 pchip 来获得一维插值）

于 2012-08-19T14:48:49.660 回答
0

我不相信有任何方法可以将缺少完整副本的东西传递给 LinearNDInterpolator （因为也没有用于三个维度的常规网格的函数）。所以唯一避免创建完整数组的地方是在创建这个点数组的过程中，我不知道你现在是怎么做的，所以在这方面它可能已经很有效了，但我想它可能不值得避免这。

其他然后 np.mgrid+reshape 可能像这样的东西可能是一个选项（对于 n 维也不难写）：

# Create broadcastest versions of xi, yi and zi
# np.broadcast_arrays does not allocate the full arrays
xi, yi, zi = np.broadcast_arrays(xi[:,None,None], yi[:,None,None], zi[:,None,None])

# then you could use .flat to fill a point array:
points = np.empty((xi.size, 3), dtype=xi.dtype)
points[:,0] = xi.flat
points[:,1] = yi.flat
points[:,2] = zi.flat

与.repeat函数相反，这里创建的临时数组并不比原来xi的等数组大。

于 2012-08-16T16:48:54.603 回答