IT TIP

3D 벡터의 회전?

itqueen 2020. 11. 25. 21:47
반응형

3D 벡터의 회전?


파이썬 목록과 각도로 두 개의 벡터가 있습니다. 예 :

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

축을 중심으로 v 벡터를 회전 할 때 결과 벡터를 얻는 가장 좋은 / 가장 쉬운 방법은 무엇입니까?

축 벡터가 가리키는 관찰자에 대해 회전은 시계 반대 방향으로 표시되어야합니다. 이것을 오른손 법칙 이라고합니다.


http://vpython.org/contents/docs/visual/VisualIntro.html을 살펴보십시오 .

vector메서드가 있는 클래스를 제공합니다 A.rotate(theta,B). 또한 rotate(A,theta,B)에서 메서드를 호출하지 않으려는 경우 도우미 함수를 제공합니다 A.

http://vpython.org/contents/docs/visual/vector.html


은 Using 오일러 - 로드리게스 공식을 :

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를 따라 단위 벡터 일 , 즉 A = 축 / 표준 (축)
A = I ×는 A가 연결된 스큐 대칭 행렬 과 단위 행렬의 예 외적

다음 M = EXP (θ의 A)는 회전 행렬이다.

from numpy import cross, eye, dot
from scipy.linalg import expm, norm

def M(axis, theta):
    return expm(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]

expm (여기 코드) 는 지수의 테일러 시리즈를 계산
\sum_{k=0}^{20} \frac{1}{k!} (θ A)^k하므로 시간이 많이 걸리지 만 읽기 쉽고 안전합니다. 할 회전이 적고 벡터가 많은 경우 좋은 방법이 될 수 있습니다.


속도가 필요한 경우 unutbu의 코드를 scipy의 weave.inline에 래핑하고 이미 존재하는 행렬을 매개 변수로 전달하면 실행 시간이 20 배 감소한다는 점을 언급하고 싶었습니다.

코드 (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]])

타이밍:

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

엄청나게 빠른 쿼터니언을 사용하는 우아한 방법이 있습니다. 적절하게 벡터화 된 numpy 배열을 사용하여 초당 천만 회전을 계산할 수 있습니다. 여기에있는 numpy에 대한 쿼터니언 확장에 의존 합니다 .

쿼터니언 이론 : 쿼터니언은 일반적으로 q = w + xi + yj + zk'i', 'j', 'k'가 허수 차원 인 1 개의 실수 차원과 3 개의 가상 차원이있는 숫자입니다 . 단위 복소수 'c'가 모든 2d 회전을으로 나타낼 수있는 것처럼 c=exp(i * theta)단위 쿼터니언 'q'는 모든 3d 회전을으로 나타낼 수 있습니다 q=exp(p). 여기서 'p'는 축과 각도로 설정된 순수한 가상 쿼터니언입니다.

축과 각도를 회전축에 의해 가상 차원이 지정되고 크기가 라디안 단위의 회전 각도의 절반으로 지정되는 쿼터니언으로 변환하는 것으로 시작합니다. 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와 회전 축 모두에 대해 실수 성분 w = 0으로 4 개의 요소로 구성된 numpy 배열이 구성됩니다 rot_axis. 그런 다음 축 각도 표현은 정규화 한 다음 원하는 각도의 절반을 곱하여 구성됩니다 theta. 반 각도가 필요한 이유는 여기참조 하십시오 .

지금 사원 수를 생성 vqlog라이브러리를 사용하여, 단위 회전 쿼터니언을 얻을 q지수를 취함으로써.

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

마지막으로 벡터의 회전은 다음 작업으로 계산됩니다.

v_prime = q * vec * np.conjugate(q)

print(v_prime) # quaternion(0.0, 2.7491163, 4.7718093, 1.9162971)

이제 실제 요소를 버리고 회전 된 벡터를 갖게됩니다!

v_prime_vec = v_prime.imag # [2.74911638 4.77180932 1.91629719] as a numpy array

이 방법은 쿼터니언 곱을 q = q1 * q2 * q3 * q4 * ... * qn으로 계산 한 다음 벡터 만 회전 할 수 있으므로 여러 연속 회전을 통해 벡터를 회전해야하는 경우 특히 효율적입니다. v '= q * v * conj (q)를 사용하여 맨 끝에'q '로.

이 방법은 단순히 explog함수 를 사용하여 축 각도 <---> 3d 회전 연산자 사이의 원활한 변환을 제공합니다 (예 log(q)는 축 각도 표현을 반환합니다!). 쿼터니언 곱셈 등이 작동하는 방식에 대한 자세한 설명은 여기를 참조 하십시오.


저는 Python {2,3}을위한 상당히 완벽한 3D 수학 라이브러리를 만들었습니다. 여전히 Cython을 사용하지 않지만 numpy의 효율성에 크게 의존합니다. pip로 여기에서 찾을 수 있습니다.

python[3] -m pip install math3d

