Orbit determination

The basic orbit determination task is to determine the classical orbital elements or Keplerian elements, ${\displaystyle a,e,i,\Omega ,\omega ,\nu }$, from the orbital state vectors [${\displaystyle {\vec {r}},{\vec {v}}}$], of an orbiting body with respect to the reference frame of its central body. The central bodies are the sources of the gravitational forces, like the Sun, Earth, Moon and other planets. The orbiting bodies, on the other hand, include planets around the Sun, artificial satellites around the Earth, and spacecraft around planets. Newton's laws of motion well explain the trajectory of an orbiting body, known as Keplerian orbit.

The steps of orbit determination from one state vector are summarized as follows:

• Compute the specific angular momentum ${\displaystyle {\vec {h}}}$ of the orbiting body from its state vector:
  $${\displaystyle {\vec {h}}={\vec {r}}\times {\vec {v}}=\left|{\vec {h}}\right|{\vec {k}}=h{\vec {k}},}$$

  
  where ${\displaystyle {\vec {k}}}$ is the unit vector of the z-axis of the orbital plane. The specific angular momentum is a constant vector for an orbiting body, with its direction perpendicular to the orbital plane of the orbiting body.
• Compute the ascending node vector ${\displaystyle {\vec {n}}}$ from ${\displaystyle {\vec {h}}}$, with ${\displaystyle {\vec {K}}}$ representing the unit vector of the Z-axis of the reference plane, which is perpendicular to the reference plane of the central body:
  $${\displaystyle {\vec {n}}={\vec {K}}\times {\vec {h}}.}$$

  
  The ascending node vector is a vector pointing from the central body to the ascending node of the orbital plane of the orbiting body. Since the line of ascending node is the line of intersection between the orbital plane and the reference plane, it is perpendicular to both the normal vectors of the reference plane (${\displaystyle {\vec {K}}}$) and the orbital plane (${\displaystyle {\vec {k}}}$ or ${\displaystyle {\vec {h}}}$). Therefore, the ascending node vector can be defined by the cross product of these two vectors.
• Compute the eccentricity vector ${\displaystyle {\vec {e}}}$ of the orbit. The eccentricity vector has the magnitude of the eccentricity, ${\displaystyle e}$, of the orbit, and points to the direction of the periapsis of the orbit. This direction is often defined as the x-axis of the orbital plane and has a unit vector ${\displaystyle {\vec {i}}}$. According to the law of motion, it can be expressed as:
  $${\displaystyle {\begin{aligned}{\vec {e}}&={{\vec {v}}\times {\vec {h}} \over {\mu }}-{{\vec {r}} \over {\left|{\vec {r}}\right|}}=e{\vec {i}}\\&=\left({{\left|{\vec {v}}\right|}^{2} \over {\mu }}-{1 \over {\left|{\vec {r}}\right|}}\right){\vec {r}}-{{\vec {r}}\cdot {\vec {v}} \over {\mu }}{\vec {v}}\\&={\frac {1}{\mu }}\left[\left({{\left|{\vec {v}}\right|}^{2}}-{\mu \over {\left|{\vec {r}}\right|}}\right){\vec {r}}-{({\vec {r}}\cdot {\vec {v}})}{\vec {v}}\right]\end{aligned}}}$$

  
  $${\displaystyle e=\left|{\vec {e}}\right|}$$

  
  where ${\displaystyle \mu =GM}$ is the standard gravitational parameter for the central body of mass ${\displaystyle M}$, and ${\displaystyle G}$ is the universal gravitational constant.
• Compute the semi-latus rectum ${\displaystyle p}$ of the orbit, and its semi-major axis ${\displaystyle a}$ (if it is not a parabolic orbit, where ${\displaystyle e=1}$ and ${\displaystyle a}$ is undefined or defined as infinity):
  $${\displaystyle p={\frac {h^{2}}{\mu }}=a(1-e^{2})}$$

  
  $${\displaystyle a={\frac {p}{1-e^{2}}},}$$

  
  (if ${\displaystyle e\neq 1}$).
• Compute the inclination ${\displaystyle i}$ of the orbital plane with respect to the reference plane:
  $${\displaystyle {\begin{aligned}\cos(i)&={\frac {{\vec {K}}\cdot {\vec {h}}}{h}}={\frac {h_{K}}{h}}\\\Rightarrow i&=\arccos \left({\frac {{\vec {K}}\cdot {\vec {h}}}{h}}\right),&i\in [0,180^{\circ }],\end{aligned}}}$$

  
  where ${\displaystyle h_{K}}$ is the Z-coordinate of ${\displaystyle {\vec {h}}}$ when it is projected to the reference frame.
