Understanding result of combining 3D transformations with `pytransform3d`

Question

I feel like I am missing something very basic here. Question is specific to pytransform3d. I am most likely missing a basic understanding of 3D transformations and/or the pytransform3d library.

Here is the setup:

Coordinate system O is inside the cube, and coordinate systems A and B are on the corners as shown.

I have two transformation matrices O2A and O2B (named in the style used in pytransform3d). Since I don't know how to visualize rotation matrices, I convert them into intrinsic Euler XYZ rotations. Using the translations and the intrinsic Euler XYZ rotations, I can visually confirm that the transformation matrices O2A and O2B are correct.

Now, I want to calculate the transformation matrix A2B. Using numpy, I come up with A2B = inv(O2A) @ O2B which I understand and am happy with. Again, just to verify, I can convert into intrinsic Euler XYZ rotations and check visually on the graph.

Now, according to pytransform3d docs of concat and invert_transform, this should be the same as A2B = pytr.concat(pytr.invert_transform(O2A), O2B) but it is not. If I look at the translations and intrinsic Euler XYZ rotations of this result, they are not what I would expect. Looking at this docs page, it seems like it should be the same but I am having difficulty understanding why.

Question: In what sense is the resulting A2B is a transformation from A to B? Perhaps I need to further post-process the result before extracting the translations and intrinsic Euler XYZ rotations?

Below I have some basic Python snippet to show the values:

import numpy as np
from math import degrees, radians
from numpy.linalg import inv
import pytransform3d.rotations as pyrot
import pytransform3d.transformations as pytr

def print_rp(X2Y):
    r = list(map(degrees, pyrot.intrinsic_euler_xyz_from_active_matrix(X2Y[0:3,0:3])))
    p = X2Y[0:3,3:4].ravel().tolist()
    print(f"intrinsic Euler XYZ rotations {r}")
    print(f"translations {p}")

print("# O2A")
O2A = pytr.transform_from(
    R = np.block([[1,0,0],[0,1,0],[0,0,1]]),
    p = [0.05,0.05,-0.05]
)
print_rp(O2A)

print("# O2B")
O2B = pytr.transform_from(
    R = np.block([[-1,-0,-0],[0,-1,0],[0,0,1]]),
    p = [0.05,-0.05,-0.05],
)
print_rp(O2B)

print("# A2B (should work but doesn't)")
A2B = pytr.concat(pytr.invert_transform(O2A), O2B)
print_rp(A2B)

print("# A2B (works but do not know why)")
A2B = pytr.concat(O2B, pytr.invert_transform(O2A))
print_rp(A2B)

print("# A2B (works as expected)")
A2B = inv(O2A) @ O2B
print_rp(A2B)

And the output:

# O2A
intrinsic Euler XYZ rotations [0.0, 0.0, 0.0]
translations [0.05, 0.05, -0.05]
# O2B
intrinsic Euler XYZ rotations [0.0, 0.0, 180.0]
translations [0.05, -0.05, -0.05]
# A2B (should work but doesn't)
intrinsic Euler XYZ rotations [0.0, 0.0, 180.0]
translations [0.1, 0.0, 0.0]
# A2B (works but do not know why)
intrinsic Euler XYZ rotations [0.0, 0.0, 180.0]
translations [0.0, -0.1, 0.0]
# A2B (works as expected)
intrinsic Euler XYZ rotations [0.0, 0.0, 180.0]
translations [0.0, -0.1, 0.0]