Как быстро посчитать длину вектора в пространстве
Продолжаем закрывать прорехи в программном коде "изделия", до которых раньше руки не доходили. Координаты объекта мы уже измерили, и в декартовых координатах выдаём без проблем. Но надо ещё и дальность посчитать, причём "честно" возводить всё в квадрат и извлекать корень совсем не хочется.
Ранее я рассматривал, как оценить длину вектора на плоскости ( раз, два), но там задача была чуть другая. Не было никаких ограничений на величину каждой компоненты, при этом большой точности вычислений не требовалось. Здесь же не плоскость, а пространство, но зато мы знаем, что оно ограничено полем зрения нашего изделия, плюс-минус 12 градусов по "горизонтали" и "вертикали". Как видно, это позволяет добиться неплохой точности, используя совсем простые формулы, просто взвешенное сложение модулей!
Даже простое отбрасывание значений Y,Z даёт приличный результат. Ошибка примерно квадратично растёт по мере роста Y,Z, достигая максимума в углу нашего поля зрения, а именно: 4,32%.
Точность можно повысить примерно вдвое, просто смасштабировав результат! Людям, воспитанным на рядах Тэйлора, это может показаться жутким кощунством, а вот по Чебышёву/Ремезу это вполне естественно. Чтобы минимизировать максимальную ошибку приближённой формулы (т.е подбирать формулу по МИНИМАКСУ), надо, чтобы одна и та же максимальная ошибка достигалась на краях диапазона и на всех максимумах/минимумах, с обязательным чередованием знака.
Здесь пока нет максимумов и минимумов, так что в центре поля зрения мы даём результат, завышенный на 2,2%, а по мере отдаления в углы поля зрения, какое-то время считаем точно, а потом начинаем давать результат, заниженный вплоть до 2,2%. До сих пор немного удивляет, как можно ошибиться в самом примитивном случае, когда Y=Z=0, но что делать.
Когда мы развлекались на плоскости, мы аппроксимировали Эвклидову метрику взвешенной суммой более простых для вычисления метрик: метрикой городских кварталов (Манхэттенская) и метрикой Чебышёва. Но здесь, поскольку заведомо X больше по модулю, чем Y или Z, Чебышёвская метрика превращается просто в |X|. Так что всё сводится попросту к |X| с одним весом, очень близким к единице, и |Y|+|Z| с другим весом.
И следующая по сложности формула - если мы оставим |X| с единичным весом, что кажется наиболее естественным. Единственный коэффициент подбирается так, чтобы относительная ошибка "в плюс" уравнялась с относительной ошибкой "в минус", вот что мы получаем:
Меня, честно говоря, удивило, что кривая ошибки "по диагонали" оказалась практически смасштабированной версией ошибки "по горизонтали". Я как-то думал, что где-то там выползет корень из двух, но нет. Когда вклад от компонент Y,Z мал в сравнении с X, двойка под корнем спокойно остаётся двойкой, т.е:
^2}\approx|X|(1+k/2(Y/X)^2))
Если мы хорошо подобрали линейное приближение для k=1 (т.е "работа по горизонтали", когда одна из компонент нулевая), то это же приближение столь же хорошо сработает и для k=2.
Следующий шаг - "отпустить" коэффициент при |X|: внеся ошибку в районе нуля, её можно будет снизить во всей остальной области. И действительно, так и получается:
Строго по Чебышёву: ошибка на краях совпадает с ошибкой в максимуме, при этом знак ошибки чередуется. Ошибка не более 0,55% - для нашего применения вполне достаточно. Вообще, эта опция - чисто "для галочки". В ТЗ написано, что мы должны измерять дальность, но на самом деле заказчик всё равно будет брать декартовы координаты, чтобы легко можно было переместить начало отсчёта куда надо. В сферических координатах это делается куда хуже, фактически переходом назад в декартовы! В общем, такая формула вполне позволяет вписаться в требования ТЗ, а большево нам и не надо!
Это, разумеется, не единственный способ. Так-то я заранее знаю, что значение X будет в диапазоне от 1 до 2 ("псевдо-плавающая точка, где экспонента общая для всего вектора). Я пробовал другие приближения, более похожие на разложение в ряд Тэйлора, но с подбором коэффициентов по минимаксу, однако, получалось гораздо хуже, единицы процентов. У рассмотренных выше формул хорошая особенность - когда масштабируешь вектор, точно также масштабируется и результат, что в нашем случае "решает".
И ещё не надо сбрасывать со счетов мой вариант Cordic'а (см. atan(y/x) на QuatCore), там мы серией поворотов притаскиваем вектор поближе к оси X, причём по ходу дела подсчитываем, на какой угол довернули. Так можно одновременно узнать и длину (теперь уже длина = X без зазрения совести, когда углы фактически устранены), и угол поворота, т.е тот самый арктангенс. Но там речь шла о векторе на плоскости, так что метод всё равно надо приспосабливать для пространства - сначала выкрутить на плоскость, а потом уже стянуть в прямую. Так что пока склоняюсь к вычислению углов и дальности "по отдельности", а там видно будет.
Ранее я рассматривал, как оценить длину вектора на плоскости ( раз, два), но там задача была чуть другая. Не было никаких ограничений на величину каждой компоненты, при этом большой точности вычислений не требовалось. Здесь же не плоскость, а пространство, но зато мы знаем, что оно ограничено полем зрения нашего изделия, плюс-минус 12 градусов по "горизонтали" и "вертикали". Как видно, это позволяет добиться неплохой точности, используя совсем простые формулы, просто взвешенное сложение модулей!
Даже простое отбрасывание значений Y,Z даёт приличный результат. Ошибка примерно квадратично растёт по мере роста Y,Z, достигая максимума в углу нашего поля зрения, а именно: 4,32%.
Точность можно повысить примерно вдвое, просто смасштабировав результат! Людям, воспитанным на рядах Тэйлора, это может показаться жутким кощунством, а вот по Чебышёву/Ремезу это вполне естественно. Чтобы минимизировать максимальную ошибку приближённой формулы (т.е подбирать формулу по МИНИМАКСУ), надо, чтобы одна и та же максимальная ошибка достигалась на краях диапазона и на всех максимумах/минимумах, с обязательным чередованием знака.
Здесь пока нет максимумов и минимумов, так что в центре поля зрения мы даём результат, завышенный на 2,2%, а по мере отдаления в углы поля зрения, какое-то время считаем точно, а потом начинаем давать результат, заниженный вплоть до 2,2%. До сих пор немного удивляет, как можно ошибиться в самом примитивном случае, когда Y=Z=0, но что делать.
Когда мы развлекались на плоскости, мы аппроксимировали Эвклидову метрику взвешенной суммой более простых для вычисления метрик: метрикой городских кварталов (Манхэттенская) и метрикой Чебышёва. Но здесь, поскольку заведомо X больше по модулю, чем Y или Z, Чебышёвская метрика превращается просто в |X|. Так что всё сводится попросту к |X| с одним весом, очень близким к единице, и |Y|+|Z| с другим весом.
И следующая по сложности формула - если мы оставим |X| с единичным весом, что кажется наиболее естественным. Единственный коэффициент подбирается так, чтобы относительная ошибка "в плюс" уравнялась с относительной ошибкой "в минус", вот что мы получаем:
Меня, честно говоря, удивило, что кривая ошибки "по диагонали" оказалась практически смасштабированной версией ошибки "по горизонтали". Я как-то думал, что где-то там выползет корень из двух, но нет. Когда вклад от компонент Y,Z мал в сравнении с X, двойка под корнем спокойно остаётся двойкой, т.е:
Если мы хорошо подобрали линейное приближение для k=1 (т.е "работа по горизонтали", когда одна из компонент нулевая), то это же приближение столь же хорошо сработает и для k=2.
Следующий шаг - "отпустить" коэффициент при |X|: внеся ошибку в районе нуля, её можно будет снизить во всей остальной области. И действительно, так и получается:
Строго по Чебышёву: ошибка на краях совпадает с ошибкой в максимуме, при этом знак ошибки чередуется. Ошибка не более 0,55% - для нашего применения вполне достаточно. Вообще, эта опция - чисто "для галочки". В ТЗ написано, что мы должны измерять дальность, но на самом деле заказчик всё равно будет брать декартовы координаты, чтобы легко можно было переместить начало отсчёта куда надо. В сферических координатах это делается куда хуже, фактически переходом назад в декартовы! В общем, такая формула вполне позволяет вписаться в требования ТЗ, а большево нам и не надо!
Это, разумеется, не единственный способ. Так-то я заранее знаю, что значение X будет в диапазоне от 1 до 2 ("псевдо-плавающая точка, где экспонента общая для всего вектора). Я пробовал другие приближения, более похожие на разложение в ряд Тэйлора, но с подбором коэффициентов по минимаксу, однако, получалось гораздо хуже, единицы процентов. У рассмотренных выше формул хорошая особенность - когда масштабируешь вектор, точно также масштабируется и результат, что в нашем случае "решает".
И ещё не надо сбрасывать со счетов мой вариант Cordic'а (см. atan(y/x) на QuatCore), там мы серией поворотов притаскиваем вектор поближе к оси X, причём по ходу дела подсчитываем, на какой угол довернули. Так можно одновременно узнать и длину (теперь уже длина = X без зазрения совести, когда углы фактически устранены), и угол поворота, т.е тот самый арктангенс. Но там речь шла о векторе на плоскости, так что метод всё равно надо приспосабливать для пространства - сначала выкрутить на плоскость, а потом уже стянуть в прямую. Так что пока склоняюсь к вычислению углов и дальности "по отдельности", а там видно будет.