logo标签列表

多个四元数的“平均值”是什么？

mathquaternions
70

我正在尝试在我的OpenGL程序中将骨架动画的矩阵转换为四元数，但是我遇到了一个问题：

给定一些单位四元数，我需要获得一个四元数，用它来变换一个向量，使得变换后的向量是每个四元数分别变换后向量的平均值。(对于矩阵，我只需要将矩阵相加并除以矩阵数量即可)

5
这不是一个可解决的问题。四元数无法编码三维空间中任意线性变换，它们只能编码正交变换。但是多个正交变换的“平均值”（即：对于每个向量，取其所有变换并求平均）通常不是正交的，因此不能通过四元数获得。可能存在一些更微妙的“平均”概念可以保持正交性，但你得不到你所想要的平均值。 - darij grinberg
1
所有关于四元数的可能算术运算都在这个页面中（http://www.euclideanspace.com/maths/algebra/realNormedAlgebra/quaternions/index.htm）。如果你在那里找不到，那么...祝你好运。 - Jav_Rock
你有没有考虑过查看开源骨骼动画代码，而不是重复造轮子？例如：http://home.gna.org/cal3d/ - JWWalker
相关：https://math.stackexchange.com/questions/61146/averaging-quaternions - Felipe G. Nievinski
16个回答
1

使用Python>=3.6中的Numpy最简单的Markley解决方案实现如下：

import numpy as np
def q_average(Q, W=None):
    if W is not None:
        Q *= W[:, None]
    eigvals, eigvecs = np.linalg.eig(Q.T@Q)
    return eigvecs[:, eigvals.argmax()]

其中Q的大小为N x 4。得到的四元数已经被归一化。

在这种情况下，默认情况下权重相等为1。否则，您可以提供一个权重列表，大小为N（每个四元数一个权重）。

就是这样。

"limit" 版本是由 PEP465 引入的 @ 操作符引起的。https://www.python.org/dev/peps/pep-0465。@ 操作符表示矩阵乘法。 - Zailux
注意，您加权的是 Q 而不是 Q*Q'。因此，您的权重是算法不同版本的平方。 - minorlogic
1
使用四元数可以做同样的事情，但需要进行一些小的修正： 1. 如果该四元数与之前求和结果的点积为负，则在求平均值之前对其取反。 2. 如果您的库使用单位四元数，则在求平均值结束时对平均四元数进行归一化。
平均四元数将近似表示平均旋转（最大误差约为5度）。
警告：如果旋转差异太大，则不同方向的平均矩阵可能会出现错误。
1

这是我的PyTorch实现，用于计算平均四元数。

def mean(Q, weights=None):
    if weights is None:
        weights = torch.ones(len(Q), device=torch.device("cuda:0")) / len(Q)
    A = torch.zeros((4, 4), device=torch.device("cuda:0"))
    weight_sum = torch.sum(weights)

    oriented_Q = ((Q[:, 0:1] > 0).float() - 0.5) * 2 * Q
    A = torch.einsum("bi,bk->bik", (oriented_Q, oriented_Q))
    A = torch.sum(torch.einsum("bij,b->bij", (A, weights)), 0)
    A /= weight_sum

    q_avg = torch.linalg.eigh(A)[1][:, -1]
    if q_avg[0] < 0:
        return -q_avg
    return q_avg

我使用了Markley, F. Landis等人的“平均四元数”算法，该算法在http://tbirdal.blogspot.com/2019/10/i-allocate-this-post-to-providing.html中介绍。

0
我想对这些解决方案发表评论，但由于我没有50个goo goo星星或其他需要允许评论的东西，我想指出在Markley方法中寻找最大特征值时，我发现幂迭代（https://en.wikipedia.org/wiki/Power_iteration）比使用小旋转的奇异值分解更稳健。至少在Eigen的JacobiSVD中，它经常崩溃。
这是一个回答还是一个评论？在只能放回答的地方。 - starball
这并没有回答问题。一旦你拥有足够的声望，你就可以评论任何帖子；相反，提供不需要提问者澄清的答案。- 来自审查 - Yaroslavm
0
-1

有另一种基于slerp的方法，当平均四元数时，这可能是您实际想要的。

让我们首先将其与基于特征分析的平均值进行比较：

考虑两个四元数A和B的平均值，它们对应于绕X轴旋转0°和90°度，权重为w_A = 2和w_B = 1。 预期的加权平均值应该对应于绕X轴旋转30°。 基于slerp的加权平均值将返回预期值。 基于特征分析的加权平均值将返回一个26.56°的旋转。 enter image description here

基于特征分析的方法将返回最小化相应旋转矩阵的Frobenius范数的四元数。相反，基于slerp的方法将在四元数的3D旋转空间中计算平均值。

import math
import numpy as np
import quaternion # (pip install numpy-quaternion)
d2r = math.pi/180
r2d = 180/math.pi

def recover_angle(mat):
    return np.arctan2(mat[1,0], mat[0,0])

# ground truth
angles = np.array([0,90])
weights = np.array([2,1])
mean_angle = np.sum(angles*(weights/weights.sum()))
quaternion_mean = quaternion.from_euler_angles(mean_angle*d2r,0,0)

# eigen analysis
Q = quaternion.as_float_array(
    [
        (weight**0.5) * quaternion.from_euler_angles(0,0,angle*d2r) 
        for angle,weight 
        in zip(angles,weights)
    ]
).T
QQT = Q @ Q.T
eigval,eigvec = np.linalg.eig(QQT)
quaternion_mean_eig = quaternion.from_float_array( eigvec[:,np.argmax(eigval)] )

# slerp based
def slerp_avg(quaternions, weights):
    # welford's mean in terms of linear mix operations:
    toqt = quaternion.from_float_array
    mix = lambda a,b,k: quaternion.slerp(toqt(a),toqt(b),0,1,k)
    wmean, wsum, num = quaternions[0], weights[0], len(weights)
    for i in range(1,num):
        wsum += weights[i]
        k = (weights[i]/wsum)
        # wmean := wmean*(1-k) + quaternions[i]*k
        wmean = mix(wmean,quaternions[i],k) 
    return wmean

quaternion_mean_slerp = slerp_avg(quaternion.as_float_array(
    [quaternion.from_euler_angles(0,0,angle*d2r) for angle in angles]
), weights)

mat_mean = quaternion.as_rotation_matrix(quaternion_mean)
mat_mean_eig = quaternion.as_rotation_matrix(quaternion_mean_eig)
mat_mean_slerp = quaternion.as_rotation_matrix(quaternion_mean_slerp)

print("expected", recover_angle(mat_mean)*r2d)
print("eigen", recover_angle(mat_mean_eig)*r2d)
print("slerp", recover_angle(mat_mean_slerp)*r2d)

输出：

expected 29.999999999999996
eigen 26.565051177077994
slerp 30.00000000000001

您可以在此处找到一个零依赖的C++库： https://github.com/xaedes/average_affine_transform_mat/

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文，
原文链接