Search code examples
pythonrotationvtkmayavirotational-matrices

Correct conversion from rotation matrix to (pitch, roll, yaw) for Mayavi/Vtk


My problem

I want to rotate a mayavi.mlab.imshow object with a 3x3 rotation matrix. The only method I could find for rotating this object is through setting the object's actor.orientation to [pitch, roll, yaw] (this order is inherited from vtk). My only problem is that I cannot find a way to convert a rotation matrix to the parameters requested by mayavi.

How do I rotate an object in mayavi with a rotation matrix or what transformation should I use to obtain the correct (pitch, roll and yaw) used by Mayavi/Vtk?

A near solution

I have found some code over here to transform a rotation matrix to different types of Euler angles (according to the order of rotation). Correct me at this point if I am wrong in assuming Euler angles is equivalent to pitch, roll, yaw. I have tried all the different conversions, but failed to find a correct one.

Trying every combination

I tested all the different transformations by rotating x, y and z vectors with my rotation matrix and testing the paramaters on the mayavi.mlab.imshow object. I used all the available transforms on both R and transpose(R) along with using the Euler the parameters in all 9 available orders, but could not find a correct combination:

import pylab as pl
import cameraTools #my own lib
from mayavi import mlab

im = pl.imread('dice.png', format='png')[:,:,0]*255 #1 color channel
rot = pl.r_[30, 80, 230]

R_orig = cameraTools.composeRotation(*(rot*pl.pi/180))
RList = [R_orig, R_orig.T]
for ii, inOrder in enumerate(['sxyz','sxzx','syxz','szxz','rzyx','rxzx','rzxy','rzxz','sxyx','syzx','syxy','szyx','rxyx','rxzy','ryxy','rxyz','sxzy','syzy','szxy','szyz','ryzx','ryzy','ryxz','rzyz']):
    tries = 0
    for outOrder in [[0,1,2], [0,2,1], [1, 0, 2], [1, 2, 0], [2, 0, 1], [2, 1, 0]]:
        for R in  RList:
            for vector, color in zip([[800, 0, 0], [0, 800, 0], [0, 0, 800]],
                                     [(1., 0., 0.), (0., 1., 0.), (0., 0., 1.)]):
                c = pl.c_[[0, tries*1000, ii*1000]]
                if ii == 0 and tries == 0: vector = pl.r_[vector]*5 #point out the first vector
                lin = R_orig.dot(pl.c_[[0,0,0], vector]) + c
                mlab.plot3d(*lin,
                            color = color,
                            tube_radius=5)

            lin3D    = mlab.imshow(im, colormap="gray")
            rxyz = pl.array(cameraTools.euler_from_matrix(R, inOrder))*180/pi
            i,j,k = outOrder
            lin3D.actor.orientation = [rxyz[i], rxyz[j], rxyz[k]]
            lin3D.actor.position    = c.flatten()
            tries +=1


mlab.draw()
mlab.show() 

Output from Mayavi, the top left item is the starting point.

Mayavi output


Solution

  • Sorry, it seems I did not concentrate hard enough. The answer is in row 3, column 5 with syxz input order and '[1,0,2]' output order. I now use the following function to convert a Rotation matrix to the required Euler angles:

    def rotationToVtk(R):
        '''
        Concert a rotation matrix into the Mayavi/Vtk rotation paramaters (pitch, roll, yaw)
        '''
        def euler_from_matrix(matrix):
            """Return Euler angles (syxz) from rotation matrix for specified axis sequence.
            :Author:
              `Christoph Gohlke <http://www.lfd.uci.edu/~gohlke/>`_
    
            full library with coplete set of euler triplets (combinations of  s/r x-y-z) at
                http://www.lfd.uci.edu/~gohlke/code/transformations.py.html
    
            Note that many Euler angle triplets can describe one matrix.
            """
            # epsilon for testing whether a number is close to zero
            _EPS = np.finfo(float).eps * 5.0
    
            # axis sequences for Euler angles
            _NEXT_AXIS = [1, 2, 0, 1]
            firstaxis, parity, repetition, frame = (1, 1, 0, 0) # ''
    
            i = firstaxis
            j = _NEXT_AXIS[i+parity]
            k = _NEXT_AXIS[i-parity+1]
    
            M = np.array(matrix, dtype='float', copy=False)[:3, :3]
            if repetition:
                sy = np.sqrt(M[i, j]*M[i, j] + M[i, k]*M[i, k])
                if sy > _EPS:
                    ax = np.arctan2( M[i, j],  M[i, k])
                    ay = np.arctan2( sy,       M[i, i])
                    az = np.arctan2( M[j, i], -M[k, i])
                else:
                    ax = np.arctan2(-M[j, k],  M[j, j])
                    ay = np.arctan2( sy,       M[i, i])
                    az = 0.0
            else:
                cy = np.sqrt(M[i, i]*M[i, i] + M[j, i]*M[j, i])
                if cy > _EPS:
                    ax = np.arctan2( M[k, j],  M[k, k])
                    ay = np.arctan2(-M[k, i],  cy)
                    az = np.arctan2( M[j, i],  M[i, i])
                else:
                    ax = np.arctan2(-M[j, k],  M[j, j])
                    ay = np.arctan2(-M[k, i],  cy)
                    az = 0.0
    
            if parity:
                ax, ay, az = -ax, -ay, -az
            if frame:
                ax, az = az, ax
            return ax, ay, az
        r_yxz = pl.array(euler_from_matrix(R))*180/pi
        r_xyz = r_yxz[[1, 0, 2]]
        return r_xyz