2010-08-15 6 views
14

Учитывая точки ABC, как я мог найти угол ABC? Я делаю инструмент feehand для приложения векторного рисования и для сведения к минимуму количества генерируемых очков, я не буду добавлять точки, если угол положения мыши и последние 2 точки больше определенного порога. БлагодаряУгол между 3 точками?

что я имел:

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

Я делаю очень хорошо, у меня есть для этого нового алгоритма, но это не делает трюк. – jmasterx

+5

@abelenky: и это задает вопрос «неясно или не полезно», как именно? Возможно, вы неправильно поняли цель репутации. Нельзя позволить вам наказывать людей за то, что они пытались сделать что-то новое для них. – jalf

ответ

25

Первые предложения относительно вашего метода:

Что вы называете ac фактически cb. Но все в порядке, это действительно необходимо. Далее

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

Это ваша первая ошибка. реального скалярного произведения двух векторов:

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

Теперь

float rslt = acos(dacos); 

Здесь следует отметить, что в связи с некоторой потерей точности при расчете теоретически возможно, что dacos станет больше, чем 1 (или меньше -1). Следовательно, вы должны проверить это явно.

Плюс заметка о производительности: вы называете тяжелую функцию sqrt дважды для вычисления длины двух векторов. Затем вы делите точечный продукт на эти длины. Вместо этого вы могли бы позвонить sqrt о умножении квадратов длины обоих векторов.

И, наконец, вы должны отметить, что ваш результат является точным до sign. То есть ваш метод не будет отличать 20 ° и -20 °, так как косинус обоих одинаковый. Ваш метод даст одинаковый угол для ABC и CBA.

Один правильный метод для вычисления угла как "oslvbo" предлагает:

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

(я только что заменил atan на atan2).

Это самый простой способ, который всегда дает правильный результат. Недостатком этого метода является то, что вы на самом деле называете тяжелую тригонометрическую функцию atan2 дважды.

Предлагаю следующий метод. Это немного сложнее (требуется понимание некоторых навыков тригонометрии), однако он превосходит точки зрения производительности. Он просто вызывает функцию тригонометрии atan2. И нет вычислений квадратного корня.

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:

Недавно я работал на схожую тему. И тогда я понял, что есть лучший способ. Это на самом деле более или менее одинаковое (за кулисами). Однако это более прямое ИМХО.

Идея состоит в том, чтобы повернуть оба вектора так, чтобы первый был выровнен по (положительному) X-направлению. Очевидно, что вращение обоих векторов не влияет на угол между ними. OTOH после такого вращения нужно просто выяснить угол второго вектора относительно оси X. И это точно для atan2.

Вращение достигается путем умножения вектора на следующей матрицей:

  • Ax, Ay
  • -ay, топор

После того, как может видеть, что вектор a, умноженную на такой матрице в самом деле поворачивается к положительной оси X.

Примечание: Строго говоря, приведенная выше матрица не только вращается, но и масштабируется. Но это нормально в нашем случае, поскольку единственное, что имеет значение, - это направление вектора, а не его длина.

повернутого вектора b становится:

  • A.x * b.x + a.y * b.y = точка б
  • -a.y * b.x + A.x * б.у = крест б

Наконец, ответ может быть выражен как:

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

Отличное решение! Я беспокоюсь о том, как справиться с проблемой угла знака, пока не прочитаю ваш ответ. –

+0

Последняя функция (у которой есть «return (int) floor (alpha * 180./pi + 0.5);'), хороша и дает другой ответ к abc и cba. Работает хорошо! –

+0

Только решение, которое не давало мне точных проблем. Спасибо! – Tiago

1

Off тему? Но вы можете сделать это с законом косинусов:

Найдите расстояние между A и B (назовите это x) и расстояние между B и C (вызовите это y) и расстояние между A и C (вызов это z).

Тогда вы знаете, что г^2 = х^2 + у^2-2 * х у соз (УГОЛ ХОЧЕШЬ)

таким образом, что угол соз^-1 ((г^2 -x^2-у^2)/(2х)) = УГОЛ

+0

Вы потеряли негатив во фракции. В моем ответе должно быть «x^2 + y^2 - z^2». –

4

& beta; = агссоз ((а^2 + с^2 - Ь^2)/2ac)

где а сторона, противоположное угол α, b - противоположный угол β, а c - противоположный угол γ. Итак, β - это то, что вы назвали угол ABC.

+1

Как вы делаете квадрат точки? – jmasterx

+0

@ Мило: у вас нет - он использует 'a',' b' и 'c' как расстояния между парами очков. –

+0

oh ok thanks :) – 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

Вместо использования atan (dy/dx) лучше использовать atan2 (dy, dx) – valdo

+1

@valdo: Это правильно. Спасибо за ваше напоминание. – oslvbo

3

подход с arccos опасно, потому что мы рискуем иметь это аргумент равен, скажем, 1.0000001 и в конечном итоге с EDOMAIN ошибки. Даже подход atan опасен, поскольку он включает в себя деления, которые могут привести к делению на нулевую ошибку. Лучше использовать atan2, передавая dx и dy значениям.

2

Вот быстрый и правильный способ расчета правильного значения угла:

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

Возвращает ответы в диапазоне '[-2 * pi ... + 2 * pi]' (2 оборота). – chux

0

Вот путь OpenCV, чтобы получить угол между 3 точками (A, B, C) с B как вершиной:

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); 
} 

Основываясь на отличное решение по @valdo