Kąty między dwoma n-wymiarowymi wektorami w Pythonie
82
Muszę określić kąt (y) między dwoma n-wymiarowymi wektorami w Pythonie. Na przykład dane wejściowe mogą być dwiema listami, takimi jak: [1,2,3,4]i [6,7,8,9].
Najlepszą odpowiedzią jest @ MK83, ponieważ jest to dokładnie wyrażenie matematyczne theta = atan2 (u ^ v, uv). nawet przypadek, w którym u = [0] lub v = [0 0] jest uwzględnione, ponieważ jest to tylko moment, w którym atan2 wytworzy NaN w innych odpowiedziach NaN zostanie wyprodukowane przez / norm (u) lub / norm (v)
PilouPili
Odpowiedzi:
66
import math
defdotproduct(v1, v2):returnsum((a*b) for a, b inzip(v1, v2))
deflength(v):return math.sqrt(dotproduct(v, v))
defangle(v1, v2):return math.acos(dotproduct(v1, v2) / (length(v1) * length(v2)))
Ponadto, jeśli potrzebujesz tylko cos, sin, tan kąta, a nie samego kąta, możesz pominąć matematykę. Acos, aby uzyskać cosinus, i użyć iloczynu krzyżowego, aby uzyskać sinus.
mbeckish
10
Biorąc pod uwagę, że math.sqrt(x)jest to równoważne x**0.5i math.pow(x,y)jest równoważne z x**y, jestem zaskoczony, że przetrwały one topór redundancji używany podczas przejścia z Pythona 2.x-> 3.0. W praktyce zwykle robię tego typu rzeczy numeryczne w ramach większego procesu intensywnie obliczającego, a obsługa interpretera dla „**” przechodzącego bezpośrednio do kodu bajtowego BINARY_POWER, w porównaniu do wyszukiwania „matematyki”, do jego atrybutu „sqrt”, a następnie do boleśnie powolnego kodu bajtowego CALL_FUNCTION, mogą spowodować wymierną poprawę szybkości bez konieczności kodowania lub kosztów czytelności.
PaulMcG
5
Jak w odpowiedzi z numpy: Może się to nie powieść, jeśli pojawi się błąd zaokrągleń! Może się to zdarzyć w przypadku wektorów równoległych i antyrównoległych!
BandGap
2
Uwaga: to się nie powiedzie, jeśli wektory są identyczne (np angle((1., 1., 1.), (1., 1., 1.)).). Zobacz moją odpowiedź dla nieco bardziej poprawnej wersji.
David Wolever
2
Jeśli mówisz o powyższej implementacji, to kończy się ona niepowodzeniem z powodu błędów zaokrąglania, a nie z powodu równoległości wektorów.
Tempo
153
Uwaga : wszystkie inne odpowiedzi tutaj nie powiedzie się, jeśli dwa wektory mają też ten sam kierunek (ex, (1, 0, 0), (1, 0, 0)) lub przeciwnych kierunkach (ex, (-1, 0, 0), (1, 0, 0)).
Oto funkcja, która poprawnie obsłuży takie przypadki:
import numpy as np
defunit_vector(vector):""" Returns the unit vector of the vector. """return vector / np.linalg.norm(vector)
defangle_between(v1, v2):""" Returns the angle in radians between vectors 'v1' and 'v2'::
>>> angle_between((1, 0, 0), (0, 1, 0))
1.5707963267948966
>>> angle_between((1, 0, 0), (1, 0, 0))
0.0
>>> angle_between((1, 0, 0), (-1, 0, 0))
3.141592653589793
"""
v1_u = unit_vector(v1)
v2_u = unit_vector(v2)
return np.arccos(np.clip(np.dot(v1_u, v2_u), -1.0, 1.0))
Czy nie byłoby lepiej użyć np.isnanzamiast tego z biblioteki matematycznej? W teorii powinny być identyczne, ale w praktyce nie mam pewności. Tak czy inaczej, wyobrażam sobie, że byłoby to bezpieczniejsze.
Hooked
2
Mój numpy (wersja == 1.12.1) może używać arccosbezpośrednio i bezpiecznie. : W [140]: np.arccos (np.dot (np.array ([1,0,0]), np.array ([- 1,0,0]))) Na zewnątrz [140]: 3.1415926535897931 W [ 141]: np.arccos (np.dot (np.array ([1,0,0]), np.array ([1,0,0]))) Out [141]: 0.0
ene
2
Specjalny przypadek, w którym co najmniej jeden wektor wejściowy jest wektorem zerowym, jest pomijany, co jest problematyczne przy dzieleniu unit_vector. Jedną z możliwości jest po prostu zwrócenie wektora wejściowego w tej funkcji, gdy tak jest.
kafman
1
angle_between ((0, 0, 0), (0, 1, 0)) da jako wynik nan, a nie 90
FabioSpaghetti
2
@kafman Kąt 0-wektorów jest niezdefiniowany (w matematyce). Więc fakt, że wywołuje błąd, jest dobry.
from numpy import (array, dot, arccos, clip)
from numpy.linalg import norm
u = array([1.,2,3,4])
v = ...
c = dot(u,v)/norm(u)/norm(v) # -> cosine of the angle
angle = arccos(clip(c, -1, 1)) # if you really want the angle
Ostatnia linia może spowodować błąd, jak się dowiedziałem, z powodu błędów zaokrągleń. Tak więc, jeśli dotkniesz (u, u) / norm (u) ** 2, da to 1.0000000002, a arccos nie powiedzie się (również 'działa' dla wektorów antyrównoległych)
BandGap
Testowałem z u = [1,1,1]. u = [1,1,1,1] działa dobrze, ale każdy dodany wymiar zwraca nieco większe lub mniejsze wartości niż 1 ...
BandGap
3
Uwaga: to się nie powiedzie (ustąpi nan), gdy kierunek dwóch wektorów jest identyczny lub przeciwny. Zobacz moją odpowiedź, aby uzyskać bardziej poprawną wersję.
David Wolever
2
dodając do tego komentarz neo, ostatnia linia powinna być taka, angle = arccos(clip(c, -1, 1))aby uniknąć zaokrąglania problemów. To rozwiązuje problem @DavidWolever.
Tim Tisdall
4
Dla osób korzystających z powyższego fragmentu kodu: clipnależy dodać do listy numpy importów.
Liam Deacon
27
Inną możliwością jest użycie tylko numpyi daje to kąt wewnętrzny
import numpy as np
p0 = [3.5, 6.7]
p1 = [7.9, 8.4]
p2 = [10.8, 4.8]
'''
compute angle (in degrees) for p0p1p2 corner
Inputs:
p0,p1,p2 - points in the form of [x,y]
'''
v0 = np.array(p0) - np.array(p1)
v1 = np.array(p2) - np.array(p1)
angle = np.math.atan2(np.linalg.det([v0,v1]),np.dot(v0,v1))
print np.degrees(angle)
a oto wynik:
In [2]: p0, p1, p2 = [3.5, 6.7], [7.9, 8.4], [10.8, 4.8]
In [3]: v0 = np.array(p0) - np.array(p1)
In [4]: v1 = np.array(p2) - np.array(p1)
In [5]: v0
Out[5]: array([-4.4, -1.7])
In [6]: v1
Out[6]: array([ 2.9, -3.6])
In [7]: angle = np.math.atan2(np.linalg.det([v0,v1]),np.dot(v0,v1))
In [8]: angle
Out[8]: 1.8802197318858924
In [9]: np.degrees(angle)
Out[9]: 107.72865519428085
Możesz także określić kąt widzenia, aby obliczyć kąt przez rzutowanie:
vg.angle(vec1, vec2, look=vg.basis.z)
Lub oblicz kąt ze znakiem za pomocą rzutu:
vg.signed_angle(vec1, vec2, look=vg.basis.z)
Bibliotekę stworzyłem podczas mojego ostatniego uruchomienia, gdzie była motywowana takimi zastosowaniami: prostymi pomysłami, które są pełne lub nieprzejrzyste w NumPy.
Jeśli chcesz mieć kąty ze znakiem, musisz określić, czy dana para jest praworęczna czy leworęczna ( więcej informacji znajdziesz na wiki ).
Moje rozwiązanie to:
defunit_vector(vector):""" Returns the unit vector of the vector"""return vector / np.linalg.norm(vector)
defangle(vector1, vector2):""" Returns the angle in radians between given vectors"""
v1_u = unit_vector(vector1)
v2_u = unit_vector(vector2)
minor = np.linalg.det(
np.stack((v1_u[-2:], v2_u[-2:]))
)
if minor == 0:
raise NotImplementedError('Too odd vectors =(')
return np.sign(minor) * np.arccos(np.clip(np.dot(v1_u, v2_u), -1.0, 1.0))
Nie jest z tego powodu idealny, NotImplementedErrorale w moim przypadku działa dobrze. To zachowanie można naprawić (ponieważ ręczność jest określana dla dowolnej pary), ale wymaga to więcej kodu niż chcę i muszę napisać.
Dla nielicznych, którzy mogli (z powodu komplikacji SEO) zakończyć tutaj próbę obliczenia kąta między dwiema liniami w pythonie, podobnie jak w (x0, y0), (x1, y1)liniach geometrycznych, jest poniżej minimalne rozwiązanie (korzysta z shapelymodułu, ale można go łatwo zmodyfikować, aby nie):
from shapely.geometry import LineString
import numpy as np
ninety_degrees_rad = 90.0 * np.pi / 180.0defangle_between(line1, line2):
coords_1 = line1.coords
coords_2 = line2.coords
line1_vertical = (coords_1[1][0] - coords_1[0][0]) == 0.0
line2_vertical = (coords_2[1][0] - coords_2[0][0]) == 0.0# Vertical lines have undefined slope, but we know their angle in rads is = 90° * π/180if line1_vertical and line2_vertical:
# Perpendicular vertical linesreturn0.0if line1_vertical or line2_vertical:
# 90° - angle of non-vertical line
non_vertical_line = line2 if line1_vertical else line1
returnabs((90.0 * np.pi / 180.0) - np.arctan(slope(non_vertical_line)))
m1 = slope(line1)
m2 = slope(line2)
return np.arctan((m1 - m2)/(1 + m1*m2))
defslope(line):# Assignments made purely for readability. One could opt to just one-line return them
x0 = line.coords[0][0]
y0 = line.coords[0][1]
x1 = line.coords[1][0]
y1 = line.coords[1][1]
return (y1 - y0) / (x1 - x0)
Odpowiedzi:
import math def dotproduct(v1, v2): return sum((a*b) for a, b in zip(v1, v2)) def length(v): return math.sqrt(dotproduct(v, v)) def angle(v1, v2): return math.acos(dotproduct(v1, v2) / (length(v1) * length(v2)))
Uwaga : to się nie powiedzie, gdy wektory mają ten sam lub przeciwny kierunek. Prawidłowa implementacja jest tutaj: https://stackoverflow.com/a/13849249/71522
źródło
math.sqrt(x)
jest to równoważnex**0.5
imath.pow(x,y)
jest równoważne zx**y
, jestem zaskoczony, że przetrwały one topór redundancji używany podczas przejścia z Pythona 2.x-> 3.0. W praktyce zwykle robię tego typu rzeczy numeryczne w ramach większego procesu intensywnie obliczającego, a obsługa interpretera dla „**” przechodzącego bezpośrednio do kodu bajtowego BINARY_POWER, w porównaniu do wyszukiwania „matematyki”, do jego atrybutu „sqrt”, a następnie do boleśnie powolnego kodu bajtowego CALL_FUNCTION, mogą spowodować wymierną poprawę szybkości bez konieczności kodowania lub kosztów czytelności.angle((1., 1., 1.), (1., 1., 1.))
.). Zobacz moją odpowiedź dla nieco bardziej poprawnej wersji.Uwaga : wszystkie inne odpowiedzi tutaj nie powiedzie się, jeśli dwa wektory mają też ten sam kierunek (ex,
(1, 0, 0)
,(1, 0, 0)
) lub przeciwnych kierunkach (ex,(-1, 0, 0)
,(1, 0, 0)
).Oto funkcja, która poprawnie obsłuży takie przypadki:
import numpy as np def unit_vector(vector): """ Returns the unit vector of the vector. """ return vector / np.linalg.norm(vector) def angle_between(v1, v2): """ Returns the angle in radians between vectors 'v1' and 'v2':: >>> angle_between((1, 0, 0), (0, 1, 0)) 1.5707963267948966 >>> angle_between((1, 0, 0), (1, 0, 0)) 0.0 >>> angle_between((1, 0, 0), (-1, 0, 0)) 3.141592653589793 """ v1_u = unit_vector(v1) v2_u = unit_vector(v2) return np.arccos(np.clip(np.dot(v1_u, v2_u), -1.0, 1.0))
źródło
np.isnan
zamiast tego z biblioteki matematycznej? W teorii powinny być identyczne, ale w praktyce nie mam pewności. Tak czy inaczej, wyobrażam sobie, że byłoby to bezpieczniejsze.arccos
bezpośrednio i bezpiecznie. : W [140]: np.arccos (np.dot (np.array ([1,0,0]), np.array ([- 1,0,0]))) Na zewnątrz [140]: 3.1415926535897931 W [ 141]: np.arccos (np.dot (np.array ([1,0,0]), np.array ([1,0,0]))) Out [141]: 0.0unit_vector
. Jedną z możliwości jest po prostu zwrócenie wektora wejściowego w tej funkcji, gdy tak jest.Używając numpy (wysoce zalecane), zrobiłbyś:
from numpy import (array, dot, arccos, clip) from numpy.linalg import norm u = array([1.,2,3,4]) v = ... c = dot(u,v)/norm(u)/norm(v) # -> cosine of the angle angle = arccos(clip(c, -1, 1)) # if you really want the angle
źródło
nan
), gdy kierunek dwóch wektorów jest identyczny lub przeciwny. Zobacz moją odpowiedź, aby uzyskać bardziej poprawną wersję.angle = arccos(clip(c, -1, 1))
aby uniknąć zaokrąglania problemów. To rozwiązuje problem @DavidWolever.clip
należy dodać do listy numpy importów.Inną możliwością jest użycie tylko
numpy
i daje to kąt wewnętrznyimport numpy as np p0 = [3.5, 6.7] p1 = [7.9, 8.4] p2 = [10.8, 4.8] ''' compute angle (in degrees) for p0p1p2 corner Inputs: p0,p1,p2 - points in the form of [x,y] ''' v0 = np.array(p0) - np.array(p1) v1 = np.array(p2) - np.array(p1) angle = np.math.atan2(np.linalg.det([v0,v1]),np.dot(v0,v1)) print np.degrees(angle)
a oto wynik:
In [2]: p0, p1, p2 = [3.5, 6.7], [7.9, 8.4], [10.8, 4.8] In [3]: v0 = np.array(p0) - np.array(p1) In [4]: v1 = np.array(p2) - np.array(p1) In [5]: v0 Out[5]: array([-4.4, -1.7]) In [6]: v1 Out[6]: array([ 2.9, -3.6]) In [7]: angle = np.math.atan2(np.linalg.det([v0,v1]),np.dot(v0,v1)) In [8]: angle Out[8]: 1.8802197318858924 In [9]: np.degrees(angle) Out[9]: 107.72865519428085
źródło
Jeśli pracujesz z wektorami 3D, możesz to zrobić zwięźle za pomocą paska narzędzi vg . To lekka warstwa na wierzchu numpy.
import numpy as np import vg vec1 = np.array([1, 2, 3]) vec2 = np.array([7, 8, 9]) vg.angle(vec1, vec2)
Możesz także określić kąt widzenia, aby obliczyć kąt przez rzutowanie:
Lub oblicz kąt ze znakiem za pomocą rzutu:
Bibliotekę stworzyłem podczas mojego ostatniego uruchomienia, gdzie była motywowana takimi zastosowaniami: prostymi pomysłami, które są pełne lub nieprzejrzyste w NumPy.
źródło
Rozwiązanie Davida Wolevera jest dobre, ale
Jeśli chcesz mieć kąty ze znakiem, musisz określić, czy dana para jest praworęczna czy leworęczna ( więcej informacji znajdziesz na wiki ).
Moje rozwiązanie to:
def unit_vector(vector): """ Returns the unit vector of the vector""" return vector / np.linalg.norm(vector) def angle(vector1, vector2): """ Returns the angle in radians between given vectors""" v1_u = unit_vector(vector1) v2_u = unit_vector(vector2) minor = np.linalg.det( np.stack((v1_u[-2:], v2_u[-2:])) ) if minor == 0: raise NotImplementedError('Too odd vectors =(') return np.sign(minor) * np.arccos(np.clip(np.dot(v1_u, v2_u), -1.0, 1.0))
Nie jest z tego powodu idealny,
NotImplementedError
ale w moim przypadku działa dobrze. To zachowanie można naprawić (ponieważ ręczność jest określana dla dowolnej pary), ale wymaga to więcej kodu niż chcę i muszę napisać.źródło
Łatwy sposób na znalezienie kąta między dwoma wektorami (działa dla wektora n-wymiarowego),
Kod Pythona:
import numpy as np vector1 = [1,0,0] vector2 = [0,1,0] unit_vector1 = vector1 / np.linalg.norm(vector1) unit_vector2 = vector2 / np.linalg.norm(vector2) dot_product = np.dot(unit_vector1, unit_vector2) angle = np.arccos(dot_product) #angle in radian
źródło
Opierając się na świetnej odpowiedzi sgt pepper i dodając obsługę wyrównanych wektorów oraz dodając ponad 2x przyspieszenie za pomocą Numba
@njit(cache=True, nogil=True) def angle(vector1, vector2): """ Returns the angle in radians between given vectors""" v1_u = unit_vector(vector1) v2_u = unit_vector(vector2) minor = np.linalg.det( np.stack((v1_u[-2:], v2_u[-2:])) ) if minor == 0: sign = 1 else: sign = -np.sign(minor) dot_p = np.dot(v1_u, v2_u) dot_p = min(max(dot_p, -1.0), 1.0) return sign * np.arccos(dot_p) @njit(cache=True, nogil=True) def unit_vector(vector): """ Returns the unit vector of the vector. """ return vector / np.linalg.norm(vector) def test_angle(): def npf(x): return np.array(x, dtype=float) assert np.isclose(angle(npf((1, 1)), npf((1, 0))), pi / 4) assert np.isclose(angle(npf((1, 0)), npf((1, 1))), -pi / 4) assert np.isclose(angle(npf((0, 1)), npf((1, 0))), pi / 2) assert np.isclose(angle(npf((1, 0)), npf((0, 1))), -pi / 2) assert np.isclose(angle(npf((1, 0)), npf((1, 0))), 0) assert np.isclose(angle(npf((1, 0)), npf((-1, 0))), pi)
%%timeit
wyniki bez NumbaI z
źródło
Używanie numpy i dbanie o błędy zaokrąglania BandGap:
from numpy.linalg import norm from numpy import dot import math def angle_between(a,b): arccosInput = dot(a,b)/norm(a)/norm(b) arccosInput = 1.0 if arccosInput > 1.0 else arccosInput arccosInput = -1.0 if arccosInput < -1.0 else arccosInput return math.acos(arccosInput)
Uwaga, ta funkcja zgłosi wyjątek, jeśli jeden z wektorów ma zerową wielkość (podzielenie przez 0).
źródło
Dla nielicznych, którzy mogli (z powodu komplikacji SEO) zakończyć tutaj próbę obliczenia kąta między dwiema liniami w pythonie, podobnie jak w
(x0, y0), (x1, y1)
liniach geometrycznych, jest poniżej minimalne rozwiązanie (korzysta zshapely
modułu, ale można go łatwo zmodyfikować, aby nie):from shapely.geometry import LineString import numpy as np ninety_degrees_rad = 90.0 * np.pi / 180.0 def angle_between(line1, line2): coords_1 = line1.coords coords_2 = line2.coords line1_vertical = (coords_1[1][0] - coords_1[0][0]) == 0.0 line2_vertical = (coords_2[1][0] - coords_2[0][0]) == 0.0 # Vertical lines have undefined slope, but we know their angle in rads is = 90° * π/180 if line1_vertical and line2_vertical: # Perpendicular vertical lines return 0.0 if line1_vertical or line2_vertical: # 90° - angle of non-vertical line non_vertical_line = line2 if line1_vertical else line1 return abs((90.0 * np.pi / 180.0) - np.arctan(slope(non_vertical_line))) m1 = slope(line1) m2 = slope(line2) return np.arctan((m1 - m2)/(1 + m1*m2)) def slope(line): # Assignments made purely for readability. One could opt to just one-line return them x0 = line.coords[0][0] y0 = line.coords[0][1] x1 = line.coords[1][0] y1 = line.coords[1][1] return (y1 - y0) / (x1 - x0)
I pożytek będzie
>>> line1 = LineString([(0, 0), (0, 1)]) # vertical >>> line2 = LineString([(0, 0), (1, 0)]) # horizontal >>> angle_between(line1, line2) 1.5707963267948966 >>> np.degrees(angle_between(line1, line2)) 90.0
źródło