Нахождение применённого кватерниона

[Isaev]

1. Имеем прямоугольный параллелепипед с координатами (x1, y1, z1) и размерами (dx1, dy1, dz1) к нему был применён один из единичных кватернионов, он повернулся и он перешёл в положение (x2, y2, z2) и размеры (dx2, dy2, dz2). Каким образом это посчитать?
2. Ещё больше интересен вопрос, имея начальные и конечные координаты как вычислить кватернион поворота?

Если можно примерчик с конкретными цифрами, пожалуйста.

[arseniiv]

Isaev в сообщении #1336149 писал(а):

Каким образом это посчитать?

Что посчитать, результат поворота или какой кватернион был применён?

Вообще тут теория простейшая, см. напр. специальную отдельную страницу https://en.wikipedia.org/wiki/Quaternions_and_spatial_rotation.

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

[Isaev]

arseniiv в сообщении #1336151 писал(а):

Вообще тут теория простейшая...

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

-- 02.09.2018, 21:16 --

Начиная с простой точки:
Дана точка с координатами (6, 0, 4) и кватернион поворота (0.5, 0.5, 0.5, 0.5). Найти результирующую точку.

$$\begin{bmatrix}6 & 0 & 0 \\0 & 0 & 0 \\0 & 0 & 4\end{bmatrix} \bullet \begin{bmatrix}0 & 0 & 1 \\1 & 0 & 0 \\0 & 1 & 0\end{bmatrix} = \begin{bmatrix}0 & 0 & 6 \\0 & 0 & 0 \\0 & 4 & 0\end{bmatrix}$$

Уже тут что-то не сходится... По идее этот кватернион должен выразиться через

$$\begin{bmatrix}1 & 0 & 0 \\0 & 1 & 0 \\0 & 0 & 1\end{bmatrix}$$

и тогда всё получается правильно..., но в wiki дана именно такая формула перевода.
Почему получается это смещение строк?
И можно ли это сделать не через матрицы, а через вектора сразу?
Нигде не нашёл формулу перевода кватерниона в вектор...

[arseniiv]

Так вы вектор диагональной матрицей представили — а зачем? Это столбец должен быть. (А матрица правильная получилась. Кстати если вам важна будет производительность, её вычисление можно оптимизировать до 12 умножений и 12 сложений.) Тогда и результат будет $(4, 6, 0)^t$.

Isaev в сообщении #1336158 писал(а):

И можно ли это сделать не через матрицы, а через вектора сразу?

Можно применить кватернион к одному вектору за меньшее число операций, чем через перевод его в матрицу и её умножение на вектор (в статье выше есть формула с векторным произведением в разделе Performance comparisons / Used methods, но и то надо будет оптимизировать вычисления, иначе толку мало), но уже для двух векторов с матрицей будет быстрее, и разница тем больше, чем больше векторов нужно одинаково поворачивать.

Isaev в сообщении #1336158 писал(а):

Нигде не нашёл формулу перевода кватерниона в вектор...

А её и нет.

[Isaev]

arseniiv в сообщении #1336188 писал(а):

Так вы вектор диагональной матрицей представили — а зачем? Это столбец должен быть. (А матрица правильная получилась.

где-то в сети нашёл такой пример перевода вектора в матрицу... Теперь всё проясняется, спасибо.

arseniiv в сообщении #1336188 писал(а):

Кстати если вам важна будет производительность, её вычисление можно оптимизировать до 12 умножений и 12 сложений.) Тогда и результат будет $(4, 6, 0)^t$.

производительность будет нужна, очень... я потому до последнего и не хотел связываться с матрицами, т.к. координаты целочисленные в пределах [0..6], кватернионы единичные, для них и я вообще сделал предварительно вычисленную табличку умножения и производительность не потерялась, потому результат тоже всегда целочисленный... А операции с матрицами какие-то затратные... Потому хотелось обойтись без них, но без них получается слишком много кода, который сложно отлаживать, потому решил попробовать таки их использовать и посмотреть что в производительности потеряется.

Для второго пункта, например, чтобы по начальному и конечному положению вычислить кватернион, придётся же умножать одну матрицу на обратную другой? Чтобы найти обратную матрицу, нужно делить элементы на определитель и без плавающей арифметики не обойтись? Или таки для целочисленных вычислений есть какие-то свои методы?

arseniiv в сообщении #1336188 писал(а):

