python vector rotation

python - Rotación del vector 3D?



rotation (8)

Aquí hay un método elegante que usa cuaterniones que son tremendamente rápidos; Puedo calcular 10 millones de rotaciones por segundo con matrices numpy adecuadamente vectorizadas. Se basa en la extensión quaternion para numpy que se encuentra here .

Teoría del cuaternión: un cuaternión es un número con una dimensión real y 3 dimensiones imaginarias usualmente escritas como q = w + xi + yj + zk donde ''i'', ''j'', ''k'' son dimensiones imaginarias. Así como un número complejo de unidades ''c'' puede representar todas las 2d rotaciones por c=exp(i * theta) , una unidad cuaternion ''q'' puede representar todas las rotaciones 3d por q=exp(p) , donde ''p'' es un puro quaternion imaginario establecido por su eje y ángulo.

Comenzamos convirtiendo su eje y ángulo en un cuaternión cuyas dimensiones imaginarias vienen dadas por su eje de rotación, y cuya magnitud está dada por la mitad del ángulo de rotación en radianes. Los 4 vectores de elementos (w, x, y, z) se construyen de la siguiente manera:

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)

En primer lugar, se construye una matriz numpy de 4 elementos con el componente real w = 0 para el vector rot_axis vector y el eje de rotación rot_axis . La representación del ángulo del eje se construye normalizando luego multiplicando por la mitad el ángulo deseado theta . Vea here por qué se requiere la mitad del ángulo.

Ahora crea los cuaterniones v y qlog usando la biblioteca, y obtén la unidad de rotación quaternion q tomando la exponencial.

vec = quat.quaternion(*v) qlog = quat.quaternion(*axis_angle) q = np.exp(qlog)

Finalmente, la rotación del vector se calcula mediante la siguiente operación.

v_prime = q * vec * np.conjugate(q) print(v_prime) # quaternion(0.0, 2.7491163, 4.7718093, 1.9162971)

¡Ahora solo descarta el elemento real y tienes tu vector rotado! Tenga en cuenta que este método es particularmente eficiente si tiene que rotar un vector a través de muchas rotaciones secuenciales, ya que el producto del cuaternión se puede calcular simplemente como q = q1 * q2 * q3 * q4 * ... * qn y luego el vector solo se rota por ''q'' al final usando v ''= q * v * conj (q).

Este método le proporciona una transformación perfecta entre el operador de rotación 3d del ángulo del eje <---> simplemente mediante las funciones exp y log (¡sí log(q) simplemente devuelve la representación del ángulo del eje!). Para mayor aclaración de cómo funciona la multiplicación de cuaterniones, etc., consulte here

Tengo dos vectores como listas de Python y un ángulo. P.ej:

v = [3,5,0] axis = [4,4,1] theta = 1.2 #radian

¿Cuál es la mejor / más fácil forma de obtener el vector resultante al rotar el vector v alrededor del eje?

La rotación debería parecer en sentido contrario a las agujas del reloj para un observador al que apunta el vector del eje. Esto se llama la regla de la mano derecha

Descargo de responsabilidad: soy el autor de este paquete

Si bien las clases especiales para rotaciones pueden ser convenientes, en algunos casos se necesitan matrices de rotación (por ejemplo, para trabajar con otras bibliotecas como las funciones affine_transform en scipy). Para evitar que todos implementen sus propias funciones de generación de matriz, existe un pequeño paquete de Python puro que no hace más que proporcionar funciones de generación de matriz de rotación convenientes. El paquete está en github ( mgen ) y se puede instalar a través de pip:

pip install mgen

Ejemplo de uso copiado del archivo léame:

import numpy as np np.set_printoptions(suppress=True) from mgen import rotation_around_axis from mgen import rotation_from_angles from mgen import rotation_around_x matrix = rotation_from_angles([np.pi/2, 0, 0], ''XYX'') matrix.dot([0, 1, 0]) # array([0., 0., 1.]) matrix = rotation_around_axis([1, 0, 0], np.pi/2) matrix.dot([0, 1, 0]) # array([0., 0., 1.]) matrix = rotation_around_x(np.pi/2) matrix.dot([0, 1, 0]) # array([0., 0., 1.])

Hice una biblioteca bastante completa de matemáticas 3D para Python {2,3}. Todavía no usa Cython, pero depende en gran medida de la eficiencia de numpy. Puedes encontrarlo aquí con pip:

python[3] -m pip install math3d

O eche un vistazo a mi gitweb http://git.automatics.dyndns.dk/?p=pymath3d.git y ahora también en github: https://github.com/mortlind/pymath3d .

Una vez instalado, en python puede crear el objeto de orientación que puede rotar vectores o ser parte de objetos de transformación. Por ejemplo, el siguiente fragmento de código compone una orientación que representa una rotación de 1 rad alrededor del eje [1,2,3], la aplica al vector [4,5,6] e imprime el resultado:

import math3d as m3d r = m3d.Orientation.new_axis_angle([1,2,3], 1) v = m3d.Vector(4,5,6) print(r * v)

El resultado sería

<Vector: (2.53727, 6.15234, 5.71935)>

Esto es más eficiente, por un factor de aproximadamente cuatro, por lo que puedo cronometrar, que el delineador que usa scipy publicado por BM arriba. Sin embargo, requiere la instalación de mi paquete math3d.

Solo quería mencionar que si se requiere velocidad, envolver el código de unutbu en scipy''s weave.inline y pasar una matriz ya existente como parámetro produce una disminución de 20 veces en el tiempo de ejecución.

El código (en 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*c-b*b-d*d; mat[3*1 + 2] = 2*(c*d-a*b); mat[3*2 + 0] = 2*(b*d-a*c); mat[3*2 + 1] = 2*(c*d+a*b); mat[3*2 + 2] = a*a+d*d-b*b-c*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*b-c*c-d*d, 2*(b*c-a*d), 2*(b*d+a*c)], [2*(b*c+a*d), a*a+c*c-b*b-d*d, 2*(c*d-a*b)], [2*(b*d-a*c), 2*(c*d+a*b), a*a+d*d-b*b-c*c]])

La temporización:

>>> 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]

Una línea, con funciones numpy / scipy.

Usamos lo siguiente:

deje que a sea ​​el vector unitario a lo largo del eje , es decir, a = eje / norma (eje)
y A = I × a sea ​​la matriz sesgada-simétrica asociada a a , es decir, el producto cruzado de la matriz de identidad con una

entonces M = exp (θ A) es la matriz de rotación.

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 (código aquí) calcula la serie taylor del exponencial:
/sum_{k=0}^{20} /frac{1}{k!} (θ A)^k , por lo que es tiempo caro, pero legible y seguro. Puede ser una buena forma si tiene pocas rotaciones pero muchos vectores.

Usando la fórmula de Euler-Rodrigues :

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]

Usar pyquaternion es extremadamente simple; para instalarlo, ejecuta en tu consola: importar pip; pip.main ([''install'', ''pyquaternion''])

Una vez instalada:

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)