使用 GLM 进行基于四元数的点旋转
Quaternion based point rotations using GLM
我正在尝试使用在 GLM 中实现的四元数来旋转一个点。最终目标是使用此代码创建轨道相机,但这是帮助理解代码背后动机的旁注。
为了更好地理解基于四元数的旋转,我编写了一些包含两个循环的代码。第一个循环将通过围绕 X 轴一直旋转到 90 度来逐步改变四元数的方向,第二个循环将继续围绕 Z 轴一直旋转到 90 度。每个循环执行 4 个步骤。因此每个循环都围绕各自的轴递增旋转 90 / 4 = 22.5 度。使用四元数乘法应用方向的变化,并使用欧拉角跟踪。循环应以四元数结束,该四元数会将 (0, 0, 3) 处的点旋转到 (3, 0, 0)。请注意,我不只是试图确定将执行此旋转的四元数。目标是执行一系列增量旋转。
如果我们看下图,从 C 到 I 的转换发生在第一个循环中,然后从 I 到 R 的转换发生在第二个循环中(请原谅稀疏点命名)。
v' = q * v * q^-1
其中 v 应被视为纯四元数(标量项 w 为零),q 需要是单位四元数(长度为 1)。根据我的理解,右手乘法与四元数的倒数需要将结果 v' 保持在 3D space 中,而不是以 4D 向量结束。所以v'也需要是一个纯四元数。
然后是旋转的加倍效果,其中左手与 q 的乘法贡献了所需旋转的一半,而右手与逆相乘则增加了所需旋转的另一半。
Ben Eater 和 Grant Sanderson 对四元数进行了出色的交互式可视化和解释,我将其用作交叉参考。可以查到here.
所以我们首先需要使用一个绕X轴旋转11.25度的四元数和GLM returns这个四元数求欧拉角(使用四元数表示法[w, [x, y, z]] ):
Rotation of [ 11.25, 0.00, 0.00] deg => Q: [ 0.9952, [ 0.0980, 0.0000, 0.0000]]
根据this,由于我们纯粹绕X轴旋转,我们可以通过对四元数的w分量执行acos来验证GLM计算的四元数中的旋转量:
float angle = acosf(q.w)
然后:
acos(0.9952) = 0.0980 rad / 5.6 degrees
这是所需角度的一半...这也在交互式动画的交叉检查中得到证实(请原谅四舍五入):
因此 GLM 返回的 11.25 度四元数实际上旋转了所需角度的一半...如果我们查看 GLM 代码,从欧拉角计算 w 分量会稍微复杂一些,因为旋转可能发生在周围任意旋转轴...但是欧拉角明显减半:
template <typename T, precision P>
GLM_FUNC_QUALIFIER tquat<T, P>::tquat(tvec3<T, P> const & eulerAngle)
{
tvec3<T, P> c = glm::cos(eulerAngle * T(0.5));
tvec3<T, P> s = glm::sin(eulerAngle * T(0.5));
this->w = c.x * c.y * c.z + s.x * s.y * s.z;
this->x = s.x * c.y * c.z - c.x * s.y * s.z;
this->y = c.x * s.y * c.z + s.x * c.y * s.z;
this->z = c.x * c.y * s.z - s.x * s.y * c.z;
}
我的第一个问题是为什么 GLM 将角度减半?
尽管所需的旋转角度不同,但我还是继续检查了两个循环的旋转结果。结果是……出乎意料。
如果我使用“不正确的形式”旋转(一些 OpenGL 在线教程建议)并且仅通过左手乘法旋转点(但整步为 22.5 度):
v' = q * v
我得到了我希望的结果。重点是正确地遵循所有中间步骤并从 (0, 0, 3) 到 (3, 0, 0)。此外,w 分量在所有中间步骤中均为 0。
但是如果我使用旋转的“正确形式”并通过左手乘以 q 和右手乘以 q 的倒数来旋转点(对于 11.25 度的半步来解释加倍旋转):
v' = q * v * q^-1
当第二个循环开始围绕 Z 轴旋转点时,我就开始得到错误的结果。一个小而明显的 Z 分量开始悄悄进入,旋转距离 22.5 度的完整步长仅一步之遥。这在下图中的绿点中可见。
对于两种旋转方法,旋转点的 w 分量都保持为 0...
谁能解释为什么 GLM 旋转可以正确地从左侧进行一次乘法运算?
这是为了将操作次数减少到最少的某种优化吗?
我可以在 GLM 中使用 v' = q * v
旋转来获得所有旋转的一致且正确的结果吗?
代码:
const int rotSteps = 4;
// Rotate around X axis in steps to 90deg
vec3 eulerState = vec3(0.0f);
// point we want to rotate (use vec4 to track the w component during rotations)
vec4 v = vec4(0.0f, 0.0f, 3.0f, 0.0f);
// Full Euler steps for q * v rotation
quat orientF = quat(1.0f, 0.0f, 0.0f, 0.0f);
vec3 euler = vec3(RAD(90.0f), RAD(0.0f), RAD(0.0f));
vec3 eulerStep = euler / (float)rotSteps;
quat qEulerF = quat(eulerStep); // GetRotQuat(eulerStep);
vec4 qa = ToAngularForm(qEulerF);
vec3 orientEuler = eulerAngles(qEulerF);
CLogD(TAG, "Rot Full Step Q [W, X, Y, Z]: " FMT_Q(4) " / " FMT_V3(2) "deg / " FMT_QA(2), PAR_Q(qEulerF), PAR_V3(degrees(orientEuler)), PAR_QA(qa));
// Half Euler steps for q * v * q^-1 rotation
quat orientH = quat(1.0f, 0.0f, 0.0f, 0.0f);
vec3 eulerStepH = eulerStep / 2.0f;
quat qEulerH = quat(eulerStepH); // GetRotQuat(eulerStepH);
qa = ToAngularForm(qEulerH);
orientEuler = eulerAngles(qEulerH);
CLogD(TAG, "Rot Half Step Q [W, X, Y, Z]: " FMT_Q(4) " / " FMT_V3(2) "deg / " FMT_QA(2), PAR_Q(qEulerH), PAR_V3(degrees(orientEuler)), PAR_QA(qa));
quat qEulerHI = inverse(qEulerH);
vec4 qai = ToAngularForm(qEulerHI);
orientEuler = eulerAngles(qEulerHI);
CLogD(TAG, "Rot Half Step Q^-1 [W, X, Y, Z]: " FMT_Q(4) " / " FMT_V3(2) "deg / " FMT_QA(2), PAR_Q(qEulerHI), PAR_V3(degrees(orientEuler)), PAR_QA(qai));
for (int rotStep = 1; rotStep <= rotSteps; ++rotStep)
{
// Track the absolute Euler rotation
eulerState += eulerStep;
// Rotate by incremental rotation as defined by Euler angles
orientH = qEulerH * orientH;
orientEuler = eulerAngles(orientH);
CLogI(TAG, "Rot Step %d. Curr Abs Q: " FMT_Q(4) "/" FMT_V3(2) "deg, Abs Euler: " FMT_V3(2) "deg",
rotStep, PAR_Q(orientH), PAR_V3(degrees(orientEuler)), PAR_V3(degrees(eulerState)));
// Transform the point using the correct q * v * q^-1 rotation and multiply from Left and Right
quat orientHI = inverse(orientH);
qa = ToAngularForm(orientH);
qai = ToAngularForm(orientHI);
vec4 rotV = orientH * v * orientHI;
CLogD(TAG, "Rot QL: " FMT_Q(4) " / " FMT_QA(1), PAR_Q(orientH), PAR_QA(qa));
CLogD(TAG, "Rot QR: " FMT_Q(4) " / " FMT_QA(1), PAR_Q(orientHI), PAR_QA(qai));
CLogD(TAG, "Rot LR -> " FMT_V4(1), PAR_V4(rotV));
// Transform the point using the incorrect q * v rotation and multiply from Left only
orientF = qEulerF * orientF;
qa = ToAngularForm(orientF);
rotV = orientF * v;
CLogD(TAG, "Rot QR: " FMT_Q(4) " / " FMT_QA(1), PAR_Q(orientF), PAR_QA(qa));
CLogD(TAG, "Rot L -> " FMT_V4(1), PAR_V4(rotV));
}
// Rotate for 90 degrees around the Z axis
// Full Euler steps for q * v rotation
euler = vec3(RAD(0.0f), RAD(0.0f), RAD(90.0f));
eulerStep = euler / (float)rotSteps;
qEulerF = quat(eulerStep); // GetRotQuat(eulerStep);
qa = ToAngularForm(qEulerF);
orientEuler = eulerAngles(qEulerF);
CLogD(TAG, "Rot Full Step Q [W, X, Y, Z]: " FMT_Q(4) " / " FMT_V3(2) "deg / " FMT_QA(2), PAR_Q(qEulerF), PAR_V3(degrees(orientEuler)), PAR_QA(qa));
// Half Euler steps for q * v * q^-1 rotation
eulerStepH = eulerStep / 2.0f;
qEulerH = quat(eulerStepH); // GetRotQuat(eulerStepH);
qa = ToAngularForm(qEulerH);
orientEuler = eulerAngles(qEulerH);
CLogD(TAG, "Rot Half Step Q [W, X, Y, Z]: " FMT_Q(4) " / " FMT_V3(2) "deg / " FMT_QA(2), PAR_Q(qEulerH), PAR_V3(degrees(orientEuler)), PAR_QA(qa));
qEulerHI = inverse(qEulerH);
qai = ToAngularForm(qEulerHI);
orientEuler = eulerAngles(qEulerHI);
CLogD(TAG, "Rot Half Step Q^-1 [W, X, Y, Z]: " FMT_Q(4) " / " FMT_V3(2) "deg / " FMT_QA(2), PAR_Q(qEulerHI), PAR_V3(degrees(orientEuler)), PAR_QA(qai));
for (int rotStep = 1; rotStep <= rotSteps; ++rotStep)
{
// Track the absolute Euler rotation
eulerState += eulerStep;
// Rotate by incremental rotation as defined by Euler angles
orientH = qEulerH * orientH;
orientEuler = eulerAngles(orientH);
CLogI(TAG, "Rot Step %d. Curr Abs Q: " FMT_Q(4) "/" FMT_V3(2) "deg, Abs Euler: " FMT_V3(2) "deg",
rotStep, PAR_Q(orientH), PAR_V3(degrees(orientEuler)), PAR_V3(degrees(eulerState)));
// Transform the point using the correct q * v * q^-1 rotation and multiply from Left and Right
quat orientHI = inverse(orientH);
qa = ToAngularForm(orientH);
qai = ToAngularForm(orientHI);
vec4 rotV = orientH * v * orientHI;
CLogD(TAG, "Rot QL: " FMT_Q(4) " / " FMT_QA(1), PAR_Q(orientH), PAR_QA(qa));
CLogD(TAG, "Rot QR: " FMT_Q(4) " / " FMT_QA(1), PAR_Q(orientHI), PAR_QA(qai));
CLogD(TAG, "Rot LR -> " FMT_V4(1), PAR_V4(rotV));
// Transform the point using the incorrect q * v rotation and multiply from Left only
orientF = qEulerF * orientF;
qa = ToAngularForm(orientF);
rotV = orientF * v;
CLogD(TAG, "Rot QR: " FMT_Q(4) " / " FMT_QA(1), PAR_Q(orientF), PAR_QA(qa));
CLogD(TAG, "Rot L -> " FMT_V4(1), PAR_V4(rotV));
}
输出:
Rot Full Step Q [W, X, Y, Z]: [ 0.9808, [ 0.1951, 0.0000, 0.0000]] / [ 22.50, -0.00, 0.00]deg / cos( 11.25) + sin( 11.25)( 1.00i + 0.00j + 0.00k)
Rot Half Step Q [W, X, Y, Z]: [ 0.9952, [ 0.0980, 0.0000, 0.0000]] / [ 11.25, -0.00, 0.00]deg / cos( 5.63) + sin( 5.63)( 1.00i + 0.00j + 0.00k)
Rot Half Step Q^-1 [W, X, Y, Z]: [ 0.9952, [-0.0980, -0.0000, -0.0000]] / [-11.25, -0.00, 0.00]deg / cos( 5.63) + sin( 5.63)(-1.00i + -0.00j + -0.00k)
Rot Step 1. Curr Abs Q: [ 0.9952, [ 0.0980, 0.0000, 0.0000]]/[ 11.25, -0.00, 0.00]deg, Abs Euler: [ 22.50, 0.00, 0.00]deg
Rot QL: [ 0.9952, [ 0.0980, 0.0000, 0.0000]] / cos( 5.6) + sin( 5.6)( 1.0i + 0.0j + 0.0k)
Rot QR: [ 0.9952, [-0.0980, -0.0000, -0.0000]] / cos( 5.6) + sin( 5.6)(-1.0i + -0.0j + -0.0k)
Rot LR -> [ 0.0, -1.1, 2.8, 0.0]
Rot QR: [ 0.9808, [ 0.1951, 0.0000, 0.0000]] / cos( 11.3) + sin( 11.3)( 1.0i + 0.0j + 0.0k)
Rot L -> [ 0.0, -1.1, 2.8, 0.0]
Rot Step 2. Curr Abs Q: [ 0.9808, [ 0.1951, 0.0000, 0.0000]]/[ 22.50, -0.00, 0.00]deg, Abs Euler: [ 45.00, 0.00, 0.00]deg
Rot QL: [ 0.9808, [ 0.1951, 0.0000, 0.0000]] / cos( 11.3) + sin( 11.3)( 1.0i + 0.0j + 0.0k)
Rot QR: [ 0.9808, [-0.1951, -0.0000, -0.0000]] / cos( 11.2) + sin( 11.2)(-1.0i + -0.0j + -0.0k)
Rot LR -> [ 0.0, -2.1, 2.1, 0.0]
Rot QR: [ 0.9239, [ 0.3827, 0.0000, 0.0000]] / cos( 22.5) + sin( 22.5)( 1.0i + 0.0j + 0.0k)
Rot L -> [ 0.0, -2.1, 2.1, 0.0]
Rot Step 3. Curr Abs Q: [ 0.9569, [ 0.2903, 0.0000, 0.0000]]/[ 33.75, -0.00, 0.00]deg, Abs Euler: [ 67.50, 0.00, 0.00]deg
Rot QL: [ 0.9569, [ 0.2903, 0.0000, 0.0000]] / cos( 16.9) + sin( 16.9)( 1.0i + 0.0j + 0.0k)
Rot QR: [ 0.9569, [-0.2903, -0.0000, -0.0000]] / cos( 16.9) + sin( 16.9)(-1.0i + -0.0j + -0.0k)
Rot LR -> [ 0.0, -2.8, 1.1, 0.0]
Rot QR: [ 0.8315, [ 0.5556, 0.0000, 0.0000]] / cos( 33.8) + sin( 33.8)( 1.0i + 0.0j + 0.0k)
Rot L -> [ 0.0, -2.8, 1.1, 0.0]
Rot Step 4. Curr Abs Q: [ 0.9239, [ 0.3827, 0.0000, 0.0000]]/[ 45.00, -0.00, 0.00]deg, Abs Euler: [ 90.00, 0.00, 0.00]deg
Rot QL: [ 0.9239, [ 0.3827, 0.0000, 0.0000]] / cos( 22.5) + sin( 22.5)( 1.0i + 0.0j + 0.0k)
Rot QR: [ 0.9239, [-0.3827, -0.0000, -0.0000]] / cos( 22.5) + sin( 22.5)(-1.0i + -0.0j + -0.0k)
Rot LR -> [ 0.0, -3.0, 0.0, 0.0]
Rot QR: [ 0.7071, [ 0.7071, 0.0000, 0.0000]] / cos( 45.0) + sin( 45.0)( 1.0i + 0.0j + 0.0k)
Rot L -> [ 0.0, -3.0, 0.0, 0.0]
Rot Full Step Q [W, X, Y, Z]: [ 0.9808, [ 0.0000, 0.0000, 0.1951]] / [ 0.00, -0.00, 22.50]deg / cos( 11.25) + sin( 11.25)( 0.00i + 0.00j + 1.00k)
Rot Half Step Q [W, X, Y, Z]: [ 0.9952, [ 0.0000, 0.0000, 0.0980]] / [ 0.00, -0.00, 11.25]deg / cos( 5.63) + sin( 5.63)( 0.00i + 0.00j + 1.00k)
Rot Half Step Q^-1 [W, X, Y, Z]: [ 0.9952, [-0.0000, -0.0000, -0.0980]] / [ 0.00, -0.00, -11.25]deg / cos( 5.63) + sin( 5.63)(-0.00i + -0.00j + -1.00k)
Rot Step 1. Curr Abs Q: [ 0.9194, [ 0.3808, 0.0375, 0.0906]]/[ 45.00, 0.00, 11.25]deg, Abs Euler: [ 90.00, 0.00, 22.50]deg
Rot QL: [ 0.9194, [ 0.3808, 0.0375, 0.0906]] / cos( 23.2) + sin( 23.2)( 1.0i + 0.1j + 0.2k)
Rot QR: [ 0.9194, [-0.3808, -0.0375, -0.0906]] / cos( 23.2) + sin( 23.2)(-1.0i + -0.1j + -0.2k)
Rot LR -> [ 1.0, -2.8, 0.0, 0.0]
Rot QR: [ 0.6935, [ 0.6935, 0.1379, 0.1379]] / cos( 46.1) + sin( 46.1)( 1.0i + 0.2j + 0.2k)
Rot L -> [ 1.1, -2.8, 0.0, 0.0]
Rot Step 2. Curr Abs Q: [ 0.9061, [ 0.3753, 0.0747, 0.1802]]/[ 45.00, -0.00, 22.50]deg, Abs Euler: [ 90.00, 0.00, 45.00]deg
Rot QL: [ 0.9061, [ 0.3753, 0.0747, 0.1802]] / cos( 25.0) + sin( 25.0)( 0.9i + 0.2j + 0.4k)
Rot QR: [ 0.9061, [-0.3753, -0.0747, -0.1802]] / cos( 25.0) + sin( 25.0)(-0.9i + -0.2j + -0.4k)
Rot LR -> [ 1.9, -2.4, 0.1, 0.0]
Rot QR: [ 0.6533, [ 0.6533, 0.2706, 0.2706]] / cos( 49.2) + sin( 49.2)( 0.9i + 0.4j + 0.4k)
Rot L -> [ 2.1, -2.1, 0.0, 0.0]
Rot Step 3. Curr Abs Q: [ 0.8841, [ 0.3662, 0.1111, 0.2682]]/[ 45.00, 0.00, 33.75]deg, Abs Euler: [ 90.00, 0.00, 67.50]deg
Rot QL: [ 0.8841, [ 0.3662, 0.1111, 0.2682]] / cos( 27.9) + sin( 27.9)( 0.8i + 0.2j + 0.6k)
Rot QR: [ 0.8841, [-0.3662, -0.1111, -0.2682]] / cos( 27.9) + sin( 27.9)(-0.8i + -0.2j + -0.6k)
Rot LR -> [ 2.5, -1.6, 0.3, 0.0]
Rot QR: [ 0.5879, [ 0.5879, 0.3928, 0.3928]] / cos( 54.0) + sin( 54.0)( 0.7i + 0.5j + 0.5k)
Rot L -> [ 2.8, -1.1, 0.0, 0.0]
Rot Step 4. Curr Abs Q: [ 0.8536, [ 0.3536, 0.1464, 0.3536]]/[ 45.00, 0.00, 45.00]deg, Abs Euler: [ 90.00, 0.00, 90.00]deg
Rot QL: [ 0.8536, [ 0.3536, 0.1464, 0.3536]] / cos( 31.4) + sin( 31.4)( 0.7i + 0.3j + 0.7k)
Rot QR: [ 0.8536, [-0.3536, -0.1464, -0.3536]] / cos( 31.4) + sin( 31.4)(-0.7i + -0.3j + -0.7k)
Rot LR -> [ 2.9, -0.7, 0.4, 0.0]
Rot QR: [ 0.5000, [ 0.5000, 0.5000, 0.5000]] / cos( 60.0) + sin( 60.0)( 0.6i + 0.6j + 0.6k)
Rot L -> [ 3.0, 0.0, 0.0, 0.0]
我得到了我的问题的答案和一台工作正常的轨道相机,但没有时间仔细检查示例代码现在是否正常工作 - 它应该。
第一个问题是为什么 GLM 在四元数转换过程中将角度减半,根据 extended Euler's formula 看起来是这样......它必须这样做。这部分可以再研究一下,但由于时间不够,我不得不接受它。
GLM 中的矢量旋转是使用乘法运算符实现的。这意味着当将 vec3 与四元数相乘时,不会将向量转换为四元数然后执行乘法,它会执行 vector rotation instead:
template <typename T, precision P>
GLM_FUNC_QUALIFIER tvec3<T, P> operator*(tquat<T, P> const & q, tvec3<T, P> const & v)
{
tvec3<T, P> const QuatVector(q.x, q.y, q.z);
tvec3<T, P> const uv(glm::cross(QuatVector, v));
tvec3<T, P> const uuv(glm::cross(QuatVector, uv));
return v + ((uv * q.w) + uuv) * static_cast<T>(2);
}
所以,是的,使用四元数旋转向量的正确方法是在四元数和向量之间使用乘法运算符,如下所示:
v' = q * v
或者在 C++ 中:
vec3 posOrigin;
quat rotQ;
...
vec3 posRot = rotQ * posOrigin;
这段代码实际上并没有直接进行四元数乘法。它进行旋转。就个人而言,我更喜欢 GLM 提供像 rotate(quat, vec)
这样的函数调用...但我确信运算符重载混淆是有原因的。
另请注意,操作数顺序很重要,因为向量和四元数之间的乘法定义如下:
template <typename T, precision P>
GLM_FUNC_QUALIFIER tvec3<T, P> operator*(tvec3<T, P> const & v, tquat<T, P> const & q)
{
return glm::inverse(q) * v;
}
因此将反向旋转矢量。
请注意,GLM 还实现了四元数之间的乘法,但为此需要使用两个四元数之间的乘法运算符:
template <typename T, precision P>
template <typename U>
GLM_FUNC_QUALIFIER tquat<T, P> & tquat<T, P>::operator*=(tquat<U, P> const & r)
{
tquat<T, P> const p(*this);
tquat<T, P> const q(r);
this->w = p.w * q.w - p.x * q.x - p.y * q.y - p.z * q.z;
this->x = p.w * q.x + p.x * q.w + p.y * q.z - p.z * q.y;
this->y = p.w * q.y + p.y * q.w + p.z * q.x - p.x * q.z;
this->z = p.w * q.z + p.z * q.w + p.x * q.y - p.y * q.x;
return *this;
}
由于 GLM 几乎没有我能找到的宝贵文档,因此这种运算符重载会导致错误的假设和大量的时间损失。所以我想我应该阅读 GLM 代码而不是假设它的作用...
我正在尝试使用在 GLM 中实现的四元数来旋转一个点。最终目标是使用此代码创建轨道相机,但这是帮助理解代码背后动机的旁注。
为了更好地理解基于四元数的旋转,我编写了一些包含两个循环的代码。第一个循环将通过围绕 X 轴一直旋转到 90 度来逐步改变四元数的方向,第二个循环将继续围绕 Z 轴一直旋转到 90 度。每个循环执行 4 个步骤。因此每个循环都围绕各自的轴递增旋转 90 / 4 = 22.5 度。使用四元数乘法应用方向的变化,并使用欧拉角跟踪。循环应以四元数结束,该四元数会将 (0, 0, 3) 处的点旋转到 (3, 0, 0)。请注意,我不只是试图确定将执行此旋转的四元数。目标是执行一系列增量旋转。
如果我们看下图,从 C 到 I 的转换发生在第一个循环中,然后从 I 到 R 的转换发生在第二个循环中(请原谅稀疏点命名)。
v' = q * v * q^-1
其中 v 应被视为纯四元数(标量项 w 为零),q 需要是单位四元数(长度为 1)。根据我的理解,右手乘法与四元数的倒数需要将结果 v' 保持在 3D space 中,而不是以 4D 向量结束。所以v'也需要是一个纯四元数。
然后是旋转的加倍效果,其中左手与 q 的乘法贡献了所需旋转的一半,而右手与逆相乘则增加了所需旋转的另一半。
Ben Eater 和 Grant Sanderson 对四元数进行了出色的交互式可视化和解释,我将其用作交叉参考。可以查到here.
所以我们首先需要使用一个绕X轴旋转11.25度的四元数和GLM returns这个四元数求欧拉角(使用四元数表示法[w, [x, y, z]] ):
Rotation of [ 11.25, 0.00, 0.00] deg => Q: [ 0.9952, [ 0.0980, 0.0000, 0.0000]]
根据this,由于我们纯粹绕X轴旋转,我们可以通过对四元数的w分量执行acos来验证GLM计算的四元数中的旋转量:
float angle = acosf(q.w)
然后:
acos(0.9952) = 0.0980 rad / 5.6 degrees
这是所需角度的一半...这也在交互式动画的交叉检查中得到证实(请原谅四舍五入):
因此 GLM 返回的 11.25 度四元数实际上旋转了所需角度的一半...如果我们查看 GLM 代码,从欧拉角计算 w 分量会稍微复杂一些,因为旋转可能发生在周围任意旋转轴...但是欧拉角明显减半:
template <typename T, precision P>
GLM_FUNC_QUALIFIER tquat<T, P>::tquat(tvec3<T, P> const & eulerAngle)
{
tvec3<T, P> c = glm::cos(eulerAngle * T(0.5));
tvec3<T, P> s = glm::sin(eulerAngle * T(0.5));
this->w = c.x * c.y * c.z + s.x * s.y * s.z;
this->x = s.x * c.y * c.z - c.x * s.y * s.z;
this->y = c.x * s.y * c.z + s.x * c.y * s.z;
this->z = c.x * c.y * s.z - s.x * s.y * c.z;
}
我的第一个问题是为什么 GLM 将角度减半?
尽管所需的旋转角度不同,但我还是继续检查了两个循环的旋转结果。结果是……出乎意料。
如果我使用“不正确的形式”旋转(一些 OpenGL 在线教程建议)并且仅通过左手乘法旋转点(但整步为 22.5 度):
v' = q * v
我得到了我希望的结果。重点是正确地遵循所有中间步骤并从 (0, 0, 3) 到 (3, 0, 0)。此外,w 分量在所有中间步骤中均为 0。
但是如果我使用旋转的“正确形式”并通过左手乘以 q 和右手乘以 q 的倒数来旋转点(对于 11.25 度的半步来解释加倍旋转):
v' = q * v * q^-1
当第二个循环开始围绕 Z 轴旋转点时,我就开始得到错误的结果。一个小而明显的 Z 分量开始悄悄进入,旋转距离 22.5 度的完整步长仅一步之遥。这在下图中的绿点中可见。
对于两种旋转方法,旋转点的 w 分量都保持为 0...
谁能解释为什么 GLM 旋转可以正确地从左侧进行一次乘法运算?
这是为了将操作次数减少到最少的某种优化吗?
我可以在 GLM 中使用 v' = q * v
旋转来获得所有旋转的一致且正确的结果吗?
代码:
const int rotSteps = 4;
// Rotate around X axis in steps to 90deg
vec3 eulerState = vec3(0.0f);
// point we want to rotate (use vec4 to track the w component during rotations)
vec4 v = vec4(0.0f, 0.0f, 3.0f, 0.0f);
// Full Euler steps for q * v rotation
quat orientF = quat(1.0f, 0.0f, 0.0f, 0.0f);
vec3 euler = vec3(RAD(90.0f), RAD(0.0f), RAD(0.0f));
vec3 eulerStep = euler / (float)rotSteps;
quat qEulerF = quat(eulerStep); // GetRotQuat(eulerStep);
vec4 qa = ToAngularForm(qEulerF);
vec3 orientEuler = eulerAngles(qEulerF);
CLogD(TAG, "Rot Full Step Q [W, X, Y, Z]: " FMT_Q(4) " / " FMT_V3(2) "deg / " FMT_QA(2), PAR_Q(qEulerF), PAR_V3(degrees(orientEuler)), PAR_QA(qa));
// Half Euler steps for q * v * q^-1 rotation
quat orientH = quat(1.0f, 0.0f, 0.0f, 0.0f);
vec3 eulerStepH = eulerStep / 2.0f;
quat qEulerH = quat(eulerStepH); // GetRotQuat(eulerStepH);
qa = ToAngularForm(qEulerH);
orientEuler = eulerAngles(qEulerH);
CLogD(TAG, "Rot Half Step Q [W, X, Y, Z]: " FMT_Q(4) " / " FMT_V3(2) "deg / " FMT_QA(2), PAR_Q(qEulerH), PAR_V3(degrees(orientEuler)), PAR_QA(qa));
quat qEulerHI = inverse(qEulerH);
vec4 qai = ToAngularForm(qEulerHI);
orientEuler = eulerAngles(qEulerHI);
CLogD(TAG, "Rot Half Step Q^-1 [W, X, Y, Z]: " FMT_Q(4) " / " FMT_V3(2) "deg / " FMT_QA(2), PAR_Q(qEulerHI), PAR_V3(degrees(orientEuler)), PAR_QA(qai));
for (int rotStep = 1; rotStep <= rotSteps; ++rotStep)
{
// Track the absolute Euler rotation
eulerState += eulerStep;
// Rotate by incremental rotation as defined by Euler angles
orientH = qEulerH * orientH;
orientEuler = eulerAngles(orientH);
CLogI(TAG, "Rot Step %d. Curr Abs Q: " FMT_Q(4) "/" FMT_V3(2) "deg, Abs Euler: " FMT_V3(2) "deg",
rotStep, PAR_Q(orientH), PAR_V3(degrees(orientEuler)), PAR_V3(degrees(eulerState)));
// Transform the point using the correct q * v * q^-1 rotation and multiply from Left and Right
quat orientHI = inverse(orientH);
qa = ToAngularForm(orientH);
qai = ToAngularForm(orientHI);
vec4 rotV = orientH * v * orientHI;
CLogD(TAG, "Rot QL: " FMT_Q(4) " / " FMT_QA(1), PAR_Q(orientH), PAR_QA(qa));
CLogD(TAG, "Rot QR: " FMT_Q(4) " / " FMT_QA(1), PAR_Q(orientHI), PAR_QA(qai));
CLogD(TAG, "Rot LR -> " FMT_V4(1), PAR_V4(rotV));
// Transform the point using the incorrect q * v rotation and multiply from Left only
orientF = qEulerF * orientF;
qa = ToAngularForm(orientF);
rotV = orientF * v;
CLogD(TAG, "Rot QR: " FMT_Q(4) " / " FMT_QA(1), PAR_Q(orientF), PAR_QA(qa));
CLogD(TAG, "Rot L -> " FMT_V4(1), PAR_V4(rotV));
}
// Rotate for 90 degrees around the Z axis
// Full Euler steps for q * v rotation
euler = vec3(RAD(0.0f), RAD(0.0f), RAD(90.0f));
eulerStep = euler / (float)rotSteps;
qEulerF = quat(eulerStep); // GetRotQuat(eulerStep);
qa = ToAngularForm(qEulerF);
orientEuler = eulerAngles(qEulerF);
CLogD(TAG, "Rot Full Step Q [W, X, Y, Z]: " FMT_Q(4) " / " FMT_V3(2) "deg / " FMT_QA(2), PAR_Q(qEulerF), PAR_V3(degrees(orientEuler)), PAR_QA(qa));
// Half Euler steps for q * v * q^-1 rotation
eulerStepH = eulerStep / 2.0f;
qEulerH = quat(eulerStepH); // GetRotQuat(eulerStepH);
qa = ToAngularForm(qEulerH);
orientEuler = eulerAngles(qEulerH);
CLogD(TAG, "Rot Half Step Q [W, X, Y, Z]: " FMT_Q(4) " / " FMT_V3(2) "deg / " FMT_QA(2), PAR_Q(qEulerH), PAR_V3(degrees(orientEuler)), PAR_QA(qa));
qEulerHI = inverse(qEulerH);
qai = ToAngularForm(qEulerHI);
orientEuler = eulerAngles(qEulerHI);
CLogD(TAG, "Rot Half Step Q^-1 [W, X, Y, Z]: " FMT_Q(4) " / " FMT_V3(2) "deg / " FMT_QA(2), PAR_Q(qEulerHI), PAR_V3(degrees(orientEuler)), PAR_QA(qai));
for (int rotStep = 1; rotStep <= rotSteps; ++rotStep)
{
// Track the absolute Euler rotation
eulerState += eulerStep;
// Rotate by incremental rotation as defined by Euler angles
orientH = qEulerH * orientH;
orientEuler = eulerAngles(orientH);
CLogI(TAG, "Rot Step %d. Curr Abs Q: " FMT_Q(4) "/" FMT_V3(2) "deg, Abs Euler: " FMT_V3(2) "deg",
rotStep, PAR_Q(orientH), PAR_V3(degrees(orientEuler)), PAR_V3(degrees(eulerState)));
// Transform the point using the correct q * v * q^-1 rotation and multiply from Left and Right
quat orientHI = inverse(orientH);
qa = ToAngularForm(orientH);
qai = ToAngularForm(orientHI);
vec4 rotV = orientH * v * orientHI;
CLogD(TAG, "Rot QL: " FMT_Q(4) " / " FMT_QA(1), PAR_Q(orientH), PAR_QA(qa));
CLogD(TAG, "Rot QR: " FMT_Q(4) " / " FMT_QA(1), PAR_Q(orientHI), PAR_QA(qai));
CLogD(TAG, "Rot LR -> " FMT_V4(1), PAR_V4(rotV));
// Transform the point using the incorrect q * v rotation and multiply from Left only
orientF = qEulerF * orientF;
qa = ToAngularForm(orientF);
rotV = orientF * v;
CLogD(TAG, "Rot QR: " FMT_Q(4) " / " FMT_QA(1), PAR_Q(orientF), PAR_QA(qa));
CLogD(TAG, "Rot L -> " FMT_V4(1), PAR_V4(rotV));
}
输出:
Rot Full Step Q [W, X, Y, Z]: [ 0.9808, [ 0.1951, 0.0000, 0.0000]] / [ 22.50, -0.00, 0.00]deg / cos( 11.25) + sin( 11.25)( 1.00i + 0.00j + 0.00k)
Rot Half Step Q [W, X, Y, Z]: [ 0.9952, [ 0.0980, 0.0000, 0.0000]] / [ 11.25, -0.00, 0.00]deg / cos( 5.63) + sin( 5.63)( 1.00i + 0.00j + 0.00k)
Rot Half Step Q^-1 [W, X, Y, Z]: [ 0.9952, [-0.0980, -0.0000, -0.0000]] / [-11.25, -0.00, 0.00]deg / cos( 5.63) + sin( 5.63)(-1.00i + -0.00j + -0.00k)
Rot Step 1. Curr Abs Q: [ 0.9952, [ 0.0980, 0.0000, 0.0000]]/[ 11.25, -0.00, 0.00]deg, Abs Euler: [ 22.50, 0.00, 0.00]deg
Rot QL: [ 0.9952, [ 0.0980, 0.0000, 0.0000]] / cos( 5.6) + sin( 5.6)( 1.0i + 0.0j + 0.0k)
Rot QR: [ 0.9952, [-0.0980, -0.0000, -0.0000]] / cos( 5.6) + sin( 5.6)(-1.0i + -0.0j + -0.0k)
Rot LR -> [ 0.0, -1.1, 2.8, 0.0]
Rot QR: [ 0.9808, [ 0.1951, 0.0000, 0.0000]] / cos( 11.3) + sin( 11.3)( 1.0i + 0.0j + 0.0k)
Rot L -> [ 0.0, -1.1, 2.8, 0.0]
Rot Step 2. Curr Abs Q: [ 0.9808, [ 0.1951, 0.0000, 0.0000]]/[ 22.50, -0.00, 0.00]deg, Abs Euler: [ 45.00, 0.00, 0.00]deg
Rot QL: [ 0.9808, [ 0.1951, 0.0000, 0.0000]] / cos( 11.3) + sin( 11.3)( 1.0i + 0.0j + 0.0k)
Rot QR: [ 0.9808, [-0.1951, -0.0000, -0.0000]] / cos( 11.2) + sin( 11.2)(-1.0i + -0.0j + -0.0k)
Rot LR -> [ 0.0, -2.1, 2.1, 0.0]
Rot QR: [ 0.9239, [ 0.3827, 0.0000, 0.0000]] / cos( 22.5) + sin( 22.5)( 1.0i + 0.0j + 0.0k)
Rot L -> [ 0.0, -2.1, 2.1, 0.0]
Rot Step 3. Curr Abs Q: [ 0.9569, [ 0.2903, 0.0000, 0.0000]]/[ 33.75, -0.00, 0.00]deg, Abs Euler: [ 67.50, 0.00, 0.00]deg
Rot QL: [ 0.9569, [ 0.2903, 0.0000, 0.0000]] / cos( 16.9) + sin( 16.9)( 1.0i + 0.0j + 0.0k)
Rot QR: [ 0.9569, [-0.2903, -0.0000, -0.0000]] / cos( 16.9) + sin( 16.9)(-1.0i + -0.0j + -0.0k)
Rot LR -> [ 0.0, -2.8, 1.1, 0.0]
Rot QR: [ 0.8315, [ 0.5556, 0.0000, 0.0000]] / cos( 33.8) + sin( 33.8)( 1.0i + 0.0j + 0.0k)
Rot L -> [ 0.0, -2.8, 1.1, 0.0]
Rot Step 4. Curr Abs Q: [ 0.9239, [ 0.3827, 0.0000, 0.0000]]/[ 45.00, -0.00, 0.00]deg, Abs Euler: [ 90.00, 0.00, 0.00]deg
Rot QL: [ 0.9239, [ 0.3827, 0.0000, 0.0000]] / cos( 22.5) + sin( 22.5)( 1.0i + 0.0j + 0.0k)
Rot QR: [ 0.9239, [-0.3827, -0.0000, -0.0000]] / cos( 22.5) + sin( 22.5)(-1.0i + -0.0j + -0.0k)
Rot LR -> [ 0.0, -3.0, 0.0, 0.0]
Rot QR: [ 0.7071, [ 0.7071, 0.0000, 0.0000]] / cos( 45.0) + sin( 45.0)( 1.0i + 0.0j + 0.0k)
Rot L -> [ 0.0, -3.0, 0.0, 0.0]
Rot Full Step Q [W, X, Y, Z]: [ 0.9808, [ 0.0000, 0.0000, 0.1951]] / [ 0.00, -0.00, 22.50]deg / cos( 11.25) + sin( 11.25)( 0.00i + 0.00j + 1.00k)
Rot Half Step Q [W, X, Y, Z]: [ 0.9952, [ 0.0000, 0.0000, 0.0980]] / [ 0.00, -0.00, 11.25]deg / cos( 5.63) + sin( 5.63)( 0.00i + 0.00j + 1.00k)
Rot Half Step Q^-1 [W, X, Y, Z]: [ 0.9952, [-0.0000, -0.0000, -0.0980]] / [ 0.00, -0.00, -11.25]deg / cos( 5.63) + sin( 5.63)(-0.00i + -0.00j + -1.00k)
Rot Step 1. Curr Abs Q: [ 0.9194, [ 0.3808, 0.0375, 0.0906]]/[ 45.00, 0.00, 11.25]deg, Abs Euler: [ 90.00, 0.00, 22.50]deg
Rot QL: [ 0.9194, [ 0.3808, 0.0375, 0.0906]] / cos( 23.2) + sin( 23.2)( 1.0i + 0.1j + 0.2k)
Rot QR: [ 0.9194, [-0.3808, -0.0375, -0.0906]] / cos( 23.2) + sin( 23.2)(-1.0i + -0.1j + -0.2k)
Rot LR -> [ 1.0, -2.8, 0.0, 0.0]
Rot QR: [ 0.6935, [ 0.6935, 0.1379, 0.1379]] / cos( 46.1) + sin( 46.1)( 1.0i + 0.2j + 0.2k)
Rot L -> [ 1.1, -2.8, 0.0, 0.0]
Rot Step 2. Curr Abs Q: [ 0.9061, [ 0.3753, 0.0747, 0.1802]]/[ 45.00, -0.00, 22.50]deg, Abs Euler: [ 90.00, 0.00, 45.00]deg
Rot QL: [ 0.9061, [ 0.3753, 0.0747, 0.1802]] / cos( 25.0) + sin( 25.0)( 0.9i + 0.2j + 0.4k)
Rot QR: [ 0.9061, [-0.3753, -0.0747, -0.1802]] / cos( 25.0) + sin( 25.0)(-0.9i + -0.2j + -0.4k)
Rot LR -> [ 1.9, -2.4, 0.1, 0.0]
Rot QR: [ 0.6533, [ 0.6533, 0.2706, 0.2706]] / cos( 49.2) + sin( 49.2)( 0.9i + 0.4j + 0.4k)
Rot L -> [ 2.1, -2.1, 0.0, 0.0]
Rot Step 3. Curr Abs Q: [ 0.8841, [ 0.3662, 0.1111, 0.2682]]/[ 45.00, 0.00, 33.75]deg, Abs Euler: [ 90.00, 0.00, 67.50]deg
Rot QL: [ 0.8841, [ 0.3662, 0.1111, 0.2682]] / cos( 27.9) + sin( 27.9)( 0.8i + 0.2j + 0.6k)
Rot QR: [ 0.8841, [-0.3662, -0.1111, -0.2682]] / cos( 27.9) + sin( 27.9)(-0.8i + -0.2j + -0.6k)
Rot LR -> [ 2.5, -1.6, 0.3, 0.0]
Rot QR: [ 0.5879, [ 0.5879, 0.3928, 0.3928]] / cos( 54.0) + sin( 54.0)( 0.7i + 0.5j + 0.5k)
Rot L -> [ 2.8, -1.1, 0.0, 0.0]
Rot Step 4. Curr Abs Q: [ 0.8536, [ 0.3536, 0.1464, 0.3536]]/[ 45.00, 0.00, 45.00]deg, Abs Euler: [ 90.00, 0.00, 90.00]deg
Rot QL: [ 0.8536, [ 0.3536, 0.1464, 0.3536]] / cos( 31.4) + sin( 31.4)( 0.7i + 0.3j + 0.7k)
Rot QR: [ 0.8536, [-0.3536, -0.1464, -0.3536]] / cos( 31.4) + sin( 31.4)(-0.7i + -0.3j + -0.7k)
Rot LR -> [ 2.9, -0.7, 0.4, 0.0]
Rot QR: [ 0.5000, [ 0.5000, 0.5000, 0.5000]] / cos( 60.0) + sin( 60.0)( 0.6i + 0.6j + 0.6k)
Rot L -> [ 3.0, 0.0, 0.0, 0.0]
我得到了我的问题的答案和一台工作正常的轨道相机,但没有时间仔细检查示例代码现在是否正常工作 - 它应该。
第一个问题是为什么 GLM 在四元数转换过程中将角度减半,根据 extended Euler's formula 看起来是这样......它必须这样做。这部分可以再研究一下,但由于时间不够,我不得不接受它。
GLM 中的矢量旋转是使用乘法运算符实现的。这意味着当将 vec3 与四元数相乘时,不会将向量转换为四元数然后执行乘法,它会执行 vector rotation instead:
template <typename T, precision P>
GLM_FUNC_QUALIFIER tvec3<T, P> operator*(tquat<T, P> const & q, tvec3<T, P> const & v)
{
tvec3<T, P> const QuatVector(q.x, q.y, q.z);
tvec3<T, P> const uv(glm::cross(QuatVector, v));
tvec3<T, P> const uuv(glm::cross(QuatVector, uv));
return v + ((uv * q.w) + uuv) * static_cast<T>(2);
}
所以,是的,使用四元数旋转向量的正确方法是在四元数和向量之间使用乘法运算符,如下所示:
v' = q * v
或者在 C++ 中:
vec3 posOrigin;
quat rotQ;
...
vec3 posRot = rotQ * posOrigin;
这段代码实际上并没有直接进行四元数乘法。它进行旋转。就个人而言,我更喜欢 GLM 提供像 rotate(quat, vec)
这样的函数调用...但我确信运算符重载混淆是有原因的。
另请注意,操作数顺序很重要,因为向量和四元数之间的乘法定义如下:
template <typename T, precision P>
GLM_FUNC_QUALIFIER tvec3<T, P> operator*(tvec3<T, P> const & v, tquat<T, P> const & q)
{
return glm::inverse(q) * v;
}
因此将反向旋转矢量。
请注意,GLM 还实现了四元数之间的乘法,但为此需要使用两个四元数之间的乘法运算符:
template <typename T, precision P>
template <typename U>
GLM_FUNC_QUALIFIER tquat<T, P> & tquat<T, P>::operator*=(tquat<U, P> const & r)
{
tquat<T, P> const p(*this);
tquat<T, P> const q(r);
this->w = p.w * q.w - p.x * q.x - p.y * q.y - p.z * q.z;
this->x = p.w * q.x + p.x * q.w + p.y * q.z - p.z * q.y;
this->y = p.w * q.y + p.y * q.w + p.z * q.x - p.x * q.z;
this->z = p.w * q.z + p.z * q.w + p.x * q.y - p.y * q.x;
return *this;
}
由于 GLM 几乎没有我能找到的宝贵文档,因此这种运算符重载会导致错误的假设和大量的时间损失。所以我想我应该阅读 GLM 代码而不是假设它的作用...