In principle you need two rotations and a matrix-multiplication. Call the recentered matrices Q ad S and P as T.
Then find the rotation s, which rotates S to "triangular" form of its coordinates, so that the coordinates of the first point becomes $[x_{s,1},0,0]$, of the second becomes $[x_{s,2},y_{s,2},0]$. This can be done with your ansatz of using unknowns, if you assume three rotations (for each pair of planes 1 rotation, and gives a 3x3-matrix whose 3'rd column is zero and describes thus a triangle in a plane)
Then do the same with another rotation t which rotes T in the same manner.
Then $T = S \cdot s \cdot t^{-1} = S \cdot r$ where r is the matrix for the complete rotation.
A bit of pseudocode (code in my proprietary MatMate-tool):
S = Q - Meanzl(Q) // do the translation to the origin T = P - Meanzl(P) // get the rotation-matrices which rotate some matrix to triangular shape ts = gettrans(S,"tri") tt = gettrans(T,"tri") // make one rotation-metrix. the quote-symbol means transposition tr = ts * tt' // difference should be zero CHK = T - S*tr
Example:
We assume
S and
T being centered. Let
$ \small \text{ S =} \begin{bmatrix} 14.469944&22.964690&-7.581631\\ -15.275348&5.923432&23.720255\\ 0.805404&-28.888122&-16.138624 \end{bmatrix} $ and
$ \small \text{ T =} \begin{bmatrix} 22.808501&2.515200&16.361035\\ 8.393637&-5.071089&-27.109127\\ -31.202138&2.555889&10.748092 \end{bmatrix} $ Now you can solve for a rotation in
y/z-plane, which makes the entry in
$S_{1,3}=0$. The rotation-parameters are some cos/sin-values. Apply this and you get
$\small \text{ S}^{(1)} = \begin{bmatrix} 14.469944&24.183840&0.000000\\ -15.275348&-1.811476&24.381471\\ 0.805404&-22.372364&-24.381471 \end{bmatrix}$ Now you can solve for a rotation in
x/y-plane, which makes the entry in
$S_{1,2}=0$. The rotation-parameters are some other cos/sin-values. Apply this and you get
$\small \text{S}^{(2)} = \begin{bmatrix} 28.182218&0.000000&0.000000\\ -9.397482&12.178056&24.381471\\ -18.784736&-12.178056&-24.381471 \end{bmatrix}$ After that a third set of rotation-parameters cos/sin-values make
S triangular. It looks then like this
$\small \text{ S}^{(3)} = \begin{bmatrix} 28.182218&0.000000&0.000000\\ -9.397482&27.253645&0.000000\\ -18.784736&-27.253645&0.000000 \end{bmatrix}$ Because the center of your original triangle was moved to the origin, the last column (the z-coordinates) are zero/not needed, since 3 points can always be placed in a plane.
From the three rotation with their cos/sin-values you can construct a 3x3-rotation-matrix, say s.
The same can be done using the matrix T leading to a rotation-matrix t. If S and T describe in fact the same triangle, only rotated, the results are equal: $T^{(3)}=S^{(3)}$. Then you can use the nice fact, that the inverse of a rotation-matrix is just its transpose, such that with
$\small \text{ tr =} \begin{bmatrix} 0.205215&0.860645&0.466022\\ 0.946329&-0.295966&0.129867\\ 0.249696&0.414359&-0.875190 \end{bmatrix}$ we get $ T = S \cdot tr $
In principle this can also be solved using the concept of
pseudoinverses: we demand
$ tr = S^{-1} \cdot T $ . But because
S has reduced rank the inverse means to divide by zero. Using SVD-decomposition (see wikipedia) the pseudoinverse can be computed if the inverse of the diagonal matrix of the SVD-factors is used (where zeros are simply left zero). This should lead to the same solution.