Дешевый алгоритм, чтобы найти меру угла между векторами - PullRequest
28 голосов
/ 15 сентября 2009

Найти угол между двумя векторами несложно , используя правило косинуса . Однако, поскольку я программирую для платформы с очень ограниченными ресурсами, я бы хотел избежать таких вычислений, как sqrt и arccos. Даже простые деления должны быть максимально ограничены.

К счастью, мне не нужен угол как таковой, а нужно только некоторое значение, пропорциональное указанному углу.

Так что я ищу какой-то вычислительно дешевый алгоритм для вычисления величины, связанной с углом между двумя векторами. До сих пор я не нашел чего-то, что отвечало бы всем требованиям, и при этом я не смог придумать что-то сам.

Ответы [ 9 ]

16 голосов
/ 03 февраля 2013

Если вам не нужен фактический евклидов угол, но что-то, что вы можете использовать в качестве основы для сравнения углов, тогда выбор геометрии такси может быть выбором, потому что вы можете отбросить тригонометрию и ее медленность при поддержании ТОЧНОСТИ ( или, по крайней мере, с очень незначительной потерей точности, см. ниже).

В основных современных браузерных движках коэффициент ускорения составляет от 1,44 до 15,2, а точность почти такая же, как в atan2. Расчет угла алмаза в среднем в 5,01 раза быстрее, чем atan2, а при использовании встроенного кода в Firefox 18 ускорение достигает коэффициента 15,2. Сравнение скорости: http://jsperf.com/diamond-angle-vs-atan2/2.

Код очень прост:

function DiamondAngle(y, x)
{
    if (y >= 0)
        return (x >= 0 ? y/(x+y) : 1-x/(-x+y)); 
    else
        return (x < 0 ? 2-y/(-x-y) : 3+x/(x-y)); 
}

Приведенный выше код дает вам угол между 0 и 4, а atan2 - угол между -PI и PI, как показано в следующей таблице:

enter image description here

Обратите внимание, что угол ромба всегда положителен и находится в диапазоне 0-4, в то время как atan2 дает также отрицательные радианы. Таким образом, угол алмаза более нормализован. И еще одно замечание: atan2 дает немного более точный результат, потому что длина диапазона равна 2 * pi (т.е. 6,283185307179586), тогда как в углах алмаза она равна 4. На практике это не очень важно, например. Рад 2.3000000000000001 и 2.3000000000000002 имеют оба угла ромба 1.4718731421442295, но если мы снизим точность, опустив один ноль, то рад 2.300000000000001 и 2.300000000000002 дают оба угла ромба. Эта «потеря точности» в углах алмаза настолько мала, что оказывает существенное влияние, только если расстояния огромны. Вы можете играть с конверсиями в http://jsbin.com/bewodonase/1/edit?output (Старая версия: http://jsbin.com/idoyon/1):

enter image description here

Приведенного выше кода достаточно для быстрого сравнения углов, но во многих случаях необходимо преобразовать угол ромба в радианы и наоборот. Если вы, например. иметь некоторый допуск в виде радиальных углов, а затем у вас есть петля в 100 000 раз, где этот допуск сравнивается с другими углами, сравнивать с использованием atan2 нецелесообразно. Вместо этого перед выполнением цикла вы изменяете допуск в радианах на допуск таксика (алмазные углы) и выполняете внутриконтурные сравнения с использованием допусков алмаза, и таким образом вам не нужно использовать медленные тригонометрические функции в критических по скорости частях кода (= в петли).

Код, который делает это преобразование:

function RadiansToDiamondAngle(rad)
{
  var P = {"x": Math.cos(rad), "y": Math.sin(rad) };
  return DiamondAngle(P.y, P.x);
}  

Как вы заметили, есть cos и sin. Как вы знаете, они медленные, но вам не нужно делать преобразование в цикле, но перед циклом и ускорение огромно.

И если по какой-то причине вам нужно преобразовать угол бриллианта в радианы, например. после зацикливания и сравнения углов для возврата, например. минимальный угол сравнения или что-либо в радианах, код выглядит следующим образом:

function DiamondAngleToRadians(dia)
{
  var P = DiamondAngleToPoint(dia);
  return Math.atan2(P.y,P.x);
}

function DiamondAngleToPoint(dia)
{
  return {"x": (dia < 2 ? 1-dia : dia-3), 
  "y": (dia < 3 ? ((dia > 1) ? 2-dia : dia) : dia-4)};
}

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

Этого должно быть достаточно для быстрого сравнения углов.

Конечно, есть другие методы ускорения atan2 (например, CORDIC и справочные таблицы), но AFAIK все они теряют точность и все же могут быть медленнее, чем atan2.

ПРЕДПОСЫЛКИ: Я проверил несколько методов: точечные произведения, внутренние произведения, закон косинуса, единичные окружности, справочные таблицы и т. Д., Но ничего не было достаточно в случае, когда важны как скорость, так и точность. Наконец я нашел страницу в http://www.freesteel.co.uk/wpblog/2009/06/encoding-2d-angles-without-trigonometry/, в которой были нужные функции и принципы.

Сначала я предположил, что для точного и быстрого сравнения расстояний можно использовать также расстояния такси, потому что большее расстояние в евклидово больше и в такси. Я понял, что вопреки евклидовым расстояниям угол между начальной и конечной точкой влияет на расстояние такси. Только длины вертикальных и горизонтальных векторов можно легко и быстро преобразовать между евклидовым и таксиским, но в любом другом случае необходимо учитывать угол, и тогда процесс идет слишком медленно (?).

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

