scipy.spatial.transform.Rotation.

as_euler#

Rotation.as_euler(seq, degrees=False)[source]#

Represent as Euler angles.

Any orientation can be expressed as a composition of 3 elementary rotations. Once the axis sequence has been chosen, Euler angles define the angle of rotation around each respective axis [1].

The algorithm from [2] has been used to calculate Euler angles for the rotation about a given sequence of axes.

Euler angles suffer from the problem of gimbal lock [3], where the representation loses a degree of freedom and it is not possible to determine the first and third angles uniquely. In this case, a warning is raised, and the third angle is set to zero. Note however that the returned angles still represent the correct rotation.

Parameters:
seqstring, length 3

3 characters belonging to the set {‘X’, ‘Y’, ‘Z’} for intrinsic rotations, or {‘x’, ‘y’, ‘z’} for extrinsic rotations [1]. Adjacent axes cannot be the same. Extrinsic and intrinsic rotations cannot be mixed in one function call.

degreesboolean, optional

Returned angles are in degrees if this flag is True, else they are in radians. Default is False.

Returns:
anglesndarray, shape (3,) or (N, 3)

Shape depends on shape of inputs used to initialize object. The returned angles are in the range:

  • First angle belongs to [-180, 180] degrees (both inclusive)

  • Third angle belongs to [-180, 180] degrees (both inclusive)

  • Second angle belongs to:

    • [-90, 90] degrees if all axes are different (like xyz)

    • [0, 180] degrees if first and third axes are the same (like zxz)

Notes

Array API Standard Support

as_euler has experimental support for Python Array API Standard compatible backends in addition to NumPy. Please consider testing these features by setting an environment variable SCIPY_ARRAY_API=1 and providing CuPy, PyTorch, JAX, or Dask arrays as array arguments. The following combinations of backend and device (or other capability) are supported.

Library

CPU

GPU

NumPy

n/a

CuPy

n/a

PyTorch

JAX

Dask

n/a

See Support for the array API standard for more information.

References

[2]

Bernardes E, Viollet S (2022) Quaternion to Euler angles conversion: A direct, general and computationally efficient method. PLoS ONE 17(11): e0276302. https://doi.org/10.1371/journal.pone.0276302

Examples

>>> from scipy.spatial.transform import Rotation as R
>>> import numpy as np

Represent a single rotation:

>>> r = R.from_rotvec([0, 0, np.pi/2])
>>> r.as_euler('zxy', degrees=True)
array([90.,  0.,  0.])
>>> r.as_euler('zxy', degrees=True).shape
(3,)

Represent a stack of single rotation:

>>> r = R.from_rotvec([[0, 0, np.pi/2]])
>>> r.as_euler('zxy', degrees=True)
array([[90.,  0.,  0.]])
>>> r.as_euler('zxy', degrees=True).shape
(1, 3)

Represent multiple rotations in a single object:

>>> r = R.from_rotvec([
... [0, 0, np.pi/2],
... [0, -np.pi/3, 0],
... [np.pi/4, 0, 0]])
>>> r.as_euler('zxy', degrees=True)
array([[ 90.,   0.,   0.],
       [  0.,   0., -60.],
       [  0.,  45.,   0.]])
>>> r.as_euler('zxy', degrees=True).shape
(3, 3)