• Compute the longitude of ascending node ${\displaystyle \Omega }$, which is the angle between the ascending line and the X-axis of the reference frame:
  $${\displaystyle {\begin{aligned}\cos(\Omega )&={\frac {{\vec {I}}\cdot {\vec {n}}}{n}}={\frac {n_{I}}{n}}=\cos(360-\Omega )\\\Rightarrow \Omega &=\arccos \left({\frac {{\vec {I}}\cdot {\vec {n}}}{n}}\right)=\Omega _{0},{\text{ or }}\\\Rightarrow \Omega &=360^{\circ }-\Omega _{0},{\text{ if }}n_{J}<0,\\\end{aligned}}}$$

  
  where ${\displaystyle n_{I}}$ and ${\displaystyle n_{J}}$ are the X- and Y- coordinates, respectively, of ${\displaystyle {\vec {n}}}$, in the reference frame.
  Notice that ${\displaystyle \cos(A)=\cos(-A)=\cos(360-A)=C}$, but ${\displaystyle \arccos(C)}$ is defined only in [0,180] degrees. So ${\displaystyle \arccos(C)}$ is ambiguous in that there are two angles, ${\displaystyle A}$ and ${\displaystyle 360-A}$ in [0,360], who have the same ${\displaystyle \cos }$ value. It could actually return the angle ${\displaystyle A}$ or ${\displaystyle 360-A}$. Therefore, we have to make the judgment based on the sign of the Y-coordinate of the vector in the plane where the angle is measured. In this case, ${\displaystyle n_{J}}$ can be used for such judgment.
• Compute the argument of periapsis ${\displaystyle \omega }$, which is the angle between the periapsis and the ascending line:
  $${\displaystyle {\begin{aligned}\cos(\omega )&={\frac {{\vec {n}}\cdot {\vec {e}}}{ne}}=\cos(360-\omega )\\\Rightarrow \omega &=\arccos \left({\frac {{\vec {n}}\cdot {\vec {e}}}{ne}}\right)=\omega _{0},{\text{ or }}\\\Rightarrow \omega &=360^{\circ }-\omega _{0},{\text{ if }}e_{K}<0,\\\end{aligned}}}$$

  
  where ${\displaystyle e_{K}}$ is the Z-coordinate of ${\displaystyle {\vec {e}}}$ in the reference frame.
• Compute the true anomaly ${\displaystyle \nu }$ at epoch, which is the angle between the position vector and the periapsis at the particular time ('epoch') of observation:
  $${\displaystyle {\begin{aligned}\cos(\nu )&={\frac {{\vec {e}}\cdot {\vec {r}}}{er}}=\cos(360-\nu )\\\Rightarrow \nu &=\arccos \left({\frac {{\vec {e}}\cdot {\vec {r}}}{er}}\right)=\nu _{0},{\text{ or }}\\\Rightarrow \nu &=360^{\circ }-\nu _{0},{\text{ if }}{\vec {r}}\cdot {\vec {v}}<0.\\\end{aligned}}}$$

  
  The sign of ${\displaystyle {\vec {r}}\cdot {\vec {v}}}$ can be used to check the quadrant of ${\displaystyle \nu }$ and correct the ${\displaystyle \arccos }$ angle, because it has the same sign as the fly-path angle ${\displaystyle \phi }$. And, the sign of the fly-path angle is always positive when ${\displaystyle \nu \in [0,180^{\circ }]}$, and negative when ${\displaystyle \nu \in [180^{\circ },360^{\circ }]}$.[1] Both are related by ${\displaystyle h=rv\sin(90-\phi )}$ and ${\displaystyle {\vec {r}}\cdot {\vec {v}}=rv\cos(90-\phi )=h\tan(\phi )}$.
• Optionally, we may compute the argument of latitude ${\displaystyle u=\omega +\nu }$ at epoch, which is the angle between the position vector and the ascending line at the particular time:
  $${\displaystyle {\begin{aligned}\cos(u)&={\frac {{\vec {n}}\cdot {\vec {r}}}{nr}}=\cos(360-u)\\\Rightarrow u&=\arccos \left({\frac {{\vec {n}}\cdot {\vec {r}}}{nr}}\right)=u_{0},{\text{ or }}\\\Rightarrow u&=360^{\circ }-u_{0},{\text{ if }}r_{K}<0,\\\end{aligned}}}$$

  
  where ${\displaystyle r_{K}}$ is the Z-coordinate of ${\displaystyle {\vec {r}}}$ in the reference frame.