Converting orbital elements from equatorial reference frame to ecliptic I have a bunch of orbital elements (inclination, longitude of the ascending node, mean anomaly, argument of periapsis) that are all referenced from an equatorial reference frame. But the program I am trying to use requires the orbital elements in an ecliptic frame. Is there an easy/convenient way to convert between the two frames?
I also have the RA and Dec for the orbit pole, also in the equatorial frame. I am aware of some rotation matrices that allow me to convert between equatorial and ecliptic frames—should I just apply the matrix to the RA and Dec of the orbit pole? Is there a way to turn that into orbital elements?
Or should I turn the equatorial elements into a state vector (Cartesian position and velocity) and then apply the matrix to that, and then convert the ecliptic vector into orbital elements?
What is the easiest way to go about this?
The orbital elements $\omega$, $i$ and $\Omega$ are Euler angles in the sequence $(3, 1, 3)$.
The easiest way to transform them is to convert them to a representation that's easier to manipulate, e.g. unit quaternions or a matrix, apply the required transformation then convert back to Euler angles.
A useful guide to converting between the systems is given by James Diebel (2006) &quot;Representing Attitude: Euler Angles, Unit Quaternions, and Rotation Vectors&quot;.
For example, using the formulae in §8.10, the rotation matrix for the orbital elements is given by:
$$R_{313}(\Omega,i,\omega)=R_3(\Omega) R_1(i) R_3(\omega) =
\begin{bmatrix}
c_\Omega c_\omega - s_\Omega c_i s_\omega &amp; c_\Omega s_\omega + s_\Omega c_i c_\omega &amp; s_\Omega s_i \\
-s_\Omega c_\omega - c_\Omega c_i s_\omega &amp; -s_\Omega s_\omega + c_\Omega c_i c_\omega &amp; c_\Omega s_i \\
s_i s_\omega &amp; -s_i c_\omega &amp; c_i
\end{bmatrix}$$
Where $s_\omega = \sin \omega$ and $c_\omega = \cos \omega$, etc.
You can then apply the transformation from equatorial to ecliptic coordinates (or whatever other transformation you want to do) to this matrix. The equatorial and ecliptic coordinates are related by a rotation about the x-axis by $\epsilon$, the obliquity of the ecliptic, so combine the rotation matrix using matrix multiplication. Depending on the precise details of the coordinate systems of the input elements and the ones your software expects, you may also need to invert the axes, which is another matrix multiplication.
This yields a new matrix $R'$ with elements $r_{ij}'$. You then convert this back to Euler angles:
$$
\begin{bmatrix}
\Omega' \\ i' \\ \omega'
\end{bmatrix}
=
\begin{bmatrix}
\operatorname{arctan2}(r_{13}', r_{23}') \\
\arccos(r_{33}') \\
\operatorname{arctan2}(r_{31}', -r_{32}')
\end{bmatrix}
$$
If $i' \in \{0, \pi\}$ then $r_{13}' = r_{23}' = r_{31}' = r_{32}' = 0$ and the above formula for $\Omega'$ and $\omega'$ no longer works as it requires evaluating $\operatorname{arctan2}(0, 0)$. In this case, the angles $\Omega'$ and $\omega'$ are in the same plane and therefore cannot be uniquely separated. Choosing $\omega' = 0$, the transformation from the rotation matrix to Euler angles becomes:
$$
\begin{bmatrix}
\Omega' \\ i' \\ \omega'
\end{bmatrix}
=
\begin{bmatrix}
\operatorname{arctan2}(-r_{21}', r_{11}') \\
\arccos(r_{33}') \\
0
\end{bmatrix}
$$

Note that the two-argument arctangent function in the above uses the convention $\operatorname{arctan2}(y, x)$ giving the angle between the positive x-axis and the point $(x, y)$. Some systems such as Microsoft Excel, which provides the function as ATAN2(x, y) have the arguments reversed, consult the manual for the software you are using.
Depending on whether you pre- or post-multiply the rotation matrix with a vector, the rotation matrix elements may be transposed from the ones given above, in which case you need to switch the element indices in the reverse transformation.