将来自一个纬度-经度网格的数据插值到另一个网格上？答案

【问题标题】：Interpolating data from one latitude-longitude grid onto a different one?将来自一个纬度-经度网格的数据插值到另一个网格上？
【发布时间】：2021-06-22 04:20:40
【问题描述】：

我有两个位于 lat-lon 网格上的数据数组。第一个，A，具有形状 (89, 180)。第二个，B，形状为 (94, 192)。 A 的纬度按降序排列，从 88. 到 -88。 & 经度从 0. 到 358 升序。 B 的纬度从 88.54199982 到 -88.54199982 降序 & 经度从 0. 到 358.125 升序。

我想将 B 的数据重新网格化/插值到 A 的坐标系上，以便我可以使两个数组的大小相同并计算它们之间的空间相关性。（如果这更容易，我也可以将 A 的数据重新网格化/插值到 B 的坐标系上。）我尝试了 mpl_toolkits.basemap.interp(datain, xin, yin, xout, yout)，但这需要 xout 和 yout 的大小相同。我也尝试过 scipy.interpolate.griddata，但我不知道它是如何工作的，我什至不确定这是否能满足我的需求......

【问题讨论】：

简单的线性插值怎么样？分别插入 x 坐标和 y，它应该可以工作。基本上通过网格差异缩放任何坐标。例如，x_new = x_old*(B_width/A_width)。例如，如果 B_width 为 2*A_width，则 x_new = x_old*2。显然，如果 B_width = A_width 则没有缩放，因为它们的大小完全相同。
现在，一旦您有了中间值（因为它可能不是整数），您需要围绕该值创建数据点的“混合”以获得新的数据值。事实上，逆向工作效果最好。例如，对于 B 中的每个点 (i,j)，数据值为 B[i,j]。为了得到这个值，我们需要在 A 中查找它。但是我们必须将 i,j “缩放”到 A 中，这不是整数点，所以我们不能直接使用它们来访问我们的数组。我们可以使用加权和或其他方法（三次插值）简单地舍入/下限/天花板或在 A 周围的值之间进行插值。
对不起，我不太明白你的第二条评论……你可以用代码形式回答吗？
例如，B[i,j] = A[floor(iA_width/B_width), floor(jA_height/B_height)] 是一个简单的映射。但这是不准确的，因为某些点可能位于 A 中的“网格点”之间，但我们总是将这些点“钳制”到网格点（例如，如果 i = 3.9，我们使用 i = 3，这可能不能很好地表示该值。舍入是最好的，但最好在 3(*0.1) 和 4(*0.9) 之间进行线性插值。因此，根据您的需要，加权方法更好。

标签： python interpolation latitude-longitude

【解决方案1】：

您可能需要查看pyresample 以了解此问题和其他类似的地理插值问题。它提供了多种插值方法，可以很好地处理纬度/经度数据，并包含basemap 支持。我建议使用这个包，因为您还可以创建使用 Proj4 定义定义域的 AreaDefinition 对象，然后将数据注册到 AreaDefinition。

对于你的具体问题，我会做以下（注意，插值步骤不完整，见下文）：

from pyresample.geometry import SwathDefinition
from pyresample.kd_tree import resample_nearest

def interp_b_to_a(a, b):
    '''Take in two dictionaries of arrays and interpolate the second to the first.
    The dictionaries must contain the following keys: "data", "lats", "lons"
    whose values must be numpy arrays.
    '''
    def_a = SwathDefinition(lons=a['lons'], lats=a['lats'])
    def_b = SwathDefinition(lons=b['lons'], lats=b['lats'])

    interp_dat = resample_nearest(def_b, b['data'], def_a, ...)
    new_b = {'data':interp_dat,
             'lats':copy(a['lats']),
             'lons':copy(a['lons'])
            }
    return new_b

请注意，调用resample_nearest 的插值步骤不完整。您还需要指定radius_of_influence，这是在每个点周围使用的搜索半径（以米为单位）。这取决于数据的分辨率。您可能还需要指定 nprocs 以加快处理速度，如果您使用的是屏蔽数据，则可能还需要指定 fill_value。

【讨论】：

def_a 和 def_b 是否意味着具有相同的定义？
不，def_a 将是您要插入的纬度和地段。 def_b 将是您要从中插值的纬度和地段。现在编辑答案。
我开始测试它，但不知道要为 ... 部分添加什么。最后，我没有使用此页面上的任何解决方案，因为我意识到我不需要为我试图解决的问题将纬度-经度网格相互插入。不过感谢您的回答。
仅供参考，... 表示其他关键字，如上一段中所讨论的那样。因此，如果您将来需要它，只需查看 pyresample 的文档并查看其他关键字。

【解决方案2】：

基于@Vorticity，我将其编辑为：

from pyresample.geometry import SwathDefinition
from pyresample.kd_tree import resample_nearest

def interp_b_to_a(a, b):
    def_a = SwathDefinition(lons=a['lons'], lats=a['lats'])
    def_b = SwathDefinition(lons=b['lons'], lats=b['lats'])

    interp_dat = resample_nearest(def_a, a['data'], def_b,radius_of_influence = 5000)
    new_b = {'data':interp_dat,
             'lats':b['lats'],
             'lons':b['lons']
            }
    return new_b

【讨论】：