(a1, a2, a3, X=True, Y=False, ZYX=False, deg=False)
return M
NAME:
      euler_matrix_new

    PURPOSE:
      computes the Euler matrix of an arbitrary rotation described
      by 3 Euler angles
      correct bugs present in Euler_Matrix

    CALLING SEQUENCE:
      result = euler_matrix_new (a1, a2, a3 [,X, Y, ZYX, DEG ])

    INPUTS:
      a1, a2, a3 = Euler angles, scalar
                   (in radian by default, in degree if DEG is set)
                   all the angles are measured counterclockwise
                   correspond to x, y, zyx-conventions (see Goldstein)
                   the default is x

    KEYWORD PARAMETERS:
       DEG : if set the angle are measured in degree

       X :       rotation a1 around original Z
                 rotation a2 around interm   X
                 rotation a3 around final    Z
                 DEFAULT,  classical mechanics convention

       Y :       rotation a1 around original Z
                 rotation a2 around interm   Y
                 rotation a3 around final    Z
                 quantum mechanics convention (override X)

       ZYX :     rotation a1 around original Z
                 rotation a2 around interm   Y
                 rotation a3 around final    X
                 aeronautics convention (override X)

       * these last three keywords are obviously mutually exclusive *

    OUTPUTS:
       result is a 3x3 matrix

    USAGE:
      if vec is an Nx3 array containing N 3D vectors,
      vec # euler_matrix_new(a1,a2,a3,/Y) will be the rotated vectors


    MODIFICATION HISTORY:
       March 2002, EH, Caltech, rewritting of euler_matrix

       convention   euler_matrix_new           euler_matrix
      X:       M_new(a,b,c,/X)  =  M_old(-a,-b,-c,/X) = Transpose( M_old(c, b, a,/X))
      Y:       M_new(a,b,c,/Y)  =  M_old(-a, b,-c,/Y) = Transpose( M_old(c,-b, a,/Y))
    ZYX:       M_new(a,b,c,/Z)  =  M_old(-a, b,-c,/Z)
NAME:
      euler_matrix_new
NAME : euler_matrix_new
def euler_matrix_new(a1, a2, a3, X=True, Y=False, ZYX=False, deg=False):
    """
    NAME:
      euler_matrix_new

    PURPOSE:
      computes the Euler matrix of an arbitrary rotation described
      by 3 Euler angles
      correct bugs present in Euler_Matrix

    CALLING SEQUENCE:
      result = euler_matrix_new (a1, a2, a3 [,X, Y, ZYX, DEG ])

    INPUTS:
      a1, a2, a3 = Euler angles, scalar
                   (in radian by default, in degree if DEG is set)
                   all the angles are measured counterclockwise
                   correspond to x, y, zyx-conventions (see Goldstein)
                   the default is x

    KEYWORD PARAMETERS:
       DEG : if set the angle are measured in degree

       X :       rotation a1 around original Z
                 rotation a2 around interm   X
                 rotation a3 around final    Z
                 DEFAULT,  classical mechanics convention

       Y :       rotation a1 around original Z
                 rotation a2 around interm   Y
                 rotation a3 around final    Z
                 quantum mechanics convention (override X)

       ZYX :     rotation a1 around original Z
                 rotation a2 around interm   Y
                 rotation a3 around final    X
                 aeronautics convention (override X)

       * these last three keywords are obviously mutually exclusive *

    OUTPUTS:
       result is a 3x3 matrix

    USAGE:
      if vec is an Nx3 array containing N 3D vectors,
      vec # euler_matrix_new(a1,a2,a3,/Y) will be the rotated vectors


    MODIFICATION HISTORY:
       March 2002, EH, Caltech, rewritting of euler_matrix

       convention   euler_matrix_new           euler_matrix
      X:       M_new(a,b,c,/X)  =  M_old(-a,-b,-c,/X) = Transpose( M_old(c, b, a,/X))
      Y:       M_new(a,b,c,/Y)  =  M_old(-a, b,-c,/Y) = Transpose( M_old(c,-b, a,/Y))
    ZYX:       M_new(a,b,c,/Z)  =  M_old(-a, b,-c,/Z)
    """

    t_k = 0
    if ZYX:
        t_k = t_k + 1
    # if X:   t_k = t_k + 1
    if Y:
        t_k = t_k + 1
    if t_k > 1:
        raise ValueError("Choose either X, Y or ZYX convention")

    convert = 1.0
    if deg:
        convert = np.pi / 180.0

    c1 = np.cos(a1 * convert)
    s1 = np.sin(a1 * convert)
    c2 = np.cos(a2 * convert)
    s2 = np.sin(a2 * convert)
    c3 = np.cos(a3 * convert)
    s3 = np.sin(a3 * convert)

    if ZYX:
        m1 = np.array([[c1, -s1, 0], [s1, c1, 0], [0, 0, 1]])  # around   z

        m2 = np.array([[c2, 0, s2], [0, 1, 0], [-s2, 0, c2]])  # around   y

        m3 = np.array([[1, 0, 0], [0, c3, -s3], [0, s3, c3]])  # around   x

    elif Y:
        m1 = np.array([[c1, -s1, 0], [s1, c1, 0], [0, 0, 1]])  # around   z

        m2 = np.array([[c2, 0, s2], [0, 1, 0], [-s2, 0, c2]])  # around   y

        m3 = np.array([[c3, -s3, 0], [s3, c3, 0], [0, 0, 1]])  # around   z

    else:
        m1 = np.array([[c1, -s1, 0], [s1, c1, 0], [0, 0, 1]])  # around   z

        m2 = np.array([[1, 0, 0], [0, c2, -s2], [0, s2, c2]])  # around   x

        m3 = np.array([[c3, -s3, 0], [s3, c3, 0], [0, 0, 1]])  # around   z

    M = np.dot(m3.T, np.dot(m2.T, m1.T))

    return M
