如何有效地计算 2d 点的三个向量之间的完整 (2pi) 角度
How to efficiently calculate full (2pi) angles between three vectors of 2d points
我有三个形状为 [n, 2] 的 numpy 数组,其中包含一个点列表。我们称它们为 a、b 和 c。我想找到 ab 和 bc 之间的完整角度。使用 acos 只能得到 pi 弧度,但我想要完整的 2pi 比例。我考虑过使用 atan2,但不确定如何计算 atan2 所需的 y 和 x 向量 - 我尝试使用向量范数,但它们本质上是正的。有什么办法可以完全使用 numpy 函数来提高效率吗?
单独使用 arccos
方法只能给出向量之间的绝对角度,而不是 clock-wise 或 counter-clockwise。您可以通过检查 a
与 b
的垂线的点积是否为负来增强这一点,表示 counter-clockwise 角度。
import numpy as np
def dot(a, b):
return np.sum(a * b, axis=-1)
def mag(a):
return np.sqrt(np.sum(a*a, axis=-1))
def angle(a, b):
cosab = dot(a, b) / (mag(a) * mag(b)) # cosine of angle between vectors
angle = np.arccos(cosab) # what you currently have (absolute angle)
b_t = b[:,[1,0]] * [1, -1] # perpendicular of b
is_cc = dot(a, b_t) < 0
# invert the angles for counter-clockwise rotations
angle[is_cc] = 2*np.pi - angle[is_cc]
return angle
print(angle(
np.array([[1, 0], [1, 0]]),
np.array([[0, 1], [0, -1]])
))
将打印 [pi/2, 3pi/2]
的浮点值。
此函数输出范围 [0, 2*pi]
。
不出所料,这里可以使用angle
函数。它需要一个复杂的参数 x + y i
:
这种方法的优点是比较容易得到相对角度。使用 atan2
这会有点棘手。
def get_angle(a,b,yx=False):
# make sure inputs are contiguous float
# swap x and if requested
a,b = map(np.ascontiguousarray, (a[...,::-1],b[...,::-1]) if yx else (a,b), (float,float))
# view cast to complex, prune excess dimension
A,B = (z.view(complex).reshape(z.shape[:-1]) for z in (a,b))
# to get the relative angle we must either divide
# or (probably cheaper) multiply with the conjugate
return np.angle(A.conj()*B)
a,b,c = np.random.randn(3,20,2)
# let's look at a roundtrip as a test
get_angle(a,b)+get_angle(b,c)+get_angle(c,a)
# array([ 0.00000000e+00, 1.66533454e-16, 4.44089210e-16, -2.22044605e-16,
# 0.00000000e+00, 0.00000000e+00, 0.00000000e+00, -4.44089210e-16,
# 0.00000000e+00, -1.66533454e-16, 2.22044605e-16, 0.00000000e+00,
# 0.00000000e+00, 2.22044605e-16, 6.28318531e+00, 8.32667268e-17,
# 2.22044605e-16, -6.28318531e+00, -2.22044605e-16, 6.28318531e+00])
# some zeros, some 2pi and some -2pi ==> looks ok
# Let's also check the sum of angles of triangles abc:
get_angle(a-c,b-c)+get_angle(b-a,c-a)+get_angle(c-b,a-b)
# array([-3.14159265, -3.14159265, 3.14159265, -3.14159265, -3.14159265,
# 3.14159265, -3.14159265, -3.14159265, 3.14159265, -3.14159265,
# -3.14159265, 3.14159265, -3.14159265, -3.14159265, 3.14159265,
# 3.14159265, -3.14159265, -3.14159265, 3.14159265, 3.14159265])
我有三个形状为 [n, 2] 的 numpy 数组,其中包含一个点列表。我们称它们为 a、b 和 c。我想找到 ab 和 bc 之间的完整角度。使用 acos 只能得到 pi 弧度,但我想要完整的 2pi 比例。我考虑过使用 atan2,但不确定如何计算 atan2 所需的 y 和 x 向量 - 我尝试使用向量范数,但它们本质上是正的。有什么办法可以完全使用 numpy 函数来提高效率吗?
单独使用 arccos
方法只能给出向量之间的绝对角度,而不是 clock-wise 或 counter-clockwise。您可以通过检查 a
与 b
的垂线的点积是否为负来增强这一点,表示 counter-clockwise 角度。
import numpy as np
def dot(a, b):
return np.sum(a * b, axis=-1)
def mag(a):
return np.sqrt(np.sum(a*a, axis=-1))
def angle(a, b):
cosab = dot(a, b) / (mag(a) * mag(b)) # cosine of angle between vectors
angle = np.arccos(cosab) # what you currently have (absolute angle)
b_t = b[:,[1,0]] * [1, -1] # perpendicular of b
is_cc = dot(a, b_t) < 0
# invert the angles for counter-clockwise rotations
angle[is_cc] = 2*np.pi - angle[is_cc]
return angle
print(angle(
np.array([[1, 0], [1, 0]]),
np.array([[0, 1], [0, -1]])
))
将打印 [pi/2, 3pi/2]
的浮点值。
此函数输出范围 [0, 2*pi]
。
不出所料,这里可以使用angle
函数。它需要一个复杂的参数 x + y i
:
这种方法的优点是比较容易得到相对角度。使用 atan2
这会有点棘手。
def get_angle(a,b,yx=False):
# make sure inputs are contiguous float
# swap x and if requested
a,b = map(np.ascontiguousarray, (a[...,::-1],b[...,::-1]) if yx else (a,b), (float,float))
# view cast to complex, prune excess dimension
A,B = (z.view(complex).reshape(z.shape[:-1]) for z in (a,b))
# to get the relative angle we must either divide
# or (probably cheaper) multiply with the conjugate
return np.angle(A.conj()*B)
a,b,c = np.random.randn(3,20,2)
# let's look at a roundtrip as a test
get_angle(a,b)+get_angle(b,c)+get_angle(c,a)
# array([ 0.00000000e+00, 1.66533454e-16, 4.44089210e-16, -2.22044605e-16,
# 0.00000000e+00, 0.00000000e+00, 0.00000000e+00, -4.44089210e-16,
# 0.00000000e+00, -1.66533454e-16, 2.22044605e-16, 0.00000000e+00,
# 0.00000000e+00, 2.22044605e-16, 6.28318531e+00, 8.32667268e-17,
# 2.22044605e-16, -6.28318531e+00, -2.22044605e-16, 6.28318531e+00])
# some zeros, some 2pi and some -2pi ==> looks ok
# Let's also check the sum of angles of triangles abc:
get_angle(a-c,b-c)+get_angle(b-a,c-a)+get_angle(c-b,a-b)
# array([-3.14159265, -3.14159265, 3.14159265, -3.14159265, -3.14159265,
# 3.14159265, -3.14159265, -3.14159265, 3.14159265, -3.14159265,
# -3.14159265, 3.14159265, -3.14159265, -3.14159265, 3.14159265,
# 3.14159265, -3.14159265, -3.14159265, 3.14159265, 3.14159265])