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

我有一个numpy的二维数组值。数组中的每个元素表示一个网格中的一个网格点，每个网格的边长为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))

这是缓慢（约需1.5秒）由于访问myarray中超过100倍，以提取一个单一的元素的值。有没有一种矢量方法可以在这里更好地工作？我似乎无法找出一种方法来使用掩码，因为条件是基于网格点相对于另一个的位置，而不是元素本身的值。

来源

2014-09-05 TheOx

你为什么不切片数组？ – 2014-09-05 17:21:27

我不确定你的意思。你能提供一个例子吗？ – TheOx 2014-09-05 18:40:55

这样想：对于圆上的每个点，还有3个其他点的大小完全相同，但组件上的符号不同。 – 2014-09-05 18:48:45

您的代码不可运行，并且似乎包含i < num_dx或j < num_dy（然后它绕回到数组的另一侧）的错误。但对你的变量名做一些假设，这是我该怎么做的：

# 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]

来源

2014-09-05 19:32:45

对于内部点（半径不会超出图像范围），您可以计算一个用于任何内部点的单个蒙版。开始用零的数组：

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

假设你的兴趣点是在该阵列的中心，设置的半径内下降到1（这里未示出）的每个元素。然后，

indices = np.argwhere(mask.ravel() == 1)

那么对于任何内部元素(i, j)在myarray，你会喜欢的半径范围内获取值：

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

对于边界附近的点，你会作出的mask和一套复印件在设置indices之前，将图像外部的行/列置零。

来源

2014-09-05 18:58:42 bogatron

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

回答

相关问题