또는 내 gitweb http://git.automatics.dyndns.dk/?p=pymath3d.git 및 github에서도 살펴보십시오 : https://github.com/mortlind/pymath3d .

일단 설치되면, 파이썬에서 벡터를 회전하거나 변형 객체의 일부가 될 수있는 방향 객체를 만들 수 있습니다. 예를 들어 다음 코드 스 니펫은 축 [1,2,3]을 중심으로 1rad의 회전을 나타내는 방향을 구성하고이를 벡터 [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)>

위의 BM이 게시 한 scipy를 사용하는 oneliner보다 시간을 할 수있는 한 대략 4 배 정도 더 효율적입니다. 그러나 내 math3d 패키지를 설치해야합니다.


면책 조항 : 나는이 패키지의 작성자입니다.

회전을위한 특수 클래스가 편리 할 수 ​​있지만, 어떤 경우에는 회전 행렬이 필요합니다 (예 : scipy의 affine_transform 함수와 같은 다른 라이브러리 작업을 위해). 모든 사람이 자신의 작은 행렬 생성 기능을 구현하는 것을 방지하기 위해 편리한 회전 행렬 생성 기능을 제공하는 것 이상을 수행하는 작은 순수 Python 패키지가 있습니다. 패키지는 github ( mgen )에 있으며 pip를 통해 설치할 수 있습니다.

pip install mgen

Readme에서 복사 한 사용 예 :

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

행렬은 일반 배열 일 뿐이므로이 패키지를 사용할 때 새로운 데이터 구조가 도입되지 않습니다.


쿼터니언 이론을 사용하여 해결할 수도 있습니다.

def angle_axis_quat(theta, axis):
    """
    Given an angle and an axis, it returns a quaternion.
    """
    axis = np.array(axis) / np.linalg.norm(axis)
    return np.append([np.cos(theta/2)],np.sin(theta/2) * axis)

def mult_quat(q1, q2):
    """
    Quaternion multiplication.
    """
    q3 = np.copy(q1)
    q3[0] = q1[0]*q2[0] - q1[1]*q2[1] - q1[2]*q2[2] - q1[3]*q2[3]
    q3[1] = q1[0]*q2[1] + q1[1]*q2[0] + q1[2]*q2[3] - q1[3]*q2[2]
    q3[2] = q1[0]*q2[2] - q1[1]*q2[3] + q1[2]*q2[0] + q1[3]*q2[1]
    q3[3] = q1[0]*q2[3] + q1[1]*q2[2] - q1[2]*q2[1] + q1[3]*q2[0]
    return q3

def rotate_quat(quat, vect):
    """
    Rotate a vector with the rotation defined by a quaternion.
    """
    # Transfrom vect into an quaternion 
    vect = np.append([0],vect)
    # Normalize it
    norm_vect = np.linalg.norm(vect)
    vect = vect/norm_vect
    # Computes the conjugate of quat
    quat_ = np.append(quat[0],-quat[1:])
    # The result is given by: quat * vect * quat_
    res = mult_quat(quat, mult_quat(vect,quat_)) * norm_vect
    return res[1:]

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

print(rotate_quat(angle_axis_quat(theta, axis), v))
# [2.74911638 4.77180932 1.91629719]

pyquaternion을 사용하는 것은 매우 간단합니다. 설치하려면 (Python에있는 동안) 콘솔에서 실행하십시오.

import pip;
pip.main(['install','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)

I needed to rotate a 3D model around one of the three axes {x, y, z} in which that model was embedded and this was the top result for a search of how to do this in numpy. I used the following simple function:

def rotate(X, theta, axis='x'):
  '''Rotate multidimensional array `X` `theta` degrees around axis `axis`'''
  c, s = np.cos(theta), np.sin(theta)
  if axis == 'x': return np.dot(X, np.array([
    [1.,  0,  0],
    [0 ,  c, -s],
    [0 ,  s,  c]
  ]))
  elif axis == 'y': return np.dot(X, np.array([
    [c,  0,  -s],
    [0,  1,   0],
    [s,  0,   c]
  ]))
  elif axis == 'z': return np.dot(X, np.array([
    [c, -s,  0 ],
    [s,  c,  0 ],
    [0,  0,  1.],
  ]))

Use scipy's Rotation.from_rotvec(). The argument is the rotation vector (a unit vector) multiplied by the rotation angle in rads.

from scipy.spatial.transform import Rotation
from numpy.linalg import norm


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

axis = axis / norm(axis)  # normalize the rotation vector first
rot = Rotation.from_rotvec(theta * axis)

new_v = rot.apply(v)  
print(new_v)    # results in [2.74911638 4.77180932 1.91629719]

There are several more ways to use Rotation based on what data you have about the rotation:

  • from_quat Initialized from quaternions.

  • from_dcm Initialized from direction cosine matrices.

  • from_euler Initialized from Euler angles.


Off-topic note: One line code is not necessarily better code as implied by some users.

참고URL : https://stackoverflow.com/questions/6802577/rotation-of-3d-vector

반응형