# NiBabel

### Access a cacophony of neuro-imaging file formats

# `eulerangles`

Module implementing Euler angle rotations and their conversions

See:

- https://en.wikipedia.org/wiki/Rotation_matrix

- https://en.wikipedia.org/wiki/Euler_angles

- http://mathworld.wolfram.com/EulerAngles.html

See also: Representing Attitude with Euler Angles and Quaternions: A Reference (2006) by James Diebel. A cached PDF link last found here:

http://citeseerx.ist.psu.edu/viewdoc/summary?doi=10.1.1.110.5134

Euler’s rotation theorem tells us that any rotation in 3D can be described by 3 angles. Let’s call the 3 angles the Euler angle vector and call the angles in the vector \(alpha\), \(beta\) and \(gamma\). The vector is [ \(alpha\), \(beta\). \(gamma\) ] and, in this description, the order of the parameters specifies the order in which the rotations occur (so the rotation corresponding to \(alpha\) is applied first).

In order to specify the meaning of an Euler angle vector we need to specify the axes around which each of the rotations corresponding to \(alpha\), \(beta\) and \(gamma\) will occur.

There are therefore three axes for the rotations \(alpha\), \(beta\) and \(gamma\); let’s call them \(i\) \(j\), \(k\).

Let us express the rotation \(alpha\) around axis i as a 3 by 3 rotation matrix A. Similarly \(beta\) around j becomes 3 x 3 matrix B and \(gamma\) around k becomes matrix G. Then the whole rotation expressed by the Euler angle vector [ \(alpha\), \(beta\). \(gamma\) ], R is given by:

    R = np.dot(G, np.dot(B, A))

See http://mathworld.wolfram.com/EulerAngles.html

The order \(G B A\) expresses the fact that the rotations are performed in the order of the vector (\(alpha\) around axis i = A first).

To convert a given Euler angle vector to a meaningful rotation, and a rotation matrix, we need to define:

- the axes i, j, k

- whether a rotation matrix should be applied on the left of a vector to be transformed (vectors are column vectors) or on the right (vectors are row vectors).

- whether the rotations move the axes as they are applied (intrinsic rotations) - compared the situation where the axes stay fixed and the vectors move within the axis frame (extrinsic)

- the handedness of the coordinate system

See: https://en.wikipedia.org/wiki/Rotation_matrix#Ambiguities

We are using the following conventions:

- axes i, j, k are the z, y, and x axes respectively. Thus an Euler angle vector [ \(alpha\), \(beta\). \(gamma\) ] in our convention implies a \(alpha\) radian rotation around the z axis, followed by a \(beta\) rotation around the y axis, followed by a \(gamma\) rotation around the x axis.

- the rotation matrix applies on the left, to column vectors on the right, so if R is the rotation matrix, and v is a 3 x N matrix with N column vectors, the transformed vector set vdash is given by `vdash = np.dot(R, v)`.

- extrinsic rotations - the axes are fixed, and do not move with the rotations.

- a right-handed coordinate system

The convention of rotation around `z`, followed by rotation around `y`, followed by rotation around `x`, is known (confusingly) as “xyz”, pitch-roll-yaw, Cardan angles, or Tait-Bryan angles.

`angle_axis2euler`(theta, vector[, is_normalized]) | Convert angle, axis pair to Euler angles
`euler2angle_axis`([z, y, x]) | Return angle, axis corresponding to these Euler angles
`euler2mat`([z, y, x]) | Return matrix for rotations around z, y and x axes
`euler2quat`([z, y, x]) | Return quaternion corresponding to these Euler angles
`mat2euler`(M[, cy_thresh]) | Discover Euler angle vector from 3x3 matrix
`quat2euler`(q) | Return Euler angles corresponding to quaternion q

## angle_axis2euler

nibabel.eulerangles.angle_axis2euler(theta, vector, is_normalized=False)

  Convert angle, axis pair to Euler angles

  Parameters:
    theta scalar

      angle of rotation

    vector 3 element sequence

      vector specifying axis for rotation.

    is_normalized bool, optional

      True if vector is already normalized (has norm of 1). Default False

  Returns:
    z scalar
    y scalar
    x scalar

      Rotations in radians around z, y, x axes, respectively

  Notes

  It’s possible to reduce the amount of calculation a little, by combining parts of the `angle_axis2mat` and `mat2euler` functions, but the reduction in computation is small, and the code repetition is large.

  Examples

      >>> z, y, x = angle_axis2euler(0, [1, 0, 0])
      >>> np.allclose((z, y, x), 0)
      True

## euler2angle_axis

nibabel.eulerangles.euler2angle_axis(z=0, y=0, x=0)

  Return angle, axis corresponding to these Euler angles

  Uses the z, then y, then x convention above

  Parameters:
    z scalar

      Rotation angle in radians around z-axis (performed first)

    y scalar

      Rotation angle in radians around y-axis

    x scalar

      Rotation angle in radians around x-axis (performed last)

  Returns:
    theta scalar

      angle of rotation

    vector array shape (3,)

      axis around which rotation occurs

  Examples

      >>> theta, vec = euler2angle_axis(0, 1.5, 0)
      >>> print(theta)
      1.5
      >>> np.allclose(vec, [0, 1, 0])
      True

## euler2mat

