def calc_rotation_matrix(q1, q2, ref_q1, ref_q2):
ref_nv = np.cross(ref_q1, ref_q2)
q_nv = np.cross(q1, q2)
if min(norm(ref_nv), norm(q_nv)) == 0.: # avoid 0 degree including angle
return np.identity(3)
axis = np.cross(ref_nv, q_nv)
angle = rad2deg(acos(ref_nv.dot(q_nv) / (norm(ref_nv) * norm(q_nv))))
R1 = axis_angle_to_rotation_matrix(axis, angle)
rot_ref_q1, rot_ref_q2 = R1.dot(ref_q1), R1.dot(ref_q2) # rotate ref_q1,2 plane to q1,2 plane
cos1 = max(min(q1.dot(rot_ref_q1) / (norm(rot_ref_q1) * norm(q1)), 1.), -1.) # avoid math domain error
cos2 = max(min(q2.dot(rot_ref_q2) / (norm(rot_ref_q2) * norm(q2)), 1.), -1.)
angle1 = rad2deg(acos(cos1))
angle2 = rad2deg(acos(cos2))
angle = (angle1 + angle2) / 2.
axis = np.cross(rot_ref_q1, q1)
R2 = axis_angle_to_rotation_matrix(axis, angle)
R = R2.dot(R1)
return R
评论列表
文章目录