но уже для двух векторов с матрицей будет быстрее

У меня же два вектора получается? Координат и размеров... или их можно в одну матрицу запаковать?

[arseniiv]

Isaev в сообщении #1336303 писал(а):

кватернионы единичные

Ну, если они соответствуют только поворотам, которые переставляют координаты и иногда меняют им знаки, то, видимо, можно просто запоминать перестановку (и смены знака) и применять её, но это может оказаться дольше из-за неугадываний ветки условного перехода процессором, а как сделать такое без проверки условий, в отличие от умножения матрицы на столбец, не представлю. Если у вас такие повороты, то и матрицы можно брать (и куда-нибудь записать, конечно) целочисленные из $0,\pm1$ — неужели 9 целочисленных умножений и 6 сложений — так плохо?

В принципе можно попытаться заменить умножения какими-то битовыми операциями, но присутствие $-1$ усложняет дело, и я в таких совсем уж низкоуровневых оптимизациях не разбираюсь, но кто-нибудь, возможно, что-то предложит.

Isaev в сообщении #1336303 писал(а):

У меня же два вектора получается? Координат и размеров... или их можно в одну матрицу запаковать?

Можно запаковать, приставив их друг к другу, но не имею понятия, даст ли это какой-то выигрыш. Что умеют векторные операции процессора, тоже не знаю.

Isaev в сообщении #1336303 писал(а):

Для второго пункта, например, чтобы по начальному и конечному положению вычислить кватернион, придётся же умножать одну матрицу на обратную другой? Чтобы найти обратную матрицу, нужно делить элементы на определитель и без плавающей арифметики не обойтись? Или таки для целочисленных вычислений есть какие-то свои методы?

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

[Isaev]

arseniiv в сообщении #1336410 писал(а):

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

Одно умножение матриц вместо двух, надеялся я.. нет? Мы же один и тот же кватернион применяем и к координатам и к размерам объекта

arseniiv в сообщении #1336410 писал(а):

Целых чисел в любом случае не хватит, если повороты не просто переставляют координаты и меняют им знаки.

Все вращения и перемещения исключительно в границах положительного октанта, никаких смен знаков

arseniiv в сообщении #1336410 писал(а):

А если нет, то можно попробовать найти интересующую перестановку как-то в обход, смотря непосредственно кто куда переходит. Но опять не уверен, будет ли это быстрее и таким кодом, на который приятно смотреть и легко отладить.

это сейчас примерно так и работает, выглядит ужасно (потому и подумал о перестройке), работает быстро (быстродействие терять не хотелось бы)

[arseniiv]

Хм, интересно. Ну, мне пока нечего добавить.

Isaev в сообщении #1336415 писал(а):

Одно умножение матриц вместо двух, надеялся я.. нет?

Но общее число операций-то не поменяется. Вот если использовать векторные операции процессора, наверно такое объединение осмысленно. Но лучше подождать тех, кто знает — я даже не в курсе, есть ли целочисленные, или они только для флоатов.

[Isaev]

Я в общем снова запутался в теории...
У меня трёхмерное пространство, есть вектор координат кубоида (т.е. координаты его левого нижнего ближнего угла), допустим С(0,2,1) и вектор его размеров, например S(2,1,1).
[1]. Мне его нужно повернуть на 90° вокруг оси Y, проходящей через точку O(1,1,1).
У меня есть единичный кватернион поворота вокруг оси Y Q$\left( \begin{matrix}1/\sqrt{2} \\ 0 \\ 1/\sqrt{2} \\ 0 \end{matrix} \right)$. Но я же не могу его просто применить к координатам? Мне нужна матрица смещения сначала, чтобы вернуться к началу координат так? т.е.
$ \left( \begin{matrix}1 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 \\ 0 & 0 & 1 & 0 \\ -1 & -1 & -1 & 1 \end{matrix} \right)$
потом нужно кватернион перевести в матричный вид, по этой формуле
потом сделать обратное смещение на
$ \left( \begin{matrix}1 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 \\ 0 & 0 & 1 & 0 \\ 1 & 1 & 1 & 1 \end{matrix} \right)$
и в результате я ну никак не получаю в ответе C'(1,2,1). А для размеров хотелось бы получить S'(1,1,2).
и ещё получается матрица 4х4, что меня так же путает, как вернуться к трём координатам?
Это по идее получится при идентичных манипуляциях, как с координатами или тут нужно что-то другое думать?
в каком месте я не прав?

