Finding rotation axis and angle to align two 3D vector bases
The rotation matrix, $R$, can be found relatively easy, such that,
$$ R\,\vec{a} = \vec{d}; \quad R\,\vec{b} = \vec{e}; \quad R\,\vec{c} = \vec{f}. $$
The rotation matrix/second-order-tensor can be constructed by adding three dyads,
$$ R = \vec{d}\vec{a} + \vec{e}\vec{b} + \vec{f}\vec{c}. $$
When using this notation the multiplication with a vector should be written as a dot product. From this it can also be proven that this expression for $R$ is correct. For example with mapping $\vec{a}$ to $\vec{d}$, the first dyad will return $\vec{d}\vec{a}\cdot\vec{a}$ and since $\vec{a}$ is normalized that is equal to $\vec{d}$. The remaining dot products will yield zero, because the other vectors on the right side of each dyad is orthogonal to $\vec{a}$.
Each dyad can also be written as a 3$\times$3 matrix using the vector direct product. The total 3$\times$3 rotation matrix can then be constructed by adding the three 3$\times$3 matrix from the vector direct products.
The axis of rotation is equal to the eigenvector of the eigenvalue equal to one. For finding this eigenvector you can take a look at the answer to this question. In order to find the axis of rotation you have to calculate the following matrix,
$$ R-R^T = \begin{bmatrix} 0 & \alpha & \beta \\ -\alpha & 0 & \gamma \\ -\beta & -\gamma & 0 \end{bmatrix}, $$
such that the normalized axis of rotation, $\vec{r}$, is equal to,
$$ \vec{r} = \frac{ \begin{bmatrix} -\gamma & \beta & -\alpha \end{bmatrix}^T }{\sqrt{\alpha^2 + \beta^2 + \gamma^2}}. $$
The angle of rotation, $\theta$, can be found according to,
$$ \text{Tr}(R) = 1 + 2\cos\theta, $$
where $\text{Tr}(R)$ is the trace of $R$, the sum of its diagonal elements. Due to the cosine you do have to look at the sign of the angle and direction of $\vec{r}$.