如何将二维数组从较粗的分辨率插入到较细的分辨率

Question

假设我有一个形状为(21600,43200) 的排放数据， 对应于lat和lon，即

lat = np.arange(21600)*(-0.008333333)+90
lon = np.arange(43200)*0.00833333-180

我还有一个形状为(720,1440,7)的比例因子，对应于lat、lon、day of week和

lat = np.arange(720)*0.25-90 
lon = np.arange(1440)*0.25-180

现在，我想将该因子应用于排放数据，我认为我需要将 (720,1440) 上的因子插入到 (21600,43200)。之后，我可以将插值因子与排放数据相乘，得到新的排放输出。

但我对插值方法有困难。 谁能给我一些建议？

Answer 1

这是您尝试执行的插值类型的完整示例。例如，我使用了emission 形状为(10, 20) 的数据和scale 形状为(5, 10) 的数据。它使用scipy.interpolate.RectBivariateSpline，这是在规则网格上进行插值的推荐方法：

import scipy.interpolate as sci

def latlon(res):
    return (np.arange(res)*(180/res) - 90,
            np.arange(2*res)*(360/(2*res)) - 180)

lat_fine,lon_fine = latlon(10)
emission = np.ones(10*20).reshape(10,20)

lat_coarse,lon_coarse = latlon(5)
scale = np.linspace(0, .5, num=5).reshape(-1, 1) + np.linspace(0, .5, num=10)

f = sci.RectBivariateSpline(lat_coarse, lon_coarse, scale)
scale_interp = f(lat_em, lon_em)

with np.printoptions(precision=1, suppress=True, linewidth=9999):
    print('original emission data:\n%s\n' % emission)
    print('original scale data:\n%s\n' % scale)
    print('interpolated scale data:\n%s\n' % scale_interp)
    print('scaled emission data:\n%s\n' % (emission*scale_interp))

哪个输出：

original emission data:
[[1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1.]
 [1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1.]
 [1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1.]
 [1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1.]
 [1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1.]
 [1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1.]
 [1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1.]
 [1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1.]
 [1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1.]
 [1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1.]]

original scale data:
[[0.  0.1 0.1 0.2 0.2 0.3 0.3 0.4 0.4 0.5]
 [0.1 0.2 0.2 0.3 0.3 0.4 0.5 0.5 0.6 0.6]
 [0.2 0.3 0.4 0.4 0.5 0.5 0.6 0.6 0.7 0.8]
 [0.4 0.4 0.5 0.5 0.6 0.7 0.7 0.8 0.8 0.9]
 [0.5 0.6 0.6 0.7 0.7 0.8 0.8 0.9 0.9 1. ]]

interpolated scale data:
[[0.  0.  0.1 0.1 0.1 0.1 0.2 0.2 0.2 0.2 0.3 0.3 0.3 0.4 0.4 0.4 0.4 0.5 0.5 0.5]
 [0.1 0.1 0.1 0.1 0.2 0.2 0.2 0.3 0.3 0.3 0.3 0.4 0.4 0.4 0.5 0.5 0.5 0.5 0.6 0.6]
 [0.1 0.2 0.2 0.2 0.2 0.3 0.3 0.3 0.3 0.4 0.4 0.4 0.5 0.5 0.5 0.5 0.6 0.6 0.6 0.6]
 [0.2 0.2 0.2 0.3 0.3 0.3 0.4 0.4 0.4 0.4 0.5 0.5 0.5 0.5 0.6 0.6 0.6 0.7 0.7 0.7]
 [0.3 0.3 0.3 0.3 0.4 0.4 0.4 0.4 0.5 0.5 0.5 0.6 0.6 0.6 0.6 0.7 0.7 0.7 0.8 0.8]
 [0.3 0.3 0.4 0.4 0.4 0.5 0.5 0.5 0.5 0.6 0.6 0.6 0.6 0.7 0.7 0.7 0.8 0.8 0.8 0.8]
 [0.4 0.4 0.4 0.5 0.5 0.5 0.5 0.6 0.6 0.6 0.7 0.7 0.7 0.7 0.8 0.8 0.8 0.8 0.9 0.9]
 [0.4 0.5 0.5 0.5 0.5 0.6 0.6 0.6 0.7 0.7 0.7 0.7 0.8 0.8 0.8 0.9 0.9 0.9 0.9 0.9]
 [0.5 0.5 0.6 0.6 0.6 0.6 0.7 0.7 0.7 0.7 0.8 0.8 0.8 0.9 0.9 0.9 0.9 1.  1.  1. ]
 [0.5 0.5 0.6 0.6 0.6 0.6 0.7 0.7 0.7 0.7 0.8 0.8 0.8 0.9 0.9 0.9 0.9 1.  1.  1. ]]

