scipy.spatial.transform.

RigidTransform#

class scipy.spatial.transform.RigidTransform#

Rigid transform in 3 dimensions.

This class provides an interface to initialize from and represent rigid transforms (rotation and translation) in 3D space. In different fields, this type of transform may be referred to as “pose” (especially in robotics), “extrinsic parameters”, or the “model matrix” (especially in computer graphics), but the core concept is the same: a rotation and translation describing the orientation of one 3D coordinate frame relative to another. Mathematically, these transforms belong to the Special Euclidean group SE(3), which encodes rotation (SO(3)) plus translation.

The following operations on rigid transforms are supported:

  • Application on vectors

  • Transformation composition

  • Transformation inversion

  • Transformation indexing

Note that coordinate systems must be right-handed. Because of this, this class more precisely represents proper rigid transformations in SE(3) rather than rigid transforms in E(3) more generally [1].

Indexing within a transform is supported since multiple transforms can be stored within a single RigidTransform instance.

To create RigidTransform objects use from_... methods (see examples below). RigidTransform(...) is not supposed to be instantiated directly.

For rigorous introductions to rigid transforms, see [2], [3], and [4].

Attributes:
single

Whether this instance represents a single transform.

rotation

Return the rotation component of the transform.

translation

Return the translation component of the transform.

Methods

__len__

Return the number of transforms in this object.

__getitem__

Extract transform(s) at given index(es) from this object.

__mul__

Compose this transform with the other.

__pow__

Compose this transform with itself n times.

from_matrix(cls, matrix)

Initialize from a 4x4 transformation matrix.

from_rotation(cls, rotation)

Initialize from a rotation, without a translation.

from_translation(cls, translation)

Initialize from a translation numpy array, without a rotation.

from_components(cls, translation, rotation)

Initialize a rigid transform from translation and rotation components.

from_exp_coords(cls, exp_coords)

Initialize from exponential coordinates of transform.

from_dual_quat(cls, dual_quat, *[, scalar_first])

Initialize from a unit dual quaternion.

as_matrix(self)

Return a copy of the matrix representation of the transform.

as_components(self)

Return the translation and rotation components of the transform, where the rotation is applied first, followed by the translation.

as_exp_coords(self)

Return the exponential coordinates of the transform.

as_dual_quat(self, *[, scalar_first])

Return the dual quaternion representation of the transform.

concatenate(cls, transforms)

Concatenate a sequence of RigidTransform objects into a single object.

apply(self, vector[, inverse])

Apply the transform to a vector.

inv(self)

Invert this transform.

identity(cls[, num])

Initialize an identity transform.

Notes

Added in version 1.16.0.

References

[4]

Kevin M. Lynch and Frank C. Park, “Modern Robotics: Mechanics, Planning, and Control” Chapter 3.3, 2017, Cambridge University Press. https://hades.mech.northwestern.edu/images/2/25/MR-v2.pdf#page=107.31

[5]

Paul Furgale, “Representing Robot Pose: The good, the bad, and the ugly”, June 9, 2014. https://rpg.ifi.uzh.ch/docs/teaching/2024/FurgaleTutorial.pdf

Examples

A RigidTransform instance can be initialized in any of the above formats and converted to any of the others. The underlying object is independent of the representation used for initialization.

Notation Conventions and Composition

The notation here largely follows the convention defined in [5]. When we name transforms, we read the subscripts from right to left. So tf_A_B represents a transform A <- B and can be interpreted as:

  • the coordinates and orientation of B relative to A

  • the transformation of points from B to A

  • the pose of B described in A’s coordinate system

tf_A_B
   ^ ^
   | |
   | --- from B
   |
   ----- to A

When composing transforms, the order is important. Transforms are not commutative, so in general tf_A_B * tf_B_C is not the same as tf_B_C * tf_A_B. Transforms are composed and applied to vectors right-to-left. So (tf_A_B * tf_B_C).apply(p_C) is the same as tf_A_B.apply(tf_B_C.apply(p_C)).