[arseniiv]

Так, ну если вы пользуетесь матрицами 4×4, чтобы делать параллельные переносы и прочие дополнительные штуки (типа центральной проекции на какой-нибудь экран), то четвёртой координатой точки надо брать единицу, и потом, если всё было корректно, она останется единицей, и её можно проигнорировать в конце, взяв опять три первые координаты. (Ещё если понадобится преобразовывать векторы, им надо брать четвёртой координатой ноль; опять же его сохранение в конце — маркер того, что код правильный.)

А вот несовпадение может быть, например, от того, что смещение было в обратную сторону. Здесь правда всё с этим хорошо, зато есть другая причина: кто-то может умножать строку из координат точки на матрицу, а кто-то матрицу на столбец координат; те, которые вы написали — для первого случая, а во втором они дадут действительно не то: столбец $\begin{bmatrix}x & y & z & 1\end{bmatrix}^t$ переведётся второй матрицей в $\begin{bmatrix}x & y & z & x+y+z+1\end{bmatrix}^t$, то есть мы сделали, если руководствоваться правилами выше, что-то не то. Если же транспонировать эти матрицы, всё будет правильно: в этом случае выйдет $\begin{bmatrix}x+1 & y+1 & z+1 & 1\end{bmatrix}^t$. Вот транспонированные в случае умножения на столбец и берут.

[Isaev]

[1] начальные координаты: $\left(\begin{matrix} 0&0&0&0 \\ 2&0&0&0 \\ 1&0&0&0 \\ 0&0&0&0 \end{matrix} \right)$

[2] матрица смещения к O: $\left(\begin{matrix} 1&0&0&0 \\ 0&1&0&0 \\ 0&0&1&0 \\ -1&-1&-1&1 \end{matrix} \right)$

[3] кватернион Q: $\left(\begin{matrix} 1/\sqrt{2} \\ 0 \\ 1/\sqrt{2} \\ 0 \end{matrix} \right)$ = $\left(\begin{matrix} 0&0&1&0 \\ 0&1&0&0 \\ -1&0&0&0 \\ 0&0&0&0 \end{matrix} \right)$

[4] матрица смещения назад: $\left(\begin{matrix} 1&0&0&0 \\ 0&1&0&0 \\ 0&0&1&0 \\ 1&1&1&1 \end{matrix} \right)$

тут почему-то смещения никакого не происходит:
[5] = [1] * [2] = $\left(\begin{matrix} 0&0&0&0 \\ 2&0&0&0 \\ 1&0&0&0 \\ 0&0&0&0 \end{matrix} \right)$

[6] = [4] * [3] = $\left(\begin{matrix} 0&0&0&0 \\ 0&0&2&0 \\ 0&0&1&0 \\ 0&0&0&0 \end{matrix} \right)$

ну и назад нет смещения:
[7] = [6] * [4] = $\left(\begin{matrix} 0&0&0&0 \\ 0&0&2&0 \\ 0&0&1&0 \\ 0&0&0&0 \end{matrix} \right)$

-- 14.05.2019, 09:56 --

arseniiv в сообщении #1392805 писал(а):

если вы пользуетесь матрицами 4×4, чтобы делать параллельные переносы и прочие дополнительные штуки (типа центральной проекции на какой-нибудь экран)

я думал это необходимость, т.к. вращение не относительно центра координат... или в этом простом случае это можно упростить?

[arseniiv]

Isaev в сообщении #1392911 писал(а):

[1] начальные координаты: $\left(\begin{matrix} 0&0&0&0 \\ 2&0&0&0 \\ 1&0&0&0 \\ 0&0&0&0 \end{matrix} \right)$

Здесь не обязательно брать прям матрицу, можно один столбец, если что. И вы забыли единицу четвёртой координатой поставить. Она в некотором смысле отвечает за учёт параллельных переносов — если вместо неё ноль, они пройдут мимо, а если $a$, они будут в $a$ раз больше, чем надо. И опять несоответствующая «раскладка» компонент. Смотрите, или перенос получается как

$$\begin{bmatrix} 1 & 0 & 0 & a_x \\ 0 & 1 & 0 & a_y \\ 0 & 0 & 1 & a_z \\ 0 & 0 & 0 & 1 \end{bmatrix} \begin{bmatrix} x \\ y \\ z \\ 1 \end{bmatrix} = \begin{bmatrix} x + a_x \\ y + a_y \\ z + a_z \\ 1 \end{bmatrix},$$

