Python将地理半径内的纬度/经度点求和并求和到网格 [英] Python sum lat/lon points within geographic radius and sum to grid

问题描述

基本上，我的实验程序正在尝试查找在给定时间落在有效点半径(例如50公里)内的点数.我的数据在三个单独的数组中进行结构化(但是如果需要，我可以进行重组)，例如:

Basically, my experimental program is trying to find the number of points that fall within a (e.g., 50km) radius of a valid point at a given time. My data is structured (but I can restructure if need-be) in three separate arrays such:

1_LAT，1_LON，1_TIM

1_LAT,1_LON,1_TIM

其中1_LAT，1_LON，1_TIM分别包含大约250个值，分别对应于纬度，经度(十进制度)和时间.

Where 1_LAT,1_LON,1_TIM all contain roughly ~250 values corresponding to Latitude, Longitude (decimal degrees), and time respectively.

我有20组这些数组(即1_LAT，1_LON，1_TIM ... 20_LAT，20_LON，20_TIM).

I have 20 sets of these arrays (i.e., 1_LAT,1_LON,1_TIM...20_LAT,20_LON,20_TIM).

这是我想完成的事情:

1)找出落入每组特定半径内的经/纬度套数.例如，从其他19组点开始，在1_TIM的有效时间有多少点落在1_LAT，1_LON的50 km半径内.然后，我想遍历每个有效时间，以找出每个特定点和有效时间的有效半径中的点数.

1) Figure out the number of lat/lon sets that fall within a particular radius of each set. For example, how many points fall within a 50km radius of 1_LAT,1_LON at the valid time of 1_TIM from the other 19 sets of points. I would then like to iterate through each valid time to figure out the number of points in the valid radius at each specific point and valid time.

我在下面附上了一张图片，以帮助形象地描述.

I have attached a picture below to help visually describe.

黑色正方形代表LAT_1/LON_1阵列中的所有点. 蓝色方块代表LAT_n/LAT_n阵列中的所有点.

The black squares represent all the points in the LAT_1/LON_1 arrays. The blue squares represent all the points in the LAT_n/LAT_n arrays.

我想为每组经纬度数组在每个有效时间计算每个半径中的点数.最终的显示将是地理底图图像上每个网格点的密度的栅格或网格的总和(即，计数数/20).

I would like to count the number of points in each radius at each valid time for each set of lat/lon arrays. The final display would be a summed raster or meshgrid of the denisty (i.e., number of counts / 20) for each grid spot on a geographic basemap image.

我觉得KDEtree可能是实现此目标的最佳方法，但是我对此几乎没有经验.任何想法或建议，将不胜感激.

I have a feeling that a KDEtree may be the best way to accomplish this, but I have little/no experience with such. Any ideas or suggestions would be greatly appreciated.

推荐答案

您将执行以下操作...首先，将每个组的(x, y)坐标分组在单个points_x数组中:

You would do something like the following... First, group your (x, y) coordinates for each group in a single points_x array:

points_1 = np.column_stack((LAT_1, LON_1))
...
points_n = np.column_stack((LAT_n, LON_n))

将它们存储在数组列表中可能是一个好主意:

It may be a good idea to store them in a list of arrays:

points = [point_1, points_2, ..., points_n]

现在，从每组点中制作一个kdTree:

Now, make a kdTree out of each set of points:

import scipy.spatial as spsp
kdtrees = [spsp.cKdTree(p) for p in point]

您已经准备好出发了.如果现在运行以下代码:

And you are ready to go. If you now run the following code:

r = whatever_your_threshold_value_is
points_within_r = np.zeros((len(kdtrees), len(kdtrees)), dtype=np.int)
for j in xrange(len(kdtrees)):
    for k in xrange(j+1, len(kdtrees)):
        points_within_r[j, k] = kdtrees[j].count_neighbors(kdtrees[k], r, 2)
points_within_r = points_within_r + points_within_r.T

现在您应该发现points_within_r[j, k]容纳了points_j中点的半径在points_k中点的半径r内的点.

You should now find that points_within_r[j, k] holds how many points in points_j are within radius r of a point in points_k.

请记住，这里的距离是坐标的欧几里得距离，而忽略了它们所测量的是球面角的事实.

Keep in mind that distances here are the euclidean distance of the coordinates, disregarding the fact that what they measure are spherical angles.

这篇关于Python将地理半径内的纬度/经度点求和并求和到网格的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！