Python 按特定方程列出交叉识别（重叠？）

Question

我在像

这样的 2 轴数据上遇到交叉识别问题

A = array([['x0', 'y0', 'data0', 'data0'],
           ['x0', 'y0', 'data0', 'data0'],
           ['x0', 'y0', 'data0', 'data0']])

B = array([['x1', 'y1', 'data1', 'data1'],
           ['x1', 'y1', 'data1', 'data1'],
           ['x1', 'y1', 'data1', 'data1']])

我需要的是找到具有相同 position 的 2 个列表的行。 position 需要描述，因为它们的 distance 足够接近，即：

distance = acos(cos(y0)*cos(y1)*cos(x0-x1)+sin(y0)*sin(y1))
if(distance < 0.001):
    position = True

目前，我使用如下代码：

from math import *
def distance(x1,y1,x2,y2):
    a = acos(cos(y1)*cos(y2)*cos(x1-x2)+sin(y1)*sin(y2))
    if(a < 0.001):
        return True
    else:
        return False
f = open('cross-identification')
for i in range(len(A[0])):
    for j in range(len(B[0])):
        if(distance(A[0][i],A[1][i],B[0][j],B[1][j])==True):
            print(A[0][i],A[1][i],A[2][i],B[2][j],A[3][i],B[3][j],file=f)
        else:continue

几行还可以，但问题是我有海量数据，速度极慢。有什么方法可以让它更快吗？

顺便说一句，我已阅读 this，接近我想要的内容，但我无法更改它。也许我可以得到你的帮助？

Answer 1

这是一个性能问题。你有很多唾手可得的成果，可以极大地改善你的距离功能。

而不是 acos(x) < eps 使用 x > 1 - (eps**2)/2 可以节省你一个膨胀的 acos() （这有效是因为泰勒展开），显然预先计算 1 - (eps**2)/2
而不是计算 cos(y0) 等 4*len(A[0])*len(B[0]) 次，预先计算它们只需要 2*(len(A[0])+len(B[0])) 次，如果 A*B变大
在我看来预计算 cos(x1-x2) 也是一个胜利，但我不是 100% 确定

这只是稍微优化了每个距离评估。我认为更大的问题是您的双循环正在进行 N*M 评估。理想情况下，您可以使用一些智能几何算法来减少该循环。（例如，参见 Mark de Berg 等人的 Computational Geometry，Springer 2008）

但是你的距离（度量）函数太奇葩了！我看不出如何使用该度量函数将您的 (x,y) 转换为 space。我很好奇那个距离函数是从哪里来的。

Answer 2

这是一个经典的外积类型问题：比较两个向量的每个组合。

前面的回答是正确的，你的双循环正在降低效率。一个可能的答案是使用 numpy

根据输入的大小，您可以使用 meshgrid approach or (more complicated) the numpy universal functions with outer product。

对于前者，您可以尝试类似的方法：

import numpy as np
A = [[1, 1],
     [1, 3],
     [3, 1],
     [3, 3]]

B = [[1, 1],
     [3, 1],
     [1, 2],
     [1, 3],
     [3, 3]]

a_array = np.array(A)
b_array = np.array(B)

x0, x1 = np.meshgrid(a_array[:, 0], b_array[:, 0])
y0, y1 = np.meshgrid(a_array[:, 1], b_array[:, 1])

distance = np.arccos(np.cos(y0)*np.cos(y1)*np.cos(x0-x1)+np.sin(y0)*np.sin(y1))
good_pairs = np.nonzero(distance < 0.001)

print(distance)
print(good_pairs)
for b_ind, a_ind in zip(*good_pairs):
    print('Entry {} of A ({}) matches entry {} of B ({})'.format(a_ind, A[a_ind], b_ind, B[b_ind]))

Answer 3

为了不仅避免昂贵的 Haversine 公式，而且为了打开使用 KDTrees 的选项，我建议转换为欧氏坐标和距离。

def to_eucl_coords(lat, lon):
    z = np.sin(lat)
    x = np.sin(lon)*np.cos(lat)
    y = np.cos(lon)*np.cos(lat)
    return x, y, z

def to_eucl_dist(sphdist):
    return 2*np.arcsin(sphdist/2)

KDTrees 易于使用，这里是一个框架，可以帮助您入门。

from scipy.spatial import cKDTree as KDTree

eucl_1 = np.c_[to_eucl_coords(lat1, lon1)]
eucl_2 = np.c_[to_eucl_coords(lat2, lon2)]
t1, t2 = KDTree(eucl_1), KDTree(eucl2)
neighbors = t1.query_ball_tree(t2, threshold)

Python 按特定方程列出交叉识别（重叠？）

Python list cross-identification (overlap?) by specific equation

python

numpy

match

overlap