2010-08-15 9 views
14

Podane punkty ABC, jak mogę znaleźć kąt ABC? Przygotowuję narzędzie do pobierania zleceń dla aplikacji do rysowania wektorowego i aby zminimalizować liczbę punktów, które generuje, nie dodam punktów, chyba że kąt pozycji myszy i ostatnie 2 punkty jest większy niż pewien próg. DziękiKąt między 3 punktami?

co miałem:

int CGlEngineFunctions::GetAngleABC(POINTFLOAT a, POINTFLOAT b, POINTFLOAT c) 
{ 
    POINTFLOAT ab; 
    POINTFLOAT ac; 

    ab.x = b.x - a.x; 
    ab.y = b.y - a.y; 

    ac.x = b.x - c.x; 
    ac.y = b.y - c.y; 

    float dotabac = (ab.x * ab.y + ac.x * ac.y); 
    float lenab = sqrt(ab.x * ab.x + ab.y * ab.y); 
    float lenac = sqrt(ac.x * ac.x + ac.y * ac.y); 

    float dacos = dotabac/lenab/lenac; 

    float rslt = acos(dacos); 
    float rs = (rslt * 180)/3.141592; 
    RoundNumber(rs); 
    return (int)rs; 


} 
+2

robię całkiem dobrze, masz algorthm dla niego, ale to nie robi trick. – jmasterx

+5

@abelenky: a to sprawia, że ​​pytanie jest "niejasne lub nieprzydatne" jak dokładnie? Być może źle zrozumiałeś cel rep. To nie jest tam, aby pozwolić ci karać ludzi za próbowanie zrobienia czegoś, co jest dla nich nowe. – jalf

Odpowiedz

25

Pierwsze sugestie dotyczące wybranej metody:

Co nazywasz ac faktycznie cb. Ale jest ok, to jest to, co naprawdę potrzebne. Dalej,

float dotabac = (ab.x * ab.y + ac.x * ac.y); 

To jest twój pierwszy błąd. realny produkt skalarny dwóch wektorów jest:

float dotabac = (ab.x * ac.x + ab.y * ac.y); 

Teraz

float rslt = acos(dacos); 

Tutaj należy zauważyć, że z powodu jakiegoś precyzyjnego straty podczas obliczeń jest to teoretycznie możliwe, że dacos będzie większy niż 1 (lub mniejszy niż -1). Dlatego należy to wyraźnie sprawdzić.

Plus uwaga o wydajności: wywołujecie ciężką funkcję sqrt dwukrotnie do obliczenia długości dwóch wektorów. Następnie dzielisz produkt kropkowany przez te długości. Zamiast tego można nazwać sqrt na mnożenie kwadratów o długości obu wektorów.

Na koniec należy zauważyć, że wynik jest dokładny aż do sign. Oznacza to, że twoja metoda nie rozróżnia 20 ° i -20 °, ponieważ cosinus obu jest taki sam. Twoja metoda daje ten sam kąt dla ABC i CBA.

jedna poprawna metoda obliczania kąta jest jako "oslvbo" proponuje:

float angba = atan2(ab.y, ab.x); 
float angbc = atan2(cb.y, cb.x); 
float rslt = angba - angbc; 
float rs = (rslt * 180)/3.141592; 

(Właśnie zastąpiony atan przez atan2).

To najprostsza metoda, która zawsze daje prawidłowy wynik. Wadą tej metody jest dwukrotne wywołanie ciężkiej funkcji trygonometrii atan2.

Proponuję następującą metodę. Jest nieco bardziej złożony (wymaga zrozumienia pewnych umiejętności trygonometrycznych), jednak jest lepszy od punktu widzenia wydajności. Po prostu wywołuje funkcję trygonometrii atan2. I bez obliczeń pierwiastków kwadratowych.

int CGlEngineFunctions::GetAngleABC(POINTFLOAT a, POINTFLOAT b, POINTFLOAT c) 
{ 
    POINTFLOAT ab = { b.x - a.x, b.y - a.y }; 
    POINTFLOAT cb = { b.x - c.x, b.y - c.y }; 

    // dot product 
    float dot = (ab.x * cb.x + ab.y * cb.y); 

    // length square of both vectors 
    float abSqr = ab.x * ab.x + ab.y * ab.y; 
    float cbSqr = cb.x * cb.x + cb.y * cb.y; 

    // square of cosine of the needed angle  
    float cosSqr = dot * dot/abSqr/cbSqr; 

    // this is a known trigonometric equality: 
    // cos(alpha * 2) = [ cos(alpha) ]^2 * 2 - 1 
    float cos2 = 2 * cosSqr - 1; 

    // Here's the only invocation of the heavy function. 
    // It's a good idea to check explicitly if cos2 is within [-1 .. 1] range 

    const float pi = 3.141592f; 

    float alpha2 = 
     (cos2 <= -1) ? pi : 
     (cos2 >= 1) ? 0 : 
     acosf(cos2); 

    float rslt = alpha2/2; 

    float rs = rslt * 180./pi; 


    // Now revolve the ambiguities. 
    // 1. If dot product of two vectors is negative - the angle is definitely 
    // above 90 degrees. Still we have no information regarding the sign of the angle. 

    // NOTE: This ambiguity is the consequence of our method: calculating the cosine 
    // of the double angle. This allows us to get rid of calling sqrt. 

    if (dot < 0) 
     rs = 180 - rs; 

    // 2. Determine the sign. For this we'll use the Determinant of two vectors. 

    float det = (ab.x * cb.y - ab.y * cb.y); 
    if (det < 0) 
     rs = -rs; 

    return (int) floor(rs + 0.5); 


} 

