我有两个向量 u 和 v。有没有办法找到代表从 u 到 v 的旋转的四元数?
9 回答
Quaternion q;
vector a = crossproduct(v1, v2);
q.xyz = a;
q.w = sqrt((v1.Length ^ 2) * (v2.Length ^ 2)) + dotproduct(v1, v2);
不要忘记标准化 q。
理查德关于没有独特的旋转是正确的,但上面应该给出“最短的弧线”,这可能是你需要的。
中途矢量解决方案
我想出了我认为 Imbrondir 试图提出的解决方案(尽管有一个小错误,这可能是 sinisterchipmunk 无法验证它的原因)。
假设我们可以构造一个表示绕轴旋转的四元数,如下所示:
q.w == cos(angle / 2)
q.x == sin(angle / 2) * axis.x
q.y == sin(angle / 2) * axis.y
q.z == sin(angle / 2) * axis.z
并且两个归一化向量的点积和叉积是:
dot == cos(theta)
cross.x == sin(theta) * perpendicular.x
cross.y == sin(theta) * perpendicular.y
cross.z == sin(theta) * perpendicular.z
看到从u到v的旋转可以通过围绕垂直向量旋转 theta(向量之间的角度)来实现,看起来我们可以直接从点积和叉积的结果构造一个表示这种旋转的四元数; 但是,就目前而言,theta = angle / 2,这意味着这样做会导致所需旋转的两倍。
一种解决方案是计算u和v之间的中间向量,并使用u和中间向量的点积和叉积来构造一个四元数,表示u和中间向量之间角度两倍的旋转,这将我们一路带到v!
有一种特殊情况,其中u == -v并且无法计算唯一的中途向量。这是意料之中的,因为可以将我们从u带到v的无限多“最短弧”旋转,并且我们必须简单地围绕与u(或v )正交的任何向量旋转 180 度作为我们的特例解决方案。这是通过将 u 与不平行于 u 的任何其他向量进行归一化叉积来完成的。
伪代码如下(显然,在现实中,特殊情况必须考虑浮点不准确性——可能通过根据某个阈值而不是绝对值检查点积)。
另请注意,当u == v时没有特殊情况(产生身份四元数 - 请自行检查并查看)。
// N.B. the arguments are _not_ axis and angle, but rather the
// raw scalar-vector components.
Quaternion(float w, Vector3 xyz);
Quaternion get_rotation_between(Vector3 u, Vector3 v)
{
// It is important that the inputs are of equal length when
// calculating the half-way vector.
u = normalized(u);
v = normalized(v);
// Unfortunately, we have to check for when u == -v, as u + v
// in this case will be (0, 0, 0), which cannot be normalized.
if (u == -v)
{
// 180 degree rotation around any orthogonal vector
return Quaternion(0, normalized(orthogonal(u)));
}
Vector3 half = normalized(u + v);
return Quaternion(dot(u, half), cross(u, half));
}
该orthogonal
函数返回与给定向量正交的任何向量。此实现使用具有最正交基向量的叉积。
Vector3 orthogonal(Vector3 v)
{
float x = abs(v.x);
float y = abs(v.y);
float z = abs(v.z);
Vector3 other = x < y ? (x < z ? X_AXIS : Z_AXIS) : (y < z ? Y_AXIS : Z_AXIS);
return cross(v, other);
}
中途四元数解
这实际上是接受的答案中提出的解决方案,它似乎比中途矢量解决方案略快(根据我的测量,快了约 20%,尽管不要相信我的话)。我在这里添加它以防像我这样的其他人对解释感兴趣。
本质上,您可以计算导致所需旋转两倍的四元数,而不是使用中途向量计算四元数(如其他解决方案中所述),并在该和零度之间找到四元数。
正如我之前解释的,双倍旋转所需的四元数是:
q.w == dot(u, v)
q.xyz == cross(u, v)
零旋转的四元数是:
q.w == 1
q.xyz == (0, 0, 0)
计算中途四元数只是简单地对四元数求和并归一化结果,就像向量一样。但是,与向量的情况一样,四元数必须具有相同的大小,否则结果将偏向具有较大大小的四元数。
由两个向量的点积和叉积构成的四元数将具有与这些积相同的大小:length(u) * length(v)
。与其将所有四个分量都除以这个因子,不如扩大单位四元数。如果您想知道为什么接受的答案似乎使用 会使事情变得复杂sqrt(length(u) ^ 2 * length(v) ^ 2)
,那是因为向量的平方长度比长度计算得更快,所以我们可以省去一次sqrt
计算。结果是:
q.w = dot(u, v) + sqrt(length_2(u) * length_2(v))
q.xyz = cross(u, v)
然后对结果进行归一化。伪代码如下:
Quaternion get_rotation_between(Vector3 u, Vector3 v)
{
float k_cos_theta = dot(u, v);
float k = sqrt(length_2(u) * length_2(v));
if (k_cos_theta / k == -1)
{
// 180 degree rotation around any orthogonal vector
return Quaternion(0, normalized(orthogonal(u)));
}
return normalized(Quaternion(k_cos_theta + k, cross(u, v)));
}
如上所述的问题没有明确定义:给定的向量对没有唯一的旋转。例如,考虑 u = <1, 0, 0>和 v = <0, 1, 0>的情况。从 u 到 v 的一次旋转将是围绕 z 轴的pi / 2旋转。从 u 到 v 的另一个旋转将是围绕向量<1, 1, 0>的pi旋转。
为什么不使用纯四元数来表示向量?如果你先规范化它们会更好。
q 1 = (0 u x u y u z )'
q 2 = (0 v x v y v z )'
q 1 q rot = q 2
预乘以 q 1 -1
q rot = q 1 -1 q 2
其中 q 1 -1 = q 1 conj / q norm
这可以被认为是“左分割”。右除法,这不是你想要的:
q rot,right = q 2 -1 q 1
我不太擅长四元数。然而,我为此苦苦挣扎了几个小时,无法让 Polaris878 解决方案发挥作用。我尝试过对 v1 和 v2 进行预规范化。归一化 q。标准化 q.xyz。然而我还是不明白。结果仍然没有给我正确的结果。
最后,虽然我找到了一个解决方案。如果它对其他人有帮助,这是我的工作(python)代码:
def diffVectors(v1, v2):
""" Get rotation Quaternion between 2 vectors """
v1.normalize(), v2.normalize()
v = v1+v2
v.normalize()
angle = v.dot(v2)
axis = v.cross(v2)
return Quaternion( angle, *axis )
如果 v1 和 v2 像 v1 == v2 或 v1 == -v2 一样平行(有一定的公差),则必须提出一个特殊情况,我认为解决方案应该是 Quaternion(1, 0,0,0) (无旋转)或 Quaternion(0, *v1)(180 度旋转)
从算法的角度来看,最快的解决方案是伪代码。
Quaternion shortest_arc(const vector3& v1, const vector3& v2 )
{
// input vectors NOT unit
Quaternion q( cross(v1, v2), dot(v1, v2) );
// reducing to half angle
q.w += q.magnitude(); // 4 multiplication instead of 6 and more numerical stable
// handling close to 180 degree case
//... code skipped
return q.normalized(); // normalize if you need UNIT quaternion
}
确保您需要单位四元数(通常,它是插值所必需的)。
注意:非单位四元数可以用于比单位更快的一些操作。
一些答案似乎没有考虑叉积可能为0的可能性。下面的代码片段使用角轴表示:
//v1, v2 are assumed to be normalized
Vector3 axis = v1.cross(v2);
if (axis == Vector3::Zero())
axis = up();
else
axis = axis.normalized();
return toQuaternion(axis, ang);
toQuaternion
可以按如下方式实现:
static Quaternion toQuaternion(const Vector3& axis, float angle)
{
auto s = std::sin(angle / 2);
auto u = axis.normalized();
return Quaternion(std::cos(angle / 2), u.x() * s, u.y() * s, u.z() * s);
}
如果您使用的是 Eigen 库,您也可以这样做:
Quaternion::FromTwoVectors(from, to)
仅使用归一化四元数,我们可以用以下术语表达约瑟夫汤普森的答案。
令 q_v = (0, u_x, v_y, v_z) 和 q_w = (0, v_x, v_y, v_z) 并考虑
q = q_v * q_w = (-u 点 v, uxv)。
所以将 q 表示为 q(q_0, q_1, q_2, q_3) 我们有
q_r = (1 - q_0, q_1, q_2, q_3).normalize()
根据两个角度之间四元数旋转的推导,可以将向量u旋转到向量v
function fromVectors(u, v) {
d = dot(u, v)
w = cross(u, v)
return Quaternion(d + sqrt(d * d + dot(w, w)), w).normalize()
}
如果已知向量u到向量v是单位向量,则函数简化为
function fromUnitVectors(u, v) {
return Quaternion(1 + dot(u, v), cross(u, v)).normalize()
}
根据您的用例,可能需要处理点积为 1(平行向量)和 -1(指向相反方向的向量)的情况。