5

我正在将同事 IDL 代码重写为 python,并提出了一些我感到困惑的差异。根据我发现的其他 SO 问题和邮件列表线程,如果您使用scipy.ndimage.interpolation.map_coordinates并指定order=1它应该进行双线性插值。在比较 IDL 代码(在 GDL 中运行)和 python(map_coordinates)之间的结果时,我得到了不同的结果。然后我尝试使用mpl_toolkits.basemap.interp并得到与 IDL 代码相同的结果。下面是一个简单的例子,显示出了什么问题。有人可以帮我弄清楚我做错了什么map_coordinates或者order=1不是双线性的吗?

from scipy.ndimage.interpolation import map_coordinates
from mpl_toolkits.basemap import interp
import numpy

in_data = numpy.array([[ 25.89125824,  25.88840675],[ 25.90930748,  25.90640068]], dtype=numpy.float32)

map_coordinates(in_data, [[0.0],[0.125]], order=1, mode='nearest')
# map_coordinates results in "25.89090157"
interp(in_data, numpy.array([0,1]), numpy.array([0,1]), numpy.array([0.0]), numpy.array([0.125]), order=1)
# interp results in "25.89351439", same as GDL's "25.8935" when printed

我很好用interp,但我很好奇为什么map_coordinates没有返回相同的结果。我注意到map_coordinates文档没有提到双线性,它实际上是双线性的吗?我错过了什么?

4

1 回答 1

7

使用时map_coordinates,需要转置数组或将坐标更改为 (y, x) 格式,因为数组的形状是(height, width)

from scipy.ndimage.interpolation import map_coordinates
from mpl_toolkits.basemap import interp
import numpy

in_data = numpy.array([[ 25.89125824,  25.88840675],[ 25.90930748,  25.90640068]], dtype=numpy.float32)

print map_coordinates(in_data.T, [[0.0],[0.125]], order=1, mode='nearest')
print interp(in_data, numpy.array([0,1]), numpy.array([0,1]), numpy.array([0.0]), numpy.array([0.125]), order=1)

这将输出:

[ 25.89351463]
[ 25.89351439]
于 2013-02-25T01:09:04.757 回答