确定给定点半径内网格值平均值的最快方法

Question

我有一个数值的二维数组。数组中的每个元素表示一个网格点，其中每个长方体的一侧为13km。我需要确定网格上某个特定点50英里范围内所有点的平均值。在

我当前的解决方案确定一个边界框，然后使用它们的索引引用该框中数组中的项，这对于numpy来说很慢。我在想一个更快的解决方案。在

当前解决方案：

num_x = 400        #horizontal dimension of the 2D array
num_y = 300        #vertical dimension of the 2D array
num_dx = 6         #maximum number of horizontal grid points that fit within 50 miles
num_dy = 6         #same as above but for vertical (square grid)
radius_m = 80467.2 #50 miles expressed in meters
values = []        # stores the extracted values

for ix in range(-num_dx,num_dx+1):
    for jy in range(-num_dy,num_dy+1):
        # Determine distance to this point
        dist = ((ix*dx)**2+(jy*dy)**2)**0.5
        if dist <= radius_m:
            # Ensure this grid point actually exists within the grid
            if (j+jy) < num_y and (i+ix) < num_x:
                value = myarray[i+ix,j+jy]
                    if value is not masked and value >= 0:
                        values.append(float(value))

average = sum(values) / float(len(values))

由于访问myarray超过100次以提取单个元素的值，因此速度很慢（大约需要1.5秒）。有没有一种向量方法在这里更有效？我似乎想不出一种方法来使用掩码来实现这一点，因为条件是基于网格点相对于另一个网格点的位置，而不是元素本身的值。在

Answer 1

您的代码不可运行，并且似乎包含i < num_dx或j < num_dy（然后它会绕到数组的另一端）的bug。但是对变量名做一些假设，我会这样做：

# First make sure we stay in the grid
i1, i2 = max(i-num_dx, 0), min(i+num_dx+1, num_x)
j1, j2 = max(j-num_dy, 0), min(j+num_dy+1, num_y)

# Get the radius in blocks, grid should be homogeneous
radius_i = radius_m / 13000.0

# Calc distances per element by broadcasting
DX = np.arange(i1, i2) - i
DY = np.arange(j1, j2)[:, None] - j
mask = DX*DX + DY*DY <= radius_i*radius_i

# Get block of interest and apply mask
values = myarray[i1:i2, j1:j2][mask]

Answer 2

对于内部点（半径不延伸到图像外部），您只需计算一个用于任何内部点的遮罩。从零数组开始：

mask = np.zeros((2 * num_dx + 1, 2 * num_dy + 1), dtype=np.int)

假设您的兴趣点位于该数组的中心，则将半径范围内的每个元素设置为1（此处未显示）。然后

^{pr2}$

然后对于myarray中的任何内部元素(i, j)，您将得到半径范围内的值，如下所示：

values = myarray[i-num_dx: i+num_dx+1, j-num_dy: j+num_dy+1].ravel()[indices]

对于靠近边界的点，在设置mask之前，需要复制{}，并将图像外部的行/列设置为零。在