scaled emission data:
[[0.  0.  0.1 0.1 0.1 0.1 0.2 0.2 0.2 0.2 0.3 0.3 0.3 0.4 0.4 0.4 0.4 0.5 0.5 0.5]
 [0.1 0.1 0.1 0.1 0.2 0.2 0.2 0.3 0.3 0.3 0.3 0.4 0.4 0.4 0.5 0.5 0.5 0.5 0.6 0.6]
 [0.1 0.2 0.2 0.2 0.2 0.3 0.3 0.3 0.3 0.4 0.4 0.4 0.5 0.5 0.5 0.5 0.6 0.6 0.6 0.6]
 [0.2 0.2 0.2 0.3 0.3 0.3 0.4 0.4 0.4 0.4 0.5 0.5 0.5 0.5 0.6 0.6 0.6 0.7 0.7 0.7]
 [0.3 0.3 0.3 0.3 0.4 0.4 0.4 0.4 0.5 0.5 0.5 0.6 0.6 0.6 0.6 0.7 0.7 0.7 0.8 0.8]
 [0.3 0.3 0.4 0.4 0.4 0.5 0.5 0.5 0.5 0.6 0.6 0.6 0.6 0.7 0.7 0.7 0.8 0.8 0.8 0.8]
 [0.4 0.4 0.4 0.5 0.5 0.5 0.5 0.6 0.6 0.6 0.7 0.7 0.7 0.7 0.8 0.8 0.8 0.8 0.9 0.9]
 [0.4 0.5 0.5 0.5 0.5 0.6 0.6 0.6 0.7 0.7 0.7 0.7 0.8 0.8 0.8 0.9 0.9 0.9 0.9 0.9]
 [0.5 0.5 0.6 0.6 0.6 0.6 0.7 0.7 0.7 0.7 0.8 0.8 0.8 0.9 0.9 0.9 0.9 1.  1.  1. ]
 [0.5 0.5 0.6 0.6 0.6 0.6 0.7 0.7 0.7 0.7 0.8 0.8 0.8 0.9 0.9 0.9 0.9 1.  1.  1. ]]

注意事项

• scipy.interpolate 中的插值方法期望 x 和 y 都严格增加，因此您必须确保您的 emission 数据排列在一个网格中，这样：

  lat = np.arange(21600)*0.008333333 - 90
  

  代替：

  lat = np.arange(21600)*(-0.008333333) + 90
  

  就像你上面所说的那样。您可以像这样翻转emission 数据：

  emission = emission[::-1, :]
  

Answer 2

如果你只是在寻找最近邻或线性插值，你可以使用 xarray 原生的da.interp 方法：

scaling_interped = scaling_factor.interp(
    lon=emissions.lon,
    lat=emissions.lat,
    method='nearest')  # or 'linear'

请注意，这将显着增加数组的大小。假设这些是 64 位浮点数，结果将大约为 (21600*43200*7)*8/(1024**3) 或 48.7 GB。您可以通过按星期几分块数组并进行计算 out of core with dask 将内存大小减少 7 倍。

如果您想使用最近或线性以外的插值方案，请使用方法suggested by tel。