在四元数旋转上应用镜像效果的有效方法？
Posted
技术标签:
【中文标题】在四元数旋转上应用镜像效果的有效方法？【英文标题】：Efficient way to apply mirror effect on quaternion rotation? 【发布时间】：2015-12-02 23:38:13 【问题描述】：四元数表示旋转 - 它们不包含有关缩放或镜像的信息。但是，仍然可以镜像旋转的效果。
考虑在 x-y 平面上的镜像（我们也可以将其称为沿 z 轴的镜像）。在 x-y 平面上镜像的绕 x 轴的旋转将被否定。同样，围绕 y 轴旋转。但是，围绕 z 轴的旋转将保持不变。
另一个示例：围绕轴 (1,1,1) 旋转 90º，镜像在 x-y 平面中，将产生围绕 (1,1,-1) 的 -90º 旋转。为了帮助直观，如果您可以可视化轴的描述和指示旋转的圆形箭头，则镜像该可视化指示新的旋转应该是什么。
我找到了一种计算旋转镜像的方法，如下所示：
获取四元数的角轴表示。 对于每个轴 x、y 和 z。 如果沿该轴缩放为负（镜像）： 对角度和轴求反。 从修改后的角度和轴获取更新的四元数。这仅支持沿主轴 x、y 和 z 进行镜像，因为这就是我所需要的。但它适用于任意旋转。
但是，从四元数到角轴以及从角轴到四元数的转换非常昂贵。我想知道是否有办法直接在四元数本身上进行转换，但我对四元数数学的理解不足以让我自己去任何地方。
（由于计算高效方法的重要性，发布在 *** 而非数学相关论坛上。）
【问题讨论】：
我有点不确定你所说的镜像是什么意思？假设我们绕 (1,1,1) 顺时针旋转 90º（看向原点），在 x-y 平面中镜像会给出什么结果？ 在 x-y 平面镜像的轴 (1,1,1) 旋转 90º 将围绕 (1,1,-1) 旋转 -90º。为了帮助直觉，如果您可以可视化轴的描述和指示旋转的圆形箭头，则镜像该可视化指示新的旋转应该是什么。 （编辑问题以包含此内容。） 【参考方案1】：我刚刚花了相当长的时间来弄清楚这个问题的明确答案，所以我把它贴在这里记录一下。
简介
正如在其他答案中所述，镜像效果不能表示为旋转。然而，给定从坐标系 C1 到坐标系 C2 的旋转 R1to2，我们可能有兴趣在对 C1 和 C2 应用相同的镜像效果时有效地计算等效旋转（例如，我面对将在左手坐标系中给定的输入四元数转换为表示相同旋转但在右手坐标系中的四元数的问题。
就旋转矩阵而言，可以这样认为：
R_mirroredC1_to_mirroredC2 = M_mirrorC2 * R_C1_to_C2 * M_mirrorC1
这里，R_C1_to_C2 和R_mirroredC1_to_mirroredC2 都表示有效的旋转，所以在处理四元数时，如何有效地从q_C1_to_C2 计算q_mirroredC1_to_mirroredC2？
解决方案
以下假设q_C1_to_C2=[w,x,y,z]：
M_mirrorC1=M_mirrorC2=diag_3x3(-1,1,1)），则q_mirroredC1_to_mirroredC2=[w,x,-y,-z]
如果 C1 和 C2 沿 Y 轴镜像（即M_mirrorC1=M_mirrorC2=diag_3x3(1,-1,1)），则q_mirroredC1_to_mirroredC2=[w,-x,y,-z]
如果 C1 和 C2 沿 Z 轴镜像（即M_mirrorC1=M_mirrorC2=diag_3x3(1,1,-1)），则q_mirroredC1_to_mirroredC2=[w,-x,-y,z]
当考虑 C1 和 C2 的不同镜像轴时，我们有以下几点：
如果 C1 沿 X 轴镜像，C2 沿 Y 轴镜像（即M_mirrorC1=diag_3x3(-1,1,1) 和 M_mirrorC2=diag_3x3(1,-1,1)），则 q_mirroredC1_to_mirroredC2=[z,y,x,w]
如果 C1 沿 X 轴镜像，C2 沿 Z 轴镜像（即M_mirrorC1=diag_3x3(-1,1,1) 和M_mirrorC2=diag_3x3(1,1,-1)），则q_mirroredC1_to_mirroredC2=[-y,z,-w,x]
如果 C1 沿 Y 轴镜像，C2 沿 X 轴镜像（即 M_mirrorC1=diag_3x3(1,-1,1) 和 M_mirrorC2=diag_3x3(-1,1,1)），则 q_mirroredC1_to_mirroredC2=[z,-y,-x,w]
如果 C1 沿 Y 轴镜像，C2 沿 Z 轴镜像（即 M_mirrorC1=diag_3x3(1,-1,1) 和 M_mirrorC2=diag_3x3(1,1,-1)），则 q_mirroredC1_to_mirroredC2=[x,w,z,y]
如果 C1 沿 Z 轴镜像，C2 沿 X 轴镜像（即 M_mirrorC1=diag_3x3(1,1,-1) 和 M_mirrorC2=diag_3x3(-1,1,1)），则 q_mirroredC1_to_mirroredC2=[y,z,w,x]
M_mirrorC1=diag_3x3(1,1,-1) 和 M_mirrorC2=diag_3x3(1,-1,1)），则 q_mirroredC1_to_mirroredC2=[x,w,-z,-y]
测试程序
这是一个基于 OpenCV 的小型 c++ 程序来测试所有这些：
#include <opencv2/opencv.hpp>
#define CST_PI 3.1415926535897932384626433832795
// Random rotation matrix uniformly sampled from SO3 (see "Fast random rotation matrices" by J.Arvo)
cv::Matx<double,3,3> get_random_rotmat()
double theta1 = 2*CST_PI*cv::randu<double>();
double theta2 = 2*CST_PI*cv::randu<double>();
double x3 = cv::randu<double>();
cv::Matx<double,3,3> R(std::cos(theta1),std::sin(theta1),0,-std::sin(theta1),std::cos(theta1),0,0,0,1);
cv::Matx<double,3,1> v(std::cos(theta2)*std::sqrt(x3),std::sin(theta2)*std::sqrt(x3),std::sqrt(1-x3));
return -1*(cv::Matx<double,3,3>::eye()-2*v*v.t())*R;
cv::Matx<double,4,1> rotmat2quatwxyz(const cv::Matx<double,3,3> &R)
// Implementation from Ceres 1.10
const double trace = R(0,0) + R(1,1) + R(2,2);
cv::Matx<double,4,1> quat_wxyz;
if (trace >= 0.0)
double t = sqrt(trace + 1.0);
quat_wxyz(0) = 0.5 * t;
t = 0.5 / t;
quat_wxyz(1) = (R(2,1) - R(1,2)) * t;
quat_wxyz(2) = (R(0,2) - R(2,0)) * t;
quat_wxyz(3) = (R(1,0) - R(0,1)) * t;
else
int i = 0;
if (R(1, 1) > R(0, 0))
i = 1;
if (R(2, 2) > R(i, i))
i = 2;
const int j = (i + 1) % 3;
const int k = (j + 1) % 3;
double t = sqrt(R(i, i) - R(j, j) - R(k, k) + 1.0);
quat_wxyz(i + 1) = 0.5 * t;
t = 0.5 / t;
quat_wxyz(0) = (R(k,j) - R(j,k)) * t;
quat_wxyz(j + 1) = (R(j,i) + R(i,j)) * t;
quat_wxyz(k + 1) = (R(k,i) + R(i,k)) * t;
// Check that the w element is positive
if(quat_wxyz(0)<0)
quat_wxyz *= -1; // quat and -quat represent the same rotation, but to make quaternion comparison easier, we always use the one with positive w
return quat_wxyz;
cv::Matx<double,4,1> apply_quaternion_trick(const unsigned int item_permuts[4], const int sign_flips[4], const cv::Matx<double,4,1>& quat_wxyz)
// Flip the sign of the x and z components
cv::Matx<double,4,1> quat_flipped(sign_flips[0]*quat_wxyz(item_permuts[0]),sign_flips[1]*quat_wxyz(item_permuts[1]),sign_flips[2]*quat_wxyz(item_permuts[2]),sign_flips[3]*quat_wxyz(item_permuts[3]));
// Check that the w element is positive
if(quat_flipped(0)<0)
quat_flipped *= -1; // quat and -quat represent the same rotation, but to make quaternion comparison easier, we always use the one with positive w
return quat_flipped;
void detect_quaternion_trick(const cv::Matx<double,4,1> &quat_regular, const cv::Matx<double,4,1> &quat_flipped, unsigned int item_permuts[4], int sign_flips[4])
if(abs(quat_regular(0))==abs(quat_flipped(0)))
item_permuts[0]=0;
sign_flips[0] = (quat_regular(0)/quat_flipped(0)>0 ? 1 : -1);
else if(abs(quat_regular(0))==abs(quat_flipped(1)))
item_permuts[1]=0;
sign_flips[1] = (quat_regular(0)/quat_flipped(1)>0 ? 1 : -1);
else if(abs(quat_regular(0))==abs(quat_flipped(2)))
item_permuts[2]=0;
sign_flips[2] = (quat_regular(0)/quat_flipped(2)>0 ? 1 : -1);
else if(abs(quat_regular(0))==abs(quat_flipped(3)))
item_permuts[3]=0;
sign_flips[3] = (quat_regular(0)/quat_flipped(3)>0 ? 1 : -1);
if(abs(quat_regular(1))==abs(quat_flipped(0)))
item_permuts[0]=1;
sign_flips[0] = (quat_regular(1)/quat_flipped(0)>0 ? 1 : -1);
else if(abs(quat_regular(1))==abs(quat_flipped(1)))
item_permuts[1]=1;
sign_flips[1] = (quat_regular(1)/quat_flipped(1)>0 ? 1 : -1);
else if(abs(quat_regular(1))==abs(quat_flipped(2)))
item_permuts[2]=1;
sign_flips[2] = (quat_regular(1)/quat_flipped(2)>0 ? 1 : -1);
else if(abs(quat_regular(1))==abs(quat_flipped(3)))
item_permuts[3]=1;
sign_flips[3] = (quat_regular(1)/quat_flipped(3)>0 ? 1 : -1);
if(abs(quat_regular(2))==abs(quat_flipped(0)))
item_permuts[0]=2;
sign_flips[0] = (quat_regular(2)/quat_flipped(0)>0 ? 1 : -1);
else if(abs(quat_regular(2))==abs(quat_flipped(1)))
item_permuts[1]=2;
sign_flips[1] = (quat_regular(2)/quat_flipped(1)>0 ? 1 : -1);
else if(abs(quat_regular(2))==abs(quat_flipped(2)))
item_permuts[2]=2;
sign_flips[2] = (quat_regular(2)/quat_flipped(2)>0 ? 1 : -1);
else if(abs(quat_regular(2))==abs(quat_flipped(3)))
item_permuts[3]=2;
sign_flips[3] = (quat_regular(2)/quat_flipped(3)>0 ? 1 : -1);
if(abs(quat_regular(3))==abs(quat_flipped(0)))
item_permuts[0]=3;
sign_flips[0] = (quat_regular(3)/quat_flipped(0)>0 ? 1 : -1);
else if(abs(quat_regular(3))==abs(quat_flipped(1)))
item_permuts[1]=3;
sign_flips[1] = (quat_regular(3)/quat_flipped(1)>0 ? 1 : -1);
else if(abs(quat_regular(3))==abs(quat_flipped(2)))
item_permuts[2]=3;
sign_flips[2] = (quat_regular(3)/quat_flipped(2)>0 ? 1 : -1);
else if(abs(quat_regular(3))==abs(quat_flipped(3)))
item_permuts[3]=3;
sign_flips[3] = (quat_regular(3)/quat_flipped(3)>0 ? 1 : -1);
int main(int argc, char **argv)
cv::Matx<double,3,3> M_xflip(-1,0,0,0,1,0,0,0,1);
cv::Matx<double,3,3> M_yflip(1,0,0,0,-1,0,0,0,1);
cv::Matx<double,3,3> M_zflip(1,0,0,0,1,0,0,0,-1);
// Let the user choose the configuration
char im,om;
std::cout << "Enter the axis (x,y,z) along which input ref is flipped:" << std::endl;
std::cin >> im;
std::cout << "Enter the axis (x,y,z) along which output ref is flipped:" << std::endl;
std::cin >> om;
cv::Matx<double,3,3> M_iflip,M_oflip;
if(im=='x') M_iflip=M_xflip;
else if(im=='y') M_iflip=M_yflip;
else if(im=='z') M_iflip=M_zflip;
if(om=='x') M_oflip=M_xflip;
else if(om=='y') M_oflip=M_yflip;
else if(om=='z') M_oflip=M_zflip;
// Generate random quaternions until we find one where no two elements are equal
cv::Matx<double,3,3> R;
cv::Matx<double,4,1> quat_regular,quat_flipped;
do
R = get_random_rotmat();
quat_regular = rotmat2quatwxyz(R);
while(quat_regular(0)==quat_regular(1) || quat_regular(0)==quat_regular(2) || quat_regular(0)==quat_regular(3) ||
quat_regular(1)==quat_regular(2) || quat_regular(1)==quat_regular(3) ||
quat_regular(2)==quat_regular(3));
// Determine and display the appropriate quaternion trick
quat_flipped = rotmat2quatwxyz(M_oflip*R*M_iflip);
unsigned int item_permuts[4]=0,1,2,3;
int sign_flips[4]=1,1,1,1;
detect_quaternion_trick(quat_regular,quat_flipped,item_permuts,sign_flips);
char str_quat[4]='w','x','y','z';
std::cout << std::endl << "When iref is flipped along the " << im << "-axis and oref along the " << om << "-axis:" << std::endl;
std::cout << "resulting_quat=[" << (sign_flips[0]>0?"":"-") << str_quat[item_permuts[0]] << ","
<< (sign_flips[1]>0?"":"-") << str_quat[item_permuts[1]] << ","
<< (sign_flips[2]>0?"":"-") << str_quat[item_permuts[2]] << ","
<< (sign_flips[3]>0?"":"-") << str_quat[item_permuts[3]] << "], where initial_quat=[w,x,y,z]" << std::endl;
// Test this trick on several random rotation matrices
unsigned int n_errors = 0, n_tests = 10000;
std::cout << std::endl << "Performing " << n_tests << " tests on random rotation matrices:" << std::endl;
for(unsigned int i=0; i<n_tests; ++i)
// Get a random rotation matrix and the corresponding quaternion
cv::Matx<double,3,3> R = get_random_rotmat();
cv::Matx<double,4,1> quat_regular = rotmat2quatwxyz(R);
// Get the quaternion corresponding to the flipped coordinate frames, via the sign trick and via computation on rotation matrices
cv::Matx<double,4,1> quat_tricked = apply_quaternion_trick(item_permuts,sign_flips,quat_regular);
cv::Matx<double,4,1> quat_flipped = rotmat2quatwxyz(M_oflip*R*M_iflip);
// Check that both results are identical
if(cv::norm(quat_tricked-quat_flipped,cv::NORM_INF)>1e-6)
std::cout << "Error (idx=" << i << ")!"
<< "\n quat_regular=" << quat_regular.t()
<< "\n quat_tricked=" << quat_tricked.t()
<< "\n quat_flipped=" << quat_flipped.t() << std::endl;
++n_errors;
std::cout << n_errors << " errors on " << n_tests << " tests." << std::endl;
system("pause");
return 0;
【讨论】：
【参考方案2】：我做了一些进一步的分析，似乎四元数 (w, x, y, z) 的效果可以像这样反映它的效果：
通过翻转四元数的 y 和 z 元素实现沿 x 轴旋转的镜像效果。 通过翻转四元数的 x 和 z 元素实现沿 y 轴旋转的镜像效果。 通过翻转四元数的 x 和 y 元素实现沿 z 轴旋转的镜像效果。永远不需要触摸四元数的 w 元素。
不幸的是，我仍然不能很好地理解四元数，无法解释为什么会这样，但我是从轴角格式转换的实现中推导出来的，在实现了这个解决方案之后，它的工作原理和我对它进行的所有测试中的原始测试。
【讨论】：
您没有碰巧在 y 轴上遇到任何问题吧？正如你所说，我可以沿 x 轴和 z 轴镜像，但由于某种原因，当我尝试否定x 和 z 时，它不是沿 y 轴镜像，而是相当于沿 x 轴镜像，然后是 z 轴。【参考方案3】：
有一点更容易和面向程序员的方式来思考这个问题。假设您想在坐标系中反转 z 轴（即将 z 轴翻转到 -z）。现在将四元数视为滚动、俯仰和偏航方面的方向向量。当您翻转 z 轴时，请注意滚动和俯仰的符号反转但偏航的符号保持不变。
现在您可以使用以下将欧拉角转换为四元数的代码找到对四元数的净影响（我会将这段代码放到***上）：
static Quaterniond toQuaternion(double pitch, double roll, double yaw)
Quaterniond q;
double t0 = std::cos(yaw * 0.5f);
double t1 = std::sin(yaw * 0.5f);
double t2 = std::cos(roll * 0.5f);
double t3 = std::sin(roll * 0.5f);
double t4 = std::cos(pitch * 0.5f);
double t5 = std::sin(pitch * 0.5f);
q.w() = t0 * t2 * t4 + t1 * t3 * t5;
q.x() = t0 * t3 * t4 - t1 * t2 * t5;
q.y() = t0 * t2 * t5 + t1 * t3 * t4;
q.z() = t1 * t2 * t4 - t0 * t3 * t5;
return q;
使用基本三角函数sin(-x) = -sin(x) 和cos(-x) = cos(x)。将此应用于上面的代码，您可以看到 t3 和 t5 的符号将翻转。这将导致 x 和 y 的符号翻转。
所以当你反转 z 轴时，
Q'(w, x, y, z) = Q(w, -x, -y, z)
同样，您可以找出轴反转的任何其他组合并找出对四元数的影响。
PS：如果有人想知道为什么有人需要这个...我需要在上面转换来自控制无人机的 MavLink/Pixhawk 系统的四元数。源系统使用 NED 坐标系，但通常的 3D 环境（如 Unreal）使用 NEU 坐标系，需要将 z 轴转换为 -z 才能正确使用四元数。
【讨论】：
【参考方案4】：我们可以在 3D 中检查所有旋转和反射的集合，这称为Orthogonal group O(3)。它可以看作是具有行列式 +1 或 -1 的正交矩阵的集合。所有旋转的行列式为 +1，纯反射的行列式为 -1。 O(3) 的另一个成员在点 (x,y,z)->(-x,-y,-z) 中的反演在 3D 中具有 det -1，我们稍后会讨论。如果我们将组中的两个变换结合起来，则将它们的行列式相乘。因此，两个旋转组合产生另一个旋转 (+1 * +1 = +1)，一个旋转与一个反射组合产生一个反射 (+1 * -1 = -1)，两个反射组合产生一个旋转 (-1 * -1 = +1)。
我们可以将 O(3) 限制为仅具有行列式 +1 的那些以形成 Special Orthogonal Group SO(3)。这仅包含旋转。
现在单位四元数的集合是 SO(3) 的双重覆盖，这意味着两个单位四元数对应于每个旋转。准确地说，如果a+b i+c j+d k 是一个单位四元数，那么a-b i-c j-d k 表示相同的旋转，您可以将其视为围绕向量 (b,c,d) 旋转 ø 与围绕 -ø 旋转相同向量 (-b,-c,-d)。
请注意，所有单位四元数都有行列式 +1，因此没有一个对应于纯反射。这就是为什么不能使用四元数来表示反射的原因。
你可以做的是使用倒置。现在，在倒置之后的反射是旋转。例如，在 x=0 中反射并反转，与在 y=0 中反射和在 z=0 中反射相同。这与绕 x 轴旋转 180º 相同。您可以对任何反射执行相同的过程。
我们可以通过使用法向量n = (a,b,c)来定义一个通过原点的平面。该平面中向量 v(x,y,z) 的 reflection 由
给出v - 2 (v . n ) / ( n . n) n
= (x,y,z) - 2 (a x+b y+c z) / (a^2+b^2+c^2) (a,b,c)
特别是 x-y 平面具有法线 (0,0,1)，因此反射是
(x,y,z) - 2 z (0,0,1) = (x,y,-z)
Quaternions and spatial rotation 有一个很好的轴角公式中的四元数公式。
p = cos(ø/2) + (x i + y j + z k) sin(ø/2)
这是一个四元数 W + X i + Y j + Z k W=cos(ø/2), X = x sin(ø/2), Y = y sin(ø/2), Z = z sin(ø/2)
改变旋转方向将翻转半角的sin但保持cos不变，给出
p' = cos(ø/2) - (x i + y j + z k) sin(ø/2)
现在，如果我们考虑在 x-y 平面上反射相应的向量给出
q = cos(ø/2) + (x i + y j - z k) sin(ø/2)
我们可能想改变旋转方向
q' = cos(ø/2) + (- x i - y j + z k) sin(ø/2 )
= W - X i - Y j + Z k
我认为这与您的答案相对应。
我们可以将其推广到具有单位长度法线 (a,b,c) 的一般平面中的反射。令 d 为点积 (a,b,c).(x,y,z)。 (x,y,z)的反射是
(x,y,z) - 2 d (a,b,c) = (x - 2 d a, y - 2 d b, z - 2 d c)
这个的旋转四元数
q = cos(ø/2) - ((x - 2 da) i + ((y - 2 db) j + (z - 2 dc) k) sin(ø/2)
q = cos(ø/2) - (x i + y j + z k) sin(ø/2) + 2 d sin(ø/2) (a i + b j + c k)
= W - X i - Y j - Z k + 2 d (X,Y,Z).(a,b ,c) (a i + b j + c k)
【讨论】：
我不明白如何根据这个描述进行实现。您所描述的是我已经使用角度轴格式的旋转所做的事情。我不明白的是如何直接在四元数上做同样的事情，而不需要转换成角度轴格式。给定一个包含成员 a、b、c 和 d 的四元数，执行您描述的操作的伪代码是什么样的？【参考方案5】：请注意，镜像不是旋转，因此通常您不能将其烘焙成四元数（不过，我很可能误解了您的问题）。镜像变换矩阵的 3x3 分量是
M = I-2(n*nT)
其中 I 是单位 3x3 矩阵，n 是镜像平面的法线，表示为 3x1 矩阵，nT 是 n 作为 1x3 矩阵（因此 n*nT 是 3x(1x1)x3=3x3 矩阵）。
现在，如果您要“镜像”的四元数 q 是最后一个变换，那么另一侧的最后一个变换将只是 M*q（同样，这将是一个通用的 3x3 矩阵，通常不能表示为四元数)
【讨论】：
我知道镜像不能用旋转来表示，也不是我想要的镜像变换。我知道这很容易引起混淆，这就是为什么我花了前两段试图尽可能清楚地说明这一点。我想它仍然不够清楚，但不幸的是，我不知道如何使它更清楚。 如果您描述一下您需要它的用途，也许会有所帮助？ 在 Unity 游戏引擎中，场景支持对象的树状层次结构。可以根据平移、旋转和指定沿 3 个局部轴中的每一个的缩放的矢量，相对于其父级变换子级。不支持完整的 3x3 或 4x4 矩阵变换。问题是如何处理一个常见的情况，即父级由于例如比例尺而沿一个或多个轴镜像。 并且仍然以合理的方式处理子对象的旋转。我已经找到了解决问题的方法，所以我不会花时间讨论更多细节。【参考方案6】：对于任何通过网络搜索来到这里并正在寻找数学的人，那么：
反射
通过平面 ax+by+cz=0 反射点 'p'，使用四元数：
n = 0+(a,b,c) p = 0+(x,y,z)
其中 'n' 是单位双向量（如果您愿意，也可以是纯四元数）
p' = npn
那么 p' 是反射点。
如果你用第二个反射“m”作曲：
p' = mnpnm = (mn)p(mn)^*
是一个旋转。
旋转和反射按预期组合。
统一缩放
由于标量积可以交换并且可以被分解，那么如果我们有单位四元数'Q'的旋转或单位双向量'b'（或任何组合）的反射乘以某个非零比例值' s' 导致 s^2 的统一缩放。并且由于 (sqrt(s0)*sqrt(s1))^2 = s0*s1，这些统一缩放值按预期组成。
然而这一点可能并不重要，因为在代码中我们希望能够假设单位幅度值以降低运行时复杂性。
【讨论】：
以上是关于在四元数旋转上应用镜像效果的有效方法？的主要内容，如果未能解决你的问题，请参考以下文章