When composed, transforms should be ordered such that the multiplication operator is surrounded by a single frame, so the frame “cancels out” and the outside frames are left. In the example below, B cancels out and the outside frames A and C are left. Or to put it another way, A <- C is the same as A <- B <- C.

              ----------- B cancels out
              |      |
              v      v
tf_A_C = tf_A_B * tf_B_C
            ^          ^
            |          |
            ------------ to A, from C are left

When we notate vectors, we write the subscript of the frame that the vector is defined in. So p_B means the point p defined in frame B. To transform this point from frame B to coordinates in frame A, we apply the transform tf_A_B to the vector, lining things up such that the notated B frames are next to each other and “cancel out”.

           ------------ B cancels out
           |         |
           v         v
p_A = tf_A_B.apply(p_B)
         ^
         |
         -------------- A is left

Visualization

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

The following function can be used to plot transforms with Matplotlib by showing how they transform the standard x, y, z coordinate axes:

>>> import matplotlib.pyplot as plt
>>> colors = ("#FF6666", "#005533", "#1199EE")  # Colorblind-safe RGB
>>> def plot_transformed_axes(ax, tf, name=None, scale=1):
...     r = tf.rotation
...     t = tf.translation
...     loc = np.array([t, t])
...     for i, (axis, c) in enumerate(zip((ax.xaxis, ax.yaxis, ax.zaxis),
...                                       colors)):
...         axlabel = axis.axis_name
...         axis.set_label_text(axlabel)
...         axis.label.set_color(c)
...         axis.line.set_color(c)
...         axis.set_tick_params(colors=c)
...         line = np.zeros((2, 3))
...         line[1, i] = scale
...         line_rot = r.apply(line)
...         line_plot = line_rot + loc
...         ax.plot(line_plot[:, 0], line_plot[:, 1], line_plot[:, 2], c)
...         text_loc = line[1]*1.2
...         text_loc_rot = r.apply(text_loc)
...         text_plot = text_loc_rot + t
...         ax.text(*text_plot, axlabel.upper(), color=c,
...                 va="center", ha="center")
...     ax.text(*tf.translation, name, color="k", va="center", ha="center",
...             bbox={"fc": "w", "alpha": 0.8, "boxstyle": "circle"})

Defining Frames

Let’s work through an example.

First, define the “world frame” A, also called the “base frame”. All frames are the identity transform from their own perspective.

>>> tf_A = Tf.identity()

We will visualize a new frame B in A’s coordinate system. So we need to define the transform that converts coordinates from frame B to frame A (A <- B).

Physically, let’s imagine constructing B from A by:

  1. Rotating A by +90 degrees around its x-axis.

  2. Translating the rotated frame 2 units in A’s -x direction.

From A’s perspective, B is at [-2, 0, 0] and rotated +90 degrees about the x-axis, which is exactly the transform A <- B.

>>> t_A_B = np.array([-2, 0, 0])
>>> r_A_B = R.from_euler('xyz', [90, 0, 0], degrees=True)
>>> tf_A_B = Tf.from_components(t_A_B, r_A_B)

Let’s plot these frames.

>>> fig, ax = plt.subplots(subplot_kw={"projection": "3d"})
>>> plot_transformed_axes(ax, tf_A, name="tfA")     # A plotted in A
>>> plot_transformed_axes(ax, tf_A_B, name="tfAB")  # B plotted in A
>>> ax.set_title("A, B frames with respect to A")
>>> ax.set_aspect("equal")
>>> ax.figure.set_size_inches(6, 5)
>>> plt.show()
../../_images/scipy-spatial-transform-RigidTransform-1_00_00.png

Now let’s visualize a new frame C in B’s coordinate system. Let’s imagine constructing C from B by:

  1. Translating B by 2 units in its +z direction.

  2. Rotating B by +30 degrees around its z-axis.

>>> t_B_C = np.array([0, 0, 2])
>>> r_B_C = R.from_euler('xyz', [0, 0, 30], degrees=True)
>>> tf_B_C = Tf.from_components(t_B_C, r_B_C)

In order to plot these frames from a consistent perspective, we need to calculate the transform between A and C. Note that we do not make this transform directly, but instead compose intermediate transforms that let us get from C to A:

