Scipy - 3d griddata - 为什么有必要将griddata xi参数转换为元组？

Question

为什么下面的调用griddata失败？

import scipy.interpolate 
import numpy as np 
grid_vals = np.meshgrid(*([np.linspace(-1,1,200)] * 3)) 
interp_vals = scipy.interpolate.griddata(np.random.randn(50,3), np.random.randn(50), grid_vals, 'linear') 

发生以下例外： ValueError异常：在XI维数不匹配X

如果我投的XI（grid_vals）参数元组：

interp_vals = scipy.interpolate.griddata(np.random.randn(50,3), np.random.randn(50), tuple(grid_vals), 'linear') 

误差变为远。为什么？

Answer 1

的基本原因是griddata通过points = _ndim_coords_from_arrays(points)函数，它的文档读取通过这两个points和xi：在元组

Convert a tuple of coordinate arrays to a (..., ndim)-shaped array. 

和关键行动是：

p = np.broadcast_arrays(*points) 

别的，包括列表，只是转换为阵列：

points = np.asanyarray(points) 

实际插值需要最后一个带有“3d”维的数组。

因此，您的3 (200,200,200)数组列表变为(3,200,200,200)形状的数组。但是你的points数组是(50,3)。来自200的number of dimensions in xi does not match x消息结果不匹配3。

该griddata文件是明确的关于points，较少xi。但是它的示例使用(x, Y)使用来自mgrid的数组。

所以这会工作：

X, Y, Z = np.meshgrid(*([np.linspace(-1,1,200)] * 3)) 
interp_vals = scipy.interpolate.griddata(points, values, (X,Y,Z), 'linear') 

从meshgrid列表生成所需的阵列的另一种方法是使它成为一个阵列，并推出第一维

grid_vals = np.rollaxis(np.array(grid_vals),0,4) 

发生的另一种方式一个网格是np.ix_，它以元组的形式返回一个开放的网格。像这样的开放网格确实需要广播。

单一指向将与进行插值之一：

interpolate.griddata(points,values,[[[[0,0,0]]]],'linear') 
interpolate.griddata(points,values,([0],[0],[0]),'linear') 

见反应约翰4123拉请求有关于个为什么更多的讨论。

Answer 2

它的失败是因为基础SciPy的插值模块的代码能够正确处理多维元组，而不是多维列表：

https://github.com/scipy/scipy/blob/master/scipy/interpolate/interpnd.pyx#L167

我创建issue 4123和pull request 4124来解决这个问题，在至少为列表。