EDIT:

Ostatnio pracuję na pokrewny temat. I wtedy zrozumiałem, że jest lepszy sposób. W rzeczywistości jest mniej więcej tak samo (za kulisami). Jednak jest to bardziej proste IMHO.

Chodzi o to, aby obrócić oba wektory tak, aby pierwszy był ustawiony na (dodatni) kierunek X. Oczywiście obracanie obu wektorów nie wpływa na kąt między nimi. OTOH po takim obrocie musi właśnie znaleźć kąt drugiego wektora względem osi X. Właśnie do tego służy atan2.

Obrót uzyskuje się poprzez pomnożenie wektora przez następującą macierz:

  • ax, ay
  • -ay AX

Po może zobaczyć, że wektor a pomnożonej przez takiej macierzy rzeczywiście obraca się w kierunku dodatniej osi X.

Uwaga: Ściśle mówiąc powyższa macierz nie tylko się obraca, ale także skaluje. Ale w naszym przypadku jest to ok, ponieważ jedyną rzeczą, która ma znaczenie, jest kierunek wektora, a nie jego długość.

Obrócona wektor b postać:

  • a.x * B.x + a.y * b.y = kropka b
  • -a.y * B.x + a.x * b *.Y = przekroju b

Wreszcie, że odpowiedź może być wyrażona jako:

int CGlEngineFunctions::GetAngleABC(POINTFLOAT a, POINTFLOAT b, POINTFLOAT c) 
{ 
    POINTFLOAT ab = { b.x - a.x, b.y - a.y }; 
    POINTFLOAT cb = { b.x - c.x, b.y - c.y }; 

    float dot = (ab.x * cb.x + ab.y * cb.y); // dot product 
    float cross = (ab.x * cb.y - ab.y * cb.x); // cross product 

    float alpha = atan2(cross, dot); 

    return (int) floor(alpha * 180./pi + 0.5); 
} 
+0

Ładne rozwiązanie! Martwię się, jak sobie poradzić z problemem kąta znaku, dopóki nie przeczytam twojej odpowiedzi. –

+0

Ostatnia funkcja (która ma 'return (int) floor (alpha * 180./pi + 0.5);') jest dobra i daje inną odpowiedź do abc i cba. Działa dobrze! –

+0

Jedyne rozwiązanie, które nie dawało mi precyzji. Dziękuję Ci! – Tiago

1

Off temacie? Ale możesz to zrobić z prawem cosinusów:

Znajdź odległość między A i B (wywołaj to x), a odległość między B i C (wywołaj to y), a odległością między A i C (wywołanie to) ten z).

Następnie, że z^2^2 = X + Y^2-2 * x Y cos (kąt potrzeby)

zatem, że kąt^cos -1 ((Z^2 X^2-y^2)/(2xy)) = kąt

+0

Utraciłeś ujemny we frakcji. Powinno to być 'x^2 + y^2 - z^2' jak w mojej odpowiedzi. –

4

ARccOS p = ((a + c^2^2 - b^2)/2Ac) można

gdzie a jest przeciwległa kąt α, b jest przeciwstawnym kątem β, a c jest przeciwnym kątem γ. Tak więc β jest tak zwanym kątem ABC.

+1

Jak wyrównać punkt? – jmasterx

+0

@Milo: nie - używa "a", "b" i "c" jako odległości między parami punktów. –

+0

o dzięki, ok :) – jmasterx

1
float angba = atan((a.y - b.y)/(a.x - b.x)); 
float angbc = atan((c.y - b.y)/(c.x - b.y)); 
float rslt = angba - angbc; 
float rs = (rslt * 180)/3.141592; 
+2

Zamiast używać atana (dy/dx) lepiej używać atan2 (dy, dx) – valdo

+1

@valdo: Zgadza się. Dziękuję za przypomnienie. – oslvbo

3

Podejście z arccos jest niebezpieczne, ponieważ ryzykujemy mieć to argumentem równym, powiedzmy, 1.0000001 i skończyć z EDOMAIN błędu. Nawet podejście atan jest niebezpieczne, ponieważ obejmuje podziały, co może prowadzić do podziału przez zero błędów. Lepiej użyj wartości atan2, przekazując wartości dx i dy.

2

Tutaj jest szybki i prawidłowy sposób obliczania właściwą wartość kąta:

double AngleBetweenThreePoints(POINTFLOAT pointA, POINTFLOAT pointB, POINTFLOAT pointC) 
{ 
    float a = pointB.x - pointA.x; 
    float b = pointB.y - pointA.y; 
    float c = pointB.x - pointC.x; 
    float d = pointB.y - pointC.y; 

    float atanA = atan2(a, b); 
    float atanB = atan2(c, d); 

    return atanB - atanA; 
} 
+0

Powoduje to odpowiedź w zakresie "[-2 * pi ... + 2 * pi]" (2 obroty). – chux

0

Oto metoda OpenCV, aby uzyskać kąt między 3 punktami (A, B, C) a B jako wierzchołkiem:

int getAngleABC(cv::Point2d a, cv::Point2d b, cv::Point2d c) 
{ 
    cv::Point2d ab = { b.x - a.x, b.y - a.y }; 
    cv::Point2d cb = { b.x - c.x, b.y - c.y }; 

    float dot = (ab.x * cb.x + ab.y * cb.y); // dot product 
    float cross = (ab.x * cb.y - ab.y * cb.x); // cross product 

    float alpha = atan2(cross, dot); 

    return (int) floor(alpha * 180./M_PI + 0.5); 
} 

oparciu o doskonałe rozwiązanie przez @valdo