Как эффективно рассчитать полные (2pi) углы между тремя векторами 2d точек - PullRequest
Новая
       24

Как эффективно рассчитать полные (2pi) углы между тремя векторами 2d точек

0 голосов
/ 22 июня 2019

У меня есть три массива фигуры [n, 2], содержащих список точек. Давайте назовем их а, б и в. Я хочу найти полный угол между ab и bc. Использование acos дает мне только пи радианы, но я хочу полную шкалу 2pi. Я подумал об использовании atan2, но не уверен, как рассчитать векторы y и x, необходимые для atan2 - я пытался использовать векторные нормы, но они по своей сути положительны. Есть ли способ, которым я могу сделать это полностью, используя функции NumPy для эффективности?

Ответы [ 2 ]

2 голосов
/ 22 июня 2019

Использование только метода arccos дает вам только абсолютный угол между векторами, а не по часовой стрелке или против часовой стрелки. Вы можете увеличить это, проверив, отрицательно ли скалярное произведение a против перпендикуляра b, что означает угол против часовой стрелки. 

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].

0 голосов
/ 22 июня 2019

Неудивительно, что здесь можно использовать функцию 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])
...