Resolution of a rotation Can the rotation of a point about some arbitrary axis $\vec{n}$ be described by some rotation of the point about z-axis followed by some rotation about y axis?

Say we want to rotate point $(x_0, y_0, z_0)$ to $(x_1, y_1, z_1)$ where $\rho= \sqrt{x_0^2+ y_0^2+ z_0^2}= \sqrt{x_1^2+ y_1^2+ z_1^2}$.  First rotate about the z-axis so that $(x_0, y_0, z_0)$ is rotated to $(r, 0, z_0)$ where $r= \sqrt{x_0^2+ y_0^2}$.  That can be done with a matrix of the form $\begin{pmatrix}a & b & 0 \\ -b & a & 0 \\ 0 & 0 & 1\end{pmatrix}$.  We want $\begin{pmatrix}a & b & 0 \\ -b & a & 0 \\ 0 & 0 & 1\end{pmatrix}\begin{pmatrix}x_0 \\y_0 \\z_0\end{pmatrix}= \begin{pmatrix}ax_0+ by_0 \\ -bx_0+ ay_0 \\ z_0\end{pmatrix}= \begin{pmatrix} r\\ 0 \\ z_0\end{pmatrix}$.  So we have $ax_0+ by_0= r$ and $-bx_0+ ay_0= 0$ and we want to solve for a and b.  Multiply the first equation by $x_0$ to get $ax_0^2+ bx_0y_0= rx_0$.  Multiply the second equation by $y_0$ to get $-bx_0y_0+ ay_0^2= 0$.  Adding those eliminates b: $a(x_0^2+ y_0^2)= ar^2= rx_0$.  $a= \frac{x_0}{r}$.  Then $-bx_0+ ay_0= 0$ becomes $-bx_0+ \frac{x_0y_0}{r}= 0$ so $bx_0= \frac{x_0y_0}{r}$ and $b= \frac{y_0}{r}$.
The first matrix, rotating about t the z-axis, is $\begin{pmatrix}\frac{x_0}{r} & \frac{y_0}{r} & 0 \\ -\frac{y_0}{r} & \frac{x_0}{r} & 0 \\ 0 & 0 & 1\end{pmatrix}$.
Now we want to rotate around the y-axis so that $\begin{pmatrix} r \\ 0 \\ z_0\end{pmatrix}$ rotates to $\begin{pmatrix} 0 \\ 0 \\ \rho\end{pmatrix}$. Rotating about the y-axis requires a matrix of the form $\begin{pmatrix} a & 0 & b \\ 0 & 1 & 0 \\ -b & 0 & a\end{pmatrix}$.  Then we must have $\begin{pmatrix} a & 0 & b \\ 0 & 1 & 0 \\ -b & 0 & a\end{pmatrix}\begin{pmatrix} r \\ 0 \\ z_0\end{pmatrix}= \begin{pmatrix}ar+ bz_0 \\ 0 \\ -br+ az_0 \end{pmatrix}= \begin{pmatrix} 0 \\ 0 \\ \rho \end{pmatrix}$ where $\rho= \sqrt{x_0^2+ y_0^2+ z_0^2}= \sqrt{r^2+ z_0^2}$.  Now we have the two equations $ar+ bz_0= 0$ and $-br+ az_0= \rho$ that we want to solve for a and b.
Multiply the first equation by r to get $ar^2+ brz_0= 0$ and multiply the second equation by $z_0$ to get $-brz_0+ az_0^2= \rho z_0$.  Adding the two equations eliminate b: $ar^2+ az_0^2= a(r^2+ z_0^2)= a\rho^2= \rho z_0$ so $a= \frac{z_0}{\rho}$.  Now, $ar+ bz_0= 0$ becomes $\frac{r z_0}{\rho}+ bz_0= 0$ so $bz_0= -\frac{r z_0}{\rho}$ and $b= -\frac{r}{\rho}$.
The second matrix, rotating around the y-axis, is $\begin{pmatrix}\frac{z_0}{\rho} & 0 & -\frac{r}{\rho} \\ 0 & 1 & 0 \\ \frac{r}{\rho} & 0 & \frac{z_0}{\rho} \end{pmatrix}$.