nibabel.eulerangles.euler2mat(z=0, y=0, x=0)

  Return matrix for rotations around z, y and x axes

  Uses the z, then y, then x convention above

  Parameters:
    z scalar

      Rotation angle in radians around z-axis (performed first)

    y scalar

      Rotation angle in radians around y-axis

    x scalar

      Rotation angle in radians around x-axis (performed last)

  Returns:
    M array shape (3,3)

      Rotation matrix giving same rotation as for given angles

  Notes

  The direction of rotation is given by the right-hand rule (orient the thumb of the right hand along the axis around which the rotation occurs, with the end of the thumb at the positive end of the axis; curl your fingers; the direction your fingers curl is the direction of rotation). Therefore, the rotations are counterclockwise if looking along the axis of rotation from positive to negative.

  Examples

      >>> zrot = 1.3 # radians
      >>> yrot = -0.1
      >>> xrot = 0.2
      >>> M = euler2mat(zrot, yrot, xrot)
      >>> M.shape == (3, 3)
      True

  The output rotation matrix is equal to the composition of the individual rotations

      >>> M1 = euler2mat(zrot)
      >>> M2 = euler2mat(0, yrot)
      >>> M3 = euler2mat(0, 0, xrot)
      >>> composed_M = np.dot(M3, np.dot(M2, M1))
      >>> np.allclose(M, composed_M)
      True

  You can specify rotations by named arguments

      >>> np.all(M3 == euler2mat(x=xrot))
      True

  When applying M to a vector, the vector should column vector to the right of M. If the right hand side is a 2D array rather than a vector, then each column of the 2D array represents a vector.

      >>> vec = np.array([1, 0, 0]).reshape((3,1))
      >>> v2 = np.dot(M, vec)
      >>> vecs = np.array([[1, 0, 0],[0, 1, 0]]).T # giving 3x2 array
      >>> vecs2 = np.dot(M, vecs)

  Rotations are counter-clockwise.

      >>> zred = np.dot(euler2mat(z=np.pi/2), np.eye(3))
      >>> np.allclose(zred, [[0, -1, 0],[1, 0, 0], [0, 0, 1]])
      True
      >>> yred = np.dot(euler2mat(y=np.pi/2), np.eye(3))
      >>> np.allclose(yred, [[0, 0, 1],[0, 1, 0], [-1, 0, 0]])
      True
      >>> xred = np.dot(euler2mat(x=np.pi/2), np.eye(3))
      >>> np.allclose(xred, [[1, 0, 0],[0, 0, -1], [0, 1, 0]])
      True

## euler2quat

nibabel.eulerangles.euler2quat(z=0, y=0, x=0)

  Return quaternion corresponding to these Euler angles

  Uses the z, then y, then x convention above

  Parameters:
    z scalar

      Rotation angle in radians around z-axis (performed first)

    y scalar

      Rotation angle in radians around y-axis

    x scalar

      Rotation angle in radians around x-axis (performed last)

  Returns:
    quat array shape (4,)

      Quaternion in w, x, y z (real, then vector) format

  Notes

  We can derive this formula in Sympy using:

  1. Formula giving quaternion corresponding to rotation of theta radians about arbitrary axis: http://mathworld.wolfram.com/EulerParameters.html

  2. Generated formulae from 1.) for quaternions corresponding to theta radians rotations about `x, y, z` axes

  3. Apply quaternion multiplication formula - https://en.wikipedia.org/wiki/Quaternions#Hamilton_product - to formulae from 2.) to give formula for combined rotations.

## mat2euler

nibabel.eulerangles.mat2euler(M, cy_thresh=None)

  Discover Euler angle vector from 3x3 matrix

  Uses the conventions above.

  Parameters:
    M array-like, shape (3,3)
    cy_thresh None or scalar, optional

      threshold below which to give up on straightforward arctan for estimating x rotation. If None (default), estimate from precision of input.

  Returns:
    z scalar
    y scalar
    x scalar

      Rotations in radians around z, y, x axes, respectively

  Notes

  If there was no numerical error, the routine could be derived using Sympy expression for z then y then x rotation matrix, which is:

      [                       cos(y)*cos(z),                       -cos(y)*sin(z),         sin(y)],
      [cos(x)*sin(z) + cos(z)*sin(x)*sin(y), cos(x)*cos(z) - sin(x)*sin(y)*sin(z), -cos(y)*sin(x)],
      [sin(x)*sin(z) - cos(x)*cos(z)*sin(y), cos(z)*sin(x) + cos(x)*sin(y)*sin(z),  cos(x)*cos(y)]

  with the obvious derivations for z, y, and x

  > z = atan2(-r12, r11) y = asin(r13) x = atan2(-r23, r33)

  Problems arise when cos(y) is close to zero, because both of:

      z = atan2(cos(y)*sin(z), cos(y)*cos(z))
      x = atan2(cos(y)*sin(x), cos(x)*cos(y))

  will be close to atan2(0, 0), and highly unstable.

  The `cy` fix for numerical instability below is from: Graphics Gems IV, Paul Heckbert (editor), Academic Press, 1994, ISBN: 0123361559. Specifically it comes from EulerAngles.c by Ken Shoemake, and deals with the case where cos(y) is close to zero:

  See: http://www.graphicsgems.org/

  The code appears to be licensed (from the website) as “can be used without restrictions”.

## quat2euler

nibabel.eulerangles.quat2euler(q)

  Return Euler angles corresponding to quaternion q

  Parameters:
    q 4 element sequence

      w, x, y, z of quaternion

  Returns:
    z scalar

      Rotation angle in radians around z-axis (performed first)

    y scalar

      Rotation angle in radians around y-axis

    x scalar

      Rotation angle in radians around x-axis (performed last)

  Notes

  It’s possible to reduce the amount of calculation a little, by combining parts of the `quat2mat` and `mat2euler` functions, but the reduction in computation is small, and the code repetition is large.