The known datapoints are
$\left[ \begin{matrix} n_x & s_x & a_x & p_x \\ n_y & s_y & a_y & p_y \\ n_z & s_z & a_z & p_z \\ 0 & 0 & 0 & 1 \end{matrix}\right]$
where
$$\begin{align} n_x &= -\sin\phi_1\sin\phi_2 + \cos\phi_1\cos\phi_2\cos\theta_1 \\ n_y &= \phantom{-}\cos\phi_1\sin\phi_2+ \sin\phi_1\cos\phi_2\cos\theta_1\\ n_z &= -\cos\phi_2\sin\theta_1 \\[6pt] s_x &= -\sin\phi_1\cos\phi_2 - \cos\phi_1\sin\phi_2\cos\theta_1 \\ s_y &= \phantom{-}\cos\phi_1\cos\phi_2 - \sin\phi_1\sin\phi_2\cos\theta_1 \\ s_z &= \phantom{-}\sin\phi_2\sin\theta_1 \\[6pt] a_x &= \cos\phi_1\sin\theta_1 \\ a_y &= \sin\phi_1\sin\theta_1 \\ a_z &= \cos\theta_1 \\[6pt] p_x &= h_x + r_1\cos\phi_1\cos\theta_1 +r_2\left(- \sin\phi_1\sin\phi_2 + \cos\phi_1\cos\phi_2\cos\theta_1\right) + d_2\cos\phi_1\sin\theta_1 \\ p_y &= h_y + r_1\sin\phi_1\cos\theta_1 + r_2\left(\phantom{-}\cos\phi_1\sin\phi_2 + \cos\phi_2\sin\phi_1\cos\theta_1\right) + d_2\sin\phi_1\sin\theta_1 \\ p_z &= h_z - r_1\sin\theta_1 - r_2\cos\phi_2\sin\theta_1 + d_2\cos\theta_1 \end{align}$$
$px$, $py$, and $pz$ can be simplified to
$$\begin{align} p_x &= h_x + r_1\cos\phi_1\cos\theta_1 +r_2 n_x + d_2 a_x\\ p_y &= h_y + r_1\sin\phi_1\cos\theta_1 +r_2 n_y + d_2 a_y\\ p_z &= h_z - r_1\sin\theta_1 + r_2 n_z+ d_2a_z \end{align}$$
I've gotten my simplification to here but I keep finding it hard to find an elegant way to remove the other sine and cosine terms so that we're able to solve the unknown $(h_x, h_y, h_z, r_1, r_2, d_2)$ relative to 6 or more datapoints matrices.