Reputation: 2226
I have three vectors in 3D a,b,c. Now I want to calculate a rotation r that when applied to a yields a result parallel to b. Then the rotation r needs to be applied to c.
How do I do this in python? Is it possible to do this with numpy/scipy?
Upvotes: 1
Views: 1820
Reputation: 55
Using numpy:
import numpy as np
from numpy import linalg as LA
from math import pi,cos,sin,acos
def rotate(v,angle=0,ref=np.array([0,0,1]),deg=False):
'''Rotates a vector a given angle respect the
given reference. Option: deg=False (default)'''
if(abs(angle) < 1e-5):
return v
if(deg):
angle = angle*pi/180
# Define rotation reference system
ref = versor(ref) # rotation axis
# n1 & n2 are perpendicular to ref, and to themselves
n1 = versor(np.cross(ref,np.array([-ref[1],ref[2],ref[0]])))
n2 = np.cross(ref,n1)
vp = np.inner(v,ref)*ref # parallel to ref vector
vn = v-vp # perpendicular to ref vector
vn_abs = LA.norm(vn)
if(vn_abs < 1e-5):
return v
alp = acos(np.inner(vn,n1)/vn_abs) # angle between vn & n1
if(triprod(ref,n1,vn) < 0):
alp = -alp # correct if necesary
return vp+vn_abs*(n1*cos(alp+angle)+n2*sin(alp+angle))
def triprod(a,b,c):
'''Triple product of vectors: a·(b x c)'''
return np.inner(a,np.cross(b,c))
def versor(v):
'''Unitary vector in the direction of the one given as input'''
v = np.array(v)
return v/LA.norm(v)
###### Test ################################################
a = np.array([3,4,1])
b = np.array([0,-1,2])
c = np.array([1,1,5])
r = acos(np.inner(a,b)/(LA.norm(a)*LA.norm(b)))
ref = versor(np.cross(a,b))
print rotate(c,angle=r,ref=ref)
print r
print ref
Upvotes: 1
Reputation: 7430
I'll assume the "geometry library for python" already answered in the comments on the question. So once you have a transformation that takes 'a' parallel to 'b', you'll just apply it to 'c'
The vectors 'a' and 'b' uniquely define a plane. Each vector has a canonical representation as a point difference from the origin, so you have three points: the head of 'a', the head of 'b', and the origin. First compute this plane. It will have an equation in the form Ax + By + Cz = 0.
A normal vector to this plane defines both the axis of rotation and the sign convention for the direction of rotation. All you need is one normal vector to the plane, since they're all collinear. You can solve for such a vector by picking at two non-collinear vectors in the plane and taking the dot product with the normal vector. This gives a pair of linear equations in two variables that you can solve with standard methods such as Cramer's rule. In all of these manipulations, if any of A, B, or C are zero, you have a special case to handle.
The angle of the rotation is given by the cosine relation for the dot product of 'a' and 'b' and their lengths. The sign of the angle is determined by the triple product of 'a', 'b', and the normal vector. Now you've got all the data to construct a rotation matrix in one of the many canonical forms you can look up.
Upvotes: 0