我有一个相当大的1000 x 4000像素xr.DataArray从OpenDataCube查询返回，并有一个大集（> 200,000）的xy点值。我需要对数组进行采样以在每个xy点下返回一个值，并返回内插值（例如，如果该点在a 0和1.0像素之间着陆，则返回的值应为0.5）。

xr.interp让我轻松地对插值进行采样，但它返回所有x和y值的每个组合的巨大矩阵，而不仅仅是每个xy点本身的值。我已经尝试过np.diagonal只提取xy点值，但这很慢，很快会遇到内存问题，并且感觉效率低下，因为我仍然需要等待每个值的组合进行内插xr.interp。

可重复的例子

（仅使用10,000个样本点（理想情况下，我需要可以扩展到> 200,000或更多的东西）：

# Create sample array

width, height = 1000, 4000

val_array = xr.DataArray(data=np.random.randint(0, 10, size=(height, width)).astype(np.float32),

coords={'x': np.linspace(3000, 5000, width),

'y': np.linspace(-3000, -5000, height)}, dims=['y', 'x'])

# Create sample points

n = 10000

x_points = np.random.randint(3000, 5000, size=n)

y_points = np.random.randint(-5000, -3000, size=n)

目前的做法

%%timeit

# ATTEMPT 1

np.diagonal(val_array.interp(x=x_points, y=y_points).squeeze().values)

32.6 s ± 1.01 s per loop (mean ± std. dev. of 7 runs, 1 loop each)

有没有人知道更快或更高效的内存方式来实现这一目标？

解决办法：要避免完整网格，您需要引入新维度。

x = xr.DataArray(x_points, dims='z')

y = xr.DataArray(y_points, dims='z')

val_array.interp(x=x, y=y)

将沿着新的z维度为您提供一个数组：

<xarray.DataArray (z: 10000)>

array([4.368132, 2.139781, 5.693636, ..., 3.7505 , 3.713589, 2.28494 ])

Coordinates:

x (z) int64 4647 4471 4692 3942 3468 ... 3040 3993 3027 4427 3749

y (z) int64 -3744 -4074 -3634 -3289 -3221 ... -4195 -4131 -4814 -3362

Dimensions without coordinates: z

36.9 ms ± 1.25 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)