или, если всё транспонировать, как

$$\begin{bmatrix} x & y & z & 1 \end{bmatrix} \begin{bmatrix} 1 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 \\ 0 & 0 & 1 & 0 \\ a_x & a_y & a_z & 1 \end{bmatrix} = \begin{bmatrix} x + a_x & y + a_y & z + a_z & 1 \end{bmatrix}.$$

По-другому получится действительно не то.

Isaev в сообщении #1392911 писал(а):

я думал это необходимость, т.к. вращение не относительно центра координат... или в этом простом случае это можно упростить?

Смотря что значит упростить. Матрицы используются, потому что можно всё преобразование целиком, ещё и со всякими перспективными искажениями, потом задать одной такой матрицей и наприменять её ко всем нужным точкам. Если например для преобразований вида $A\vec v + \vec b$ хранить матрицу $A$ и вектор $\vec b$, то композиция преобразований $(A,\vec b)$ и $(A',\vec b')$ будет равна $(AA', A\vec b' + \vec b)$, а если затолкать это всё в одну большую матрицу, мы их просто умножаем и не паримся. В оптимизированных алгоритмах в принципе может быть и так, и так, в зависимости от того, какие векторные операции процессор поддерживает (тут я не разбираюсь).

[Isaev]

Заработал разворот, это большой шаг вперёд, огромное спасибо!

arseniiv в сообщении #1392994 писал(а):

Матрицы используются, потому что можно всё преобразование целиком, ещё и со всякими перспективными искажениями, потом задать одной такой матрицей и наприменять её ко всем нужным точкам. Если например для преобразований вида $A\vec v + \vec b$ хранить матрицу $A$ и вектор $\vec b$, то композиция преобразований $(A,\vec b)$ и $(A',\vec b')$ будет равна $(AA', A\vec b' + \vec b)$, а если затолкать это всё в одну большую матрицу, мы их просто умножаем и не паримся.

1. Тут как бы запаковать перенос туда и обратно не получится же? они взаимоисключающие... всё равно придётся их использовать отдельно друг от друга?

Вопрос номер 2.: С размером, естественно, по той же схеме не работает (и не должно по идее)... Как его посчитать?
это какая-то проекция чтоли, даже не представляю... как бы объяснить... Идея в том, что если объект имеет размеры (2, 3, 4), мы поворачиваем вокруг оси Y, при этом конечно размеры не меняются. А вот если мы теряем информацию о повороте, то перед нами объект с размерами (4, 3, 2). Как это вычислить?

И номер 3.: имея начальный и конечный вектор получить кватернион поворота. Тут как я понимаю порядок действий обратный, но делить мы матрицы не можем, потому умножаем конечную на инвертированную начальную?

[arseniiv]

Isaev в сообщении #1393178 писал(а):

1. Тут как бы запаковать перенос туда и обратно не получится же? они взаимоисключающие... всё равно придётся их использовать отдельно друг от друга?

Не, можно. Перенос на $\vec b$, поворот $R$ и перенос обратно дают $R(\vec v + \vec b) - \vec b = R\vec v + R\vec b - \vec b$, что нормально задаётся парой $(R, R\vec b - \vec b)$. Вот проективные преобразования сюда уже не влезут.

Про размеры. Ну, если повороты оставляют кирпич ориентированным по осям, то если взять векторы $\Delta x\vec\imath, \Delta y\vec\jmath, \Delta z\vec k$ и повернуть так же, новое $\Delta x'$ будет модулем того из результатов, который коллинеарен $\vec\imath$ и т. д.. Практически, если не прибегать всё-таки к явной записи перестановок (со знаками) вместо поворотов, можно это сделать, записав $\Delta x\vec\imath, \Delta y\vec\jmath, \Delta z\vec k$ друг за другом в матрицу (получится диагональная $\operatorname{diag}(\Delta x, \Delta y, \Delta z)$), действовать поворотом на неё, а потом взять максимальное по модулю число в первой строке как $\Delta x'$, максимальное по модулю во второй строке как $\Delta y'$, ну и с третьей $\Delta z'$. (Беря максимум по модулю, избежим потенциального превращения нулей во что-то маленькое. Впрочем если у вас можно сделать всё целочисленным, можно будет прямо проверять на не нулевые значения, и даже просто брать суммы элементов строки, что должно быть чуть быстрее на процессорах, угадывающих альтернативы if.)

