Python – 3Dvector的旋转
我有两个向量作为Python列表和一个angular度。 例如：
v = [3,5,0] axis = [4,4,1] theta = 1.2 #radian 绕v轴旋转vvector时，得到最终vector的最佳/最简单的方法是什么？
对于轴向量指向的观察者，旋转应该是逆时针的。 这就是所谓的右手规则
看看http://vpython.org/contents/docs/visual/VisualIntro.html 。
它提供了一个vector类，它有一个方法A.rotate(theta,B) 。 如果你不想调用A上的方法，它还提供了一个帮助函数rotate(A,theta,B) 。
使用欧拉 – 罗德里格公式 ：
import numpy as np import math def rotation_matrix(axis, theta): """ Return the rotation matrix associated with counterclockwise rotation about the given axis by theta radians. """ axis = np.asarray(axis) axis = axis/math.sqrt(np.dot(axis, axis)) a = math.cos(theta/2.0) b, c, d = -axis*math.sin(theta/2.0) aa, bb, cc, dd = a*a, b*b, c*c, d*d bc, ad, ac, ab, bd, cd = b*c, a*d, a*c, a*b, b*d, c*d return np.array([[aa+bb-cc-dd, 2*(bc+ad), 2*(bd-ac)], [2*(bc-ad), aa+cc-bb-dd, 2*(cd+ab)], [2*(bd+ac), 2*(cd-ab), aa+dd-bb-cc]]) v = [3, 5, 0] axis = [4, 4, 1] theta = 1.2 print(np.dot(rotation_matrix(axis,theta), v)) # [ 2.74911638 4.77180932 1.91629719]
单线程，具有numpy / scipy函数。
我们使用以下内容：
令a为沿轴的单位vector，即a =轴/范数（轴）
A = I×a是与a相关的斜对称matrix，即单位matrix与a的叉积那么M = exp（θA）是旋转matrix。
from numpy import cross, eye, dot from scipy.linalg import expm3, norm def M(axis, theta): return expm3(cross(eye(3), axis/norm(axis)*theta)) v, axis, theta = [3,5,0], [4,4,1], 1.2 M0 = M(axis, theta) print(dot(M0,v)) # [ 2.74911638 4.77180932 1.91629719]
expm3 （这里的代码）计算指数的泰勒级数：
\sum_{k=0}^{20} \frac{1}{k!} (θ A)^k ，所以这是时间昂贵的，但可读性和安全性。 如果你有很less的轮换，但是有很多向量，这可能是一个好的方法。
我只想提到，如果需要速度，将unutbu的代码包装到scipy的weave.inline中，并将已经存在的matrix作为parameter passing，会使运行时间减less20倍。
代码（在rotation_matrix_test.py中）：
import numpy as np import timeit from math import cos, sin, sqrt import numpy.random as nr from scipy import weave def rotation_matrix_weave(axis, theta, mat = None): if mat == None: mat = np.eye(3,3) support = "#include <math.h>" code = """ double x = sqrt(axis[0] * axis[0] + axis[1] * axis[1] + axis[2] * axis[2]); double a = cos(theta / 2.0); double b = -(axis[0] / x) * sin(theta / 2.0); double c = -(axis[1] / x) * sin(theta / 2.0); double d = -(axis[2] / x) * sin(theta / 2.0); mat[0] = a*a + b*b - c*c - d*d; mat[1] = 2 * (b*c - a*d); mat[2] = 2 * (b*d + a*c); mat[3*1 + 0] = 2*(b*c+a*d); mat[3*1 + 1] = a*a+c*cb*bd*d; mat[3*1 + 2] = 2*(c*da*b); mat[3*2 + 0] = 2*(b*da*c); mat[3*2 + 1] = 2*(c*d+a*b); mat[3*2 + 2] = a*a+d*db*bc*c; """ weave.inline(code, ['axis', 'theta', 'mat'], support_code = support, libraries = ['m']) return mat def rotation_matrix_numpy(axis, theta): mat = np.eye(3,3) axis = axis/sqrt(np.dot(axis, axis)) a = cos(theta/2.) b, c, d = -axis*sin(theta/2.) return np.array([[a*a+b*bc*cd*d, 2*(b*ca*d), 2*(b*d+a*c)], [2*(b*c+a*d), a*a+c*cb*bd*d, 2*(c*da*b)], [2*(b*da*c), 2*(c*d+a*b), a*a+d*db*bc*c]])
时机：
>>> import timeit >>> >>> setup = """ ... import numpy as np ... import numpy.random as nr ... ... from rotation_matrix_test import rotation_matrix_weave ... from rotation_matrix_test import rotation_matrix_numpy ... ... mat1 = np.eye(3,3) ... theta = nr.random() ... axis = nr.random(3) ... """ >>> >>> timeit.repeat("rotation_matrix_weave(axis, theta, mat1)", setup=setup, number=100000) [0.36641597747802734, 0.34883809089660645, 0.3459300994873047] >>> timeit.repeat("rotation_matrix_numpy(axis, theta)", setup=setup, number=100000) [7.180983066558838, 7.172032117843628, 7.180462837219238]
我为Python {2,3}做了一个相当完整的三维math库。 它仍然不使用Cython，但很大程度上依赖于numpy的效率。 你可以在这里find它的点子：
python[3] -m pip install math3d
或者看看我的gitweb http://git.automatics.dyndns.dk/?p=pymath3d.git ，现在也在github上： https ： //github.com/mortlind/pymath3d 。
一旦安装，在python中，您可以创build可以旋转vector的方向对象，或者成为变换对象的一部分。 例如，下面的代码片段组成一个方向，表示围绕轴[1，3]旋转1 rad，将其应用于vector[4,5,6]，并打印结果：
import math3d as m3d r = m3d.Orientation.new_axis_angle([1,2,3], 1) v = m3d.Vector(4,5,6) print(r * v)
输出将是
<Vector: (2.53727, 6.15234, 5.71935)>
就我所能做到的而言，这个效率比使用上面提到的scipy的oneliner效率要高四倍左右。 但是，它需要安装我的math3d软件包。
这是一个使用快速四元数的优雅方法; 我可以用适当的向量化numpy数组计算每秒1000万转。 它依赖于在这里find的numpy的四元数的扩展。
四元数理论：四元数是一个具有一个实数和三个虚数维的数，通常写成q = w + xi + yj + zk ，其中'i'，'j'，'k'是虚数维。 正如一个单位复数“c”可以用c=exp(i * theta)表示所有二次旋转，单位四元数“q”可以用q=exp(p)表示所有的3d旋转，其中'p'是一个纯粹的由您的轴和angular度设置的虚四元数。
我们首先把你的轴和angular度转换成一个四元数，它的虚轴由你的旋转轴线给出，其大小由弧度旋转angular度的一半给出。 4个元素向量(w, x, y, z)构造如下：
import numpy as np import quaternion as quat v = [3,5,0] axis = [4,4,1] theta = 1.2 #radian vector = np.array([0.] + v) rot_axis = np.array([0.] + axis) axis_angle = (theta*0.5) * rot_axis/np.linalg.norm(rot_axis)
首先，对于要被旋转的vector和旋转轴rot_axis用实分量w = 0来构造4个元素的一个numpyarrays。 轴angular表示然后通过归一化然后乘以所需angular度theta一半来构build。 看到这里为什么一半的angular度是必需的。
现在使用库创build四元数v和qlog ，并通过取指数得到单位旋转四元数q 。
vec = quat.quaternion(*v) qlog = quat.quaternion(*axis_angle) q = np.exp(qlog)
最后，vector的旋转通过以下操作来计算。
v_prime = q * vec * np.conjugate(q) print(v_prime) # quaternion(0.0, 2.7491163, 4.7718093, 1.9162971)
现在放弃真正的元素，你有你的旋转vector！ 请注意，如果必须通过许多顺序旋转来旋转vector，则此方法特别有效，因为四元数乘积可以仅计算为q = q1 * q2 * q3 * q4 * … * qn，然后仅vector旋转由'q'在最后使用v'= q * v * conj（q）。
这种方法可以简单地通过exp和log函数（yes log(q)只是返回轴angular表示法）给出一个无缝的轴angular转换运算符。 为了进一步说明四元数乘法等等是如何工作的，请看这里
使用pyquaternion非常简单， 要安装它，请在您的控制台中运行：import pip; pip.main（[ '安装'， 'pyquaternion']）
一旦安装：
from pyquaternion import Quaternion v = [3,5,0] axis = [4,4,1] theta = 1.2 #radian rotated_v = Quaternion(axis=axis,angle=theta).rotate(v)