>>> tf_A_C = tf_A_B * tf_B_C  # A <- B <- C

Now we can plot these three frames from A’s perspective.

>>> fig, ax = plt.subplots(subplot_kw={"projection": "3d"})
>>> plot_transformed_axes(ax, tf_A, name="tfA")     # A plotted in A
>>> plot_transformed_axes(ax, tf_A_B, name="tfAB")  # B plotted in A
>>> plot_transformed_axes(ax, tf_A_C, name="tfAC")  # C plotted in A
>>> ax.set_title("A, B, C frames with respect to A")
>>> ax.set_aspect("equal")
>>> ax.figure.set_size_inches(6, 5)
>>> plt.show()
../../_images/scipy-spatial-transform-RigidTransform-1_01_00.png

Transforming Vectors

Let’s transform a vector from A, to B and C. To do this, we will first invert the transforms we already have from B and C, to A.

>>> tf_B_A = tf_A_B.inv()  # B <- A
>>> tf_C_A = tf_A_C.inv()  # C <- A

Now we can define a point in A and use the above transforms to get its coordinates in B and C:

>>> p1_A = np.array([1, 0, 0])  # +1 in x_A direction
>>> p1_B = tf_B_A.apply(p1_A)
>>> p1_C = tf_C_A.apply(p1_A)
>>> print(p1_A)  # Original point 1 in A
[1 0 0]
>>> print(p1_B)  # Point 1 in B
[3. 0. 0.]
>>> print(p1_C)  # Point 1 in C
[ 2.59807621 -1.5       -2.        ]

We can also do the reverse. We define a point in C and transform it to A:

>>> p2_C = np.array([0, 1, 0])  # +1 in y_C direction
>>> p2_A = tf_A_C.apply(p2_C)
>>> print(p2_C)  # Original point 2 in C
[0 1 0]
>>> print(p2_A)  # Point 2 in A
[-2.5       -2.         0.8660254]

Plot the frames with respect to A again, but also plot these two points:

>>> fig, ax = plt.subplots(subplot_kw={"projection": "3d"})
>>> plot_transformed_axes(ax, tf_A, name="tfA")     # A plotted in A
>>> plot_transformed_axes(ax, tf_A_B, name="tfAB")  # B plotted in A
>>> plot_transformed_axes(ax, tf_A_C, name="tfAC")  # C plotted in A
>>> ax.scatter(p1_A[0], p1_A[1], p1_A[2], color=colors[0])  # +1 x_A
>>> ax.scatter(p2_A[0], p2_A[1], p2_A[2], color=colors[1])  # +1 y_C
>>> ax.set_title("A, B, C frames and points with respect to A")
>>> ax.set_aspect("equal")
>>> ax.figure.set_size_inches(6, 5)
>>> plt.show()
../../_images/scipy-spatial-transform-RigidTransform-1_02_00.png

Switching Base Frames

Up to this point, we have been visualizing frames from A’s perspective. Let’s use the transforms we defined to visualize the frames from C’s perspective.

Now C is the “base frame” or “world frame”. All frames are the identity transform from their own perspective.

>>> tf_C = Tf.identity()

We’ve already defined the transform C <- A, and can obtain C <- B by inverting the existing transform B <- C.

>>> tf_C_B = tf_B_C.inv()  # C <- B

This lets us plot everything from C’s perspective:

>>> fig, ax = plt.subplots(subplot_kw={"projection": "3d"})
>>> plot_transformed_axes(ax, tf_C, name="tfC")     # C plotted in C
>>> plot_transformed_axes(ax, tf_C_B, name="tfCB")  # B plotted in C
>>> plot_transformed_axes(ax, tf_C_A, name="tfCA")  # A plotted in C
>>> ax.scatter(p1_C[0], p1_C[1], p1_C[2], color=colors[0])
>>> ax.scatter(p2_C[0], p2_C[1], p2_C[2], color=colors[1])
>>> ax.set_title("A, B, C frames and points with respect to C")
>>> ax.set_aspect("equal")
>>> ax.figure.set_size_inches(6, 5)
>>> plt.show()
../../_images/scipy-spatial-transform-RigidTransform-1_03_00.png