def euler_matrix_new ( a1 , a2 , a3 , X = True , Y = False , ZYX = False , deg = False ) : t_k = 0 if ZYX : t_k = t_k + 1 # if X:   t_k = t_k + 1 if Y : t_k = t_k + 1 if t_k > 1 : raise ValueError ( "Choose either X, Y or ZYX convention" ) convert = 1.0 if deg : convert = np . pi / 180.0 c1 = np . cos ( a1 * convert ) s1 = np . sin ( a1 * convert ) c2 = np . cos ( a2 * convert ) s2 = np . sin ( a2 * convert ) c3 = np . cos ( a3 * convert ) s3 = np . sin ( a3 * convert ) if ZYX : m1 = np . array ( [ [ c1 , - s1 , 0 ] , [ s1 , c1 , 0 ] , [ 0 , 0 , 1 ] ] ) # around   z m2 = np . array ( [ [ c2 , 0 , s2 ] , [ 0 , 1 , 0 ] , [ - s2 , 0 , c2 ] ] ) # around   y m3 = np . array ( [ [ 1 , 0 , 0 ] , [ 0 , c3 , - s3 ] , [ 0 , s3 , c3 ] ] ) # around   x elif Y : m1 = np . array ( [ [ c1 , - s1 , 0 ] , [ s1 , c1 , 0 ] , [ 0 , 0 , 1 ] ] ) # around   z m2 = np . array ( [ [ c2 , 0 , s2 ] , [ 0 , 1 , 0 ] , [ - s2 , 0 , c2 ] ] ) # around   y m3 = np . array ( [ [ c3 , - s3 , 0 ] , [ s3 , c3 , 0 ] , [ 0 , 0 , 1 ] ] ) # around   z else : m1 = np . array ( [ [ c1 , - s1 , 0 ] , [ s1 , c1 , 0 ] , [ 0 , 0 , 1 ] ] ) # around   z m2 = np . array ( [ [ 1 , 0 , 0 ] , [ 0 , c2 , - s2 ] , [ 0 , s2 , c2 ] ] ) # around   x m3 = np . array ( [ [ c3 , - s3 , 0 ] , [ s3 , c3 , 0 ] , [ 0 , 0 , 1 ] ] ) # around   z M = np . dot ( m3 . T , np . dot ( m2 . T , m1 . T ) ) return M