<p>既然您显然在使用numpy,那么让我们实际使用numpy;它会更快。如果在使用numpy时完全避免python循环,而是使用它的矢量化数组操作,那么它几乎总是更快,而且通常更容易阅读。你知道吗</p>
<pre><code>a, b, c = np.loadtxt('test.dat', dtype='double', unpack=True)
q = 3e5 * c / 100 # why not just 3e3 * c?
x = q * np.sin(a) * np.cos(b)
y = q * np.sin(a) * np.sin(b)
z = q * np.cos(a)
</code></pre>
<p>现在,后面的示例代码并没有完成您可能希望它做的事情—注意您每次都是如何说<code>xdist = ...</code>?您正在覆盖该变量,而没有对其执行任何操作。我假设你想要每对点之间的平方欧几里德距离,使矩阵<code>dists</code>,其中<code>dists[i, j]</code>等于第<code>i</code>点和第<code>j</code>点之间的距离。你知道吗</p>
<p>简单的方法是,如果您有可用的scipy:</p>
<pre><code># stack the points into a num_pts x 3 matrix
pts = np.hstack([thing.reshape((-1, 1)) for thing in (x, y, z)])
# get squared euclidean distances in a matrix
dists = scipy.spatial.squareform(scipy.spatial.pdist(pts, 'sqeuclidean'))
</code></pre>
<p>如果你的列表很大,不使用squareform会更节省内存,但是它是一种压缩格式,要找到特定的距离对有点困难。你知道吗</p>
<p>稍微难一点,如果你不能/不想使用scipy:</p>
<pre><code>pts = np.hstack([thing.reshape((-1, 1)) for thing in (x, y, z)])
sqnorms = np.sum(pts ** 2, axis=1)
dists = sqnorms.reshape((-1, 1)) - 2 * np.dot(pts, pts.T) + sqnorms
</code></pre>
<p>它基本上实现了公式(a-b)^2=a^2-2ab+b^2,但都是向量。你知道吗</p>