Углы Эйлера и матрица вращения из двух трехмерных точек

Я пытаюсь найти углы Эйлера, которые позволяют преобразовать точку A в точку B в трехмерном пространстве.

Рассмотрим нормализованные векторы A = [1, 0, 0] и B = [0.32 0.88 -0.34].

Я понимаю, что вычисляя перекрестное произведение A × B, я получаю ось вращения. Угол между A и B определяется как tan⁻¹(||cross||, A·B), где A·B - это скалярное произведение между A и B.

Это дает мне вектор вращения rotvec = [0 0.36 0.93 1.24359531111], который равен rotvec = [A × B; angle] (перекрестное произведение нормализовано).

Теперь мой вопрос: Как мне перейти отсюда, чтобы получить углы Эйлера, соответствующие преобразованию от A к B?

В MATLAB функция vrrotvec2mat получает в качестве входных данных вектор вращения и выводит матрицу вращения. . Затем функция rotm2eul должна вернуть соответствующие углы Эйлера. Я получаю следующий результат (в радианах): [0.2456 0.3490 1.2216] в соответствии с соглашением XYZ. Однако это не ожидаемый результат.

Правильный ответ - [0 0.3490 1.2216], что соответствует повороту на 20° и 70° в Y и Z соответственно.

Когда я использую eul2rot([0 0.3490 1.2216])eul2rot, взятым из здесь) для проверки Полученная матрица вращения отличается от той, которую я получаю при использовании vrrotvec2mat(rotvec).

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

--- Python (2.7) с использованием transform3d ---

import numpy as np
import transforms3d

cross = np.cross(A, B)
dot = np.dot(A, B.transpose())
angle = math.atan2(np.linalg.norm(cross), dot)
rotation_axes = sklearn.preprocessing.normalize(cross)
rotation_m = transforms3d.axangles.axangle2mat(rotation_axes[0], angle, True)
rotation_angles = transforms3d.euler.mat2euler(rotation_m, 'sxyz')

Что мне здесь не хватает? Что мне делать вместо этого?

Спасибо


person jcampos    schedule 27.07.2018    source источник
comment
Углы Эйлера между двумя точками не уникальны. Вот почему существуют кватерионы   -  person Ander Biguri    schedule 28.07.2018
comment
@AnderBiguri Я действительно не понимаю, как кватернионы имеют отношение к вопросу. Или вы предполагаете, что это проблема XY?   -  person jodag    schedule 28.07.2018
comment
Получу ли я уникальное решение, если вместо этого использую кватернион? В конце концов, мне все еще нужно преобразовать его в углы Эйлера.   -  person jcampos    schedule 28.07.2018
comment
Нет. Единица-кватернион - это альтернативное представление вращения в R3. Она более компактна, численно стабильна и быстрее применяется в вычислениях, чем матрица вращения, но не обладает большей репрезентативной способностью. Он менее компактен, чем углы Эйлера, но не страдает проблемой блокировки кардана. Факт остается фактом: существует несколько поворотов, которые удовлетворяют вашу проблему, поэтому вам нужно будет придумать способ наложения дополнительных ограничений, чтобы получить желаемое решение.   -  person jodag    schedule 28.07.2018


Ответы (2)


Матрица вращения имеет 3 степени свободы, но ограничения вашей задачи ограничивают только 2 из этих степеней.

Это можно сделать более конкретным, рассмотрев случай, когда у нас есть матрица вращения R, которая вращается с A на B, поэтому R*A == B. Если мы затем построим другую матрицу вращения RB, которая вращается вокруг вектора B, тогда применение этого поворота к R*A не будет иметь никакого эффекта, то есть B == R*A == RB*R*A. Однако это создаст другую матрицу вращения RB*R с другими углами Эйлера.

Вот пример в MATLAB:

A = [1; 0; 0];
B = [0.32; 0.88; -0.34];

A = A / norm(A);
B = B / norm(B);

ax = cross(A, B);
ang = atan2(norm(ax), dot(A, B)); % ang = acos(dot(A, B)) works too
R = axang2rotm([ax; ang].');

ang_arbitrary = rand()*2*pi;
RB = axang2rotm([B; ang_arbitrary].');

R*A - B
RB*R*A - B

rotm2eul(R)
rotm2eul(RB*R)

Результат

ans =
   1.0e-15 *

   -0.0555
    0.1110
         0

ans =
   1.0e-15 *

    0.2220
    0.7772
   -0.2776

ans =    
    1.2220    0.3483    0.2452

ans =    
    1.2220    0.3483    0.7549
person jodag    schedule 27.07.2018
comment
Большое спасибо за Ваш ответ. Я понимаю, что в этом примере вектор, вращающийся вокруг себя, не является проблемой, но если у меня есть цепочка векторов, связанных друг с другом, где их преобразование зависит от предыдущего, у меня будет проблема, как вы проиллюстрировали выше. Правильно ли игнорировать вращение по X? Поскольку перекрестное произведение говорит мне, что нет взаимодействия между x и другими переменными? - person jcampos; 28.07.2018
comment
Я не уверен, что понимаю, что вы имеете в виду под их преобразованием в зависимости от предыдущего. Пример, который вы опубликовали, является особым случаем. Фактически, если вы измените A на [0; 1; 0], вы обнаружите, что ось вращения имеет ноль в компоненте Y, однако все 3 угла Эйлера будут различаться между R и RB*R. Тот факт, что мы наблюдаем только изменение оси X в вашем примере, является результатом соглашения, выбранного для углов Эйлера (в вашем случае XYZ). - person jodag; 28.07.2018
comment
@jcampos В случае, если ответ все еще непонятен. Полученные вами углы Эйлера будут работать в том смысле, что они повернут точку A в точку B. Они просто отличаются от тех, которые вы ожидаете, потому что есть несколько наборов углов Эйлера, которые будут работать. - person jodag; 28.07.2018

Я дам вам решение, основанное на теореме Эйлера о вращении.

Это решение дает вам только один угол, но можно вычислить и другие углы.

import numpy as np


a_vec = np.array([1, 0, 0])/np.linalg.norm(np.array([1, 0, 0]))
b_vec = np.array([0.32, 0.88, -0.34])/np.linalg.norm(np.array([0.32, 0.88, -0.34]))

cross = np.cross(a_vec, b_vec)
ab_angle = np.arccos(np.dot(a_vec,b_vec))


vx = np.array([[0,-cross[2],cross[1]],[cross[2],0,-cross[0]],[-cross[1],cross[0],0]])
R = np.identity(3)*np.cos(ab_angle) + (1-np.cos(ab_angle))*np.outer(cross,cross) + np.sin(ab_angle)*vx


validation=np.matmul(R,a_vec)

При этом в качестве векторного произведения используется общая ось вращения (в данном случае собственный вектор).

Тогда матрица R представляет собой матрицу вращения.

Это общий и очень простой способ сделать это.

person CupinaCoffee    schedule 27.07.2018