Ответ 1
Первые предложения относительно вашего метода:
То, что вы называете 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 the only invocation of the heavy function.
// It 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
для.
Вращение достигается путем умножения вектора на следующую матрицу:
- a.x, a.y
- -a.y, a.x
Как только можно увидеть, что вектор a
, умноженный на такую матрицу, действительно вращается к положительной оси X.
Примечание: Строго говоря, приведенная выше матрица не просто вращается, но и масштабируется. Но в нашем случае это нормально, так как единственное, что имеет значение, - это векторное направление, а не его длина.
Вращающийся вектор b
становится:
- a.x * b.x + a.y * b.y = a точка b
- -a.y * b.x + a.x * b.y = a крест b
Наконец, ответ может быть выражен как:
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);
}