找到位于球形区域中的所有点答案

【问题标题】：FInd all the points that lie with in a spherical region找到位于球形区域中的所有点
【发布时间】：2016-10-28 13:57:01
【问题描述】：

例如，找到下面的图片，它解释了一个简单的 2D 案例的问题。每个点的标签 (N) 和坐标 (x,y) 是已知的。我需要找到位于红色圆圈内的所有点标签

我的实际问题是在 3D 中，点分布不均匀

此处附有包含 7.25 M 点坐标的示例输入文件point file。

我尝试了以下代码

import numpy as np

C = [50,50,50]

R = 20 

centroid = np.loadtxt('centroid') #chk the file attached

def dist(x,y): return sum([(xi-yi)**2 for xi, yi in zip(x,y)])

elabels=[i+1 for i in range(len(centroid)) if dist(C,centroid[i])<=R**2]

单次搜索大约需要 10 分钟。有什么建议可以让它更快吗？

谢谢，普里蒂维

【问题讨论】：

Python 是一种脚本语言，这意味着每次启动时每一行都会转换为机器码。一种编程语言，如 C、Java、Pascal 等，在运行之前首先完全编译，运行速度要快得多。因此，我建议您开始使用非脚本编程语言。

标签： python list search

【解决方案1】：

使用numpy 时，避免对数组使用列表推导。

您的计算可以使用像这样的矢量化表达式来完成

centre = np.array((50., 50., 50.))
points = np.loadtxt('data')

distances2= np.sum((points-centre)**2, axis=1)

points 是一个N x 2 数组，points-centre 也是一个N x 2 数组， (points-centre)**2 计算每个差值元素的平方，最终np.sum(..., axis=1) 将沿轴号的平方差值的元素相加。 1，即跨列。

要过滤位置数组，可以使用布尔索引

close = points[distances2<max_dist**2]

【讨论】：

非常感谢！它在 0.4 秒内完成。实际上我需要索引而不是值，我使用 ind = [i+1 for labels in close for i in range(len(close)) if labels == points[i]] 有更好的方法吗？
我不清楚你想做什么......无论如何，索引的向量是表达式np.arange(N)[distances2<max_dist**2]的值

【解决方案2】：

您正在大量调用dist 函数。您可以尝试对其进行低级优化，并使用更有效的 timeit Python 模块进行控制。在我的机器上，我尝试了这个：

def dist(x,y):
    d0 = y[0] -x[0]
    d1 = y[1] -x[1]
    d2 = y[2] -x[2]
    return d0 * d0 + d1*d1 + d2*d2

timeit 说它快了 3 倍以上。

这个就在中间：

def dist(x,y):
    s = 0
    for i in range(len(x)):
        d = y[i] - x[i]
        s += d * d
    return s

【讨论】：