Но если у вас только повороты—симметрии куба (переставляющие декартовы оси местами, плюс-минус знаки), я бы это закодировал (что упомянуто выше) именно тем, какая перестановка и какое изменение знаков. Тогда перестановку можно применять к $(\Delta x,\Delta y,\Delta z)$ и получать новые размерности кирпича, да и к координатам точек тоже можно применять, ну а перенос делать отдельно, что уж. Где-то с кем-то вроде эти перестановки обсуждали недавно, но забыл где и с кем.

Isaev в сообщении #1393178 писал(а):

И номер 3.: имея начальный и конечный вектор получить кватернион поворота. Тут как я понимаю порядок действий обратный, но делить мы матрицы не можем, потому умножаем конечную на инвертированную начальную?

Не, тут хитрее. Во-первых, можно повернуть один вектор в другой не единственным способом. Ну, вам наверно нужен был поворот на наименьший угол — тогда всё почти всегда хорошо, если только векторы не противоположны. Тогда есть куча равноценных поворотов (на 180°), где ничего уже не уточнить. В этом случае стоит поворачивать пару (неколлинеарных) векторов в пару векторов — вот тогда результат будет однозначен (но и угол не обязательно будет наименьшим).

Определять для одиночных векторов $\vec v,\vec v'$ можно, найдя, раз тут трёхмерие, вектор, ортогональный плоскости $\langle\vec v,\vec v'\rangle$, просто взяв их векторное произведение. Заодно оно же даст и угол, на который надо будет поворачивать. По оси и углу можно найти матрицу поворота.

Для пар векторов иначе: дополняем их до правой тройки каждую, а вот уже после этого можно записать из столбцов координат первых матрицу $A$, из столбцов вторых матрицу $B$, и взять $BA^{-1}$.

(В представлении перестановками можно было бы просто посмотреть, куда переходит $(x,0,0)$, куда $(0,y,0)$, куда $(0,0,z)$, причём третий вектор по двум заданным можно опять определить векторным произведением, или сократить вычисления, просто поставив единицу в оставшейся координате и поставив там минус, если тройка получается левая.)

[Isaev]

arseniiv в сообщении #1393245 писал(а):

Перенос на $\vec b$, поворот $R$ и перенос обратно дают $R(\vec v + \vec b) - \vec b = R\vec v + R\vec b - \vec b$, что нормально задаётся парой $(R, R\vec b - \vec b)$.

Можно это в виде матрицы, если не сложно?
При переносе, кстати, в моём случае всегда $\Delta x = \Delta y = \Delta z$, думаю это тоже может всё упростить...

arseniiv в сообщении #1393245 писал(а):

Ну, если повороты оставляют кирпич ориентированным по осям

Оставляют, повороты только на единичные кватернионы, потому ориентация по осям никогда не нарушается.

arseniiv в сообщении #1393245 писал(а):

Впрочем если у вас можно сделать всё целочисленным, можно будет прямо проверять на не нулевые значения

Да, я так и обнулял эти ошибки переполнения.

arseniiv в сообщении #1393245 писал(а):

Но если у вас только повороты—симметрии куба (переставляющие декартовы оси местами, плюс-минус знаки), я бы это закодировал (что упомянуто выше) именно тем, какая перестановка и какое изменение знаков. Тогда перестановку можно применять к $(\Delta x,\Delta y,\Delta z)$ и получать новые размерности кирпича, да и к координатам точек тоже можно применять

Да, именно тот случай, только симметрии...
Кстати да, очень не плохая идея. Я на это смотрел с позиции объекта, а не поворота и получалось очень не красиво это всё хранить.

arseniiv в сообщении #1393245 писал(а):

Во-первых, можно повернуть один вектор в другой не единственным способом. Ну, вам наверно нужен был поворот на наименьший угол — тогда всё почти всегда хорошо, если только векторы не противоположны. Тогда есть куча равноценных поворотов (на 180°), где ничего уже не уточнить.

Вроде в случае единичных кватернионы всё вполне однозначно должно получаться. Можно конечно перемножить множество различных поворотов и прийти к 180, но путь нам не важен, а результат всегда однозначен.

Далее по третьему пункту запутался, можно пример с конкретными цифрами?