12 голосов
/ 15 сентября 2009

Вы пробовали алгоритм CORDIC ? Это общая схема для решения полярной & Harr; прямоугольные задачи только с добавлением / вычитанием / битовым сдвигом + таблицей, по существу, выполняющие поворот на углы вида tan -1 (2 -n ). Вы можете поменять точность со временем выполнения, изменив количество итераций.

В вашем случае возьмите один вектор в качестве фиксированной привязки, а другой скопируйте во временный вектор, который вы поворачиваете, используя кордовые углы к первому вектору (грубо раздвоенные), пока не достигнете желаемой угловой точности.

( edit: использовать знак продукта точки для определения на каждом шаге, вращаться ли вперед или назад. Хотя, если умножения достаточно дешевы, чтобы разрешить использование продукта точки, то не стоит беспокоиться о CORDIC, возможно используйте таблицу пар sin / cos для матриц вращения углов & pi; / 2 n , чтобы решить проблему с делением пополам.)

( edit: Мне нравится предложение Эрика Бейнвилла в комментариях: поверните оба вектора к нулю и отследите разницу углов).

12 голосов
/ 15 сентября 2009

В те времена, когда у нас было несколько килобайт ОЗУ и машин с ограниченными математическими возможностями, я использовал таблицы поиска и линейную интерполяцию. Основная идея проста: создать массив с таким разрешением, которое вам необходимо (больше элементов уменьшает ошибку, создаваемую интерполяцией). Затем интерполируйте между поисковыми значениями.

Вот пример в обработке ( исходная битая ссылка ).

Вы можете сделать это и с другими функциями триггера. На процессоре 6502 это позволило рассчитать полную трехмерную графику в виде каркаса с увеличением скорости на порядок.

4 голосов
/ 03 июня 2013

Здесь, на ТАК, я до сих пор не имею права комментировать (хотя у меня есть на math.se), так что это на самом деле ответ на пост Тимо об алмазных углах.

Вся концепция углов алмаза, основанная на норме L1, наиболее интересна, и если бы это было просто сравнение того, какой вектор имеет большую / меньшую w.r.t. положительной оси X было бы достаточно. Тем не менее, OP упомянул угол между двумя общими векторами, и я предполагаю, что OP хочет сравнить его с некоторым допуском для определения состояния гладкости / угла или чего-то подобного, но, к сожалению, кажется, что только с формулами, представленными на jsperf.com или freesteel.co.uk (ссылки выше), кажется, что это невозможно сделать, используя алмазные углы.

Обратите внимание на следующий вывод из моей реализации асимптотических формул:

Vectors : 50,20 and -40,40
Angle diff found by acos      : 113.199
Diff of angles found by atan2 : 113.199
Diamond minus diamond         : 1.21429
Convert that to degrees       : 105.255
Rotate same vectors by 30 deg.
Vectors : 33.3013,42.3205 and -54.641,14.641
Angle diff found by acos      : 113.199
Diff of angles found by atan2 : 113.199
Diamond minus diamond         : 1.22904
Convert that to degrees       : 106.546
Rotate same vectors by 30 deg.
Vectors : 7.67949,53.3013 and -54.641,-14.641
Angle diff found by acos      : 113.199
Diff of angles found by atan2 : 113.199
Diamond minus diamond         : 1.33726
Convert that to degrees       : 116.971

Таким образом, дело в том, что вы не можете сделать алмазную (альфа) -бриллиант (бета) и сравнить ее с некоторым допуском, в отличие от того, что вы можете сделать с выходом atan2. Если все, что вы хотите сделать, это алмаз (альфа)> алмаз (бета), то я думаю, что с бриллиантом все в порядке.

2 голосов
/ 15 сентября 2009

Решение было бы тривиальным, если бы векторы были определены / сохранены с использованием полярных координат вместо декартовых координат (или «а также» с использованием декартовых координат).

1 голос
/ 19 декабря 2010

Перекрестное произведение пропорционально углу между двумя векторами, и когда векторы нормализованы и угол мал, перекрестное произведение очень близко к фактическому углу в радианах из-за приближения малого угла.

в частности:

I1Q2-I2Q1 пропорционально углу между I1Q1 и I2Q2.

1 голос
/ 15 сентября 2009

если вам нужно вычислить квадратный корень, рассмотрите возможность использования invsqrt hack .

acos((x1*x2 + y1*y2) * invsqrt((x1*x1+y1*y1)*(x2*x2+y2*y2)));
1 голос
/ 15 сентября 2009

произведение двух векторов (x1, y1) и (x2, y2) равно

x1 * x2 + y1 * y2 

и эквивалентно произведению длин двух векторов на косинус угла между ними.

Итак, если вы сначала нормализуете два вектора (разделите координаты на длину)

Where length of V1 L1 = sqrt(x1^2 + y1^2),  
  and length of V2 L2 = sqrt(x2^2 + y2^2),

Тогда нормализованные векторы равны

(x1/L1, y1/L1),  and (x2/L2, y2/L2),  

И скалярное произведение нормализованных векторов (которое совпадает с косинусом угла между векторами) будет

 (x1*x2 + y1*y2)
 -----------------
     (L1*L2)

конечно, это может быть столь же сложно в вычислительном отношении, как вычисление косинуса

0 голосов
/ 15 сентября 2009

Точечный продукт может работать в вашем случае. Он не пропорционален углу, а «связан».

Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...