# Astrodynamics/Orbit Determination

According to basic geometry, any three points can be fit by an elipse. This means that if we have at least three observations of a satellite's position, we can determine the entire orbit. However, the process is not an easy one.

RADAR is one of the powerful tools available to observers on earth. A radar installation can determine the position and the velocity of a satellite, with the position being in the topocentric-horizon coordinate system. Once the radar observation has been made, the coordinates can be transformed into the geocentric-equatorial system as such:

$\bold{r} = \bold{r}_e + \bold{p}$

Where re is the position vector from the center of the earth to the radar installation, and p is the position vector of the satellite from the radar installation.

## Determination from a Position and Velocity Vector

A radar installation can measure the position and the velocity of the satellite in a single observation. However, this observation is typically going to be made using the polar variety of the topocentric-horizon system (azimuth and elevation). The position of the satellite is given in terms of a line-of-sight vector, ρ, and the velocity vector is the derivative of that, ρ' We can find the distance from the center of the earth (converting to the geocentric-horizon system) by:

$\bold{r} = \bold{r}_e + \bold{\rho}$

And we can find the velocity vector v by:

$\bold{v} = \bold{\rho}' + \bold{\omega}_e \times \bold{r}$

Where ωe is the angular velocity vector of the earth.

Once we have the v and r vectors in the geocentric-equatorial system, we can find all the classical orbital elements.

## Determination from 3 Position Vectors

If we have three observations, with each observation representing a position vector, we can determine the orbit of the satellite. We assume that all three position vectors are coplanar, because if they are not it is not a proper orbit and these equations will not yield proper results.

Given three position vectors, r1, r2, and r3, we can apply the gibbsian method to determine the orbit.

### Gibbsian Method

For the gibbsian method, we first construct three new vectors, N, D, and S (not to be confused the S unit vector from the topocentric-horizon system). We construct these three new vectors as follows:

$\bold{N} = r_3\bold{r}_1 \times \bold{r}_2 + r_1\bold{r}_2\times\bold{r}_3 + r_2\bold{r}_3\times\bold{r}_1$
$\bold{D} = \bold{r}_1\times\bold{r}_2 + \bold{r}_2\times\bold{r}_3 + \bold{r}_3\times\bold{r}_1$
$\bold{S} = (r_2 - r_3)\bold{r}_1 + (r_3 - r_1)\bold{r}_2 + (r_1 - r_2)\bold{r}_3$

Now, we can find our p and e quantities from these vectors:

$p = \bold{D} = \bold{N}$

N and D always have the same direction, so:

$\bold{N} \cdot \bold{D} = ND$
$p = \frac{N}{D}$
$e = \frac{S}{D}$

Now, we can determine the unit vectors from the perifocal coordinate system, P, Q, and W from the three vectors we already have:

$\bold{Q} = \frac{\bold{S}}{S}$
$\bold{W} = \frac{\bold{N}}{N}$
$\bold{P} = \bold{Q} \times \bold{W}$

Now, we have the p and e values, and we also have the necessarily perifocal unit vectors to define the orbit.

### Determining Velocity Vectors

If we have the three position vectors, we can determine three corresponding velocity vectors at each position. If we have the D, N, and S vectors, we can define a new vector B such that:

$\bold{B} = \bold{D} \times \bold{r}_i$

Where the i subscript corresponds to one of the three position vectors. We can also form a scalar value L, in order to make the calculations easier:

$L = \sqrt{\frac{\mu}{DN}}$

Now that we have our 4 vectors and the scalar L, we can find the velocity at position i as:

$\bold{v}_i = \frac{L}{r_i}\bold{B} + L\bold{S}$

## Determination from 3 Observations

Radar observations aren't always available, but optical observations typically are. Optical measurements are typically made with the right ascension and declination coordinate system, so we will have 6 values: α1, δ1, α2, δ2, α3, δ3. We can define three unit vectors Li as:

$\bold{L}_i = \begin{bmatrix}\cos(\delta_i)\cos(\alpha_i) \\ \cos(\delta_i)\sin(\alpha_i) \\ \sin(\delta_i) \end{bmatrix}$

With the L vectors, we can find the geocentric-equatorial position vector:

$\bold{r}_i = \rho_i \bold{L}_i + \bold{r}_e$

Dotting this equation with itself gives us a second equation:

$r^2 = \rho^2 + r_e^2 + 2 \rho \bold{L} \cdot \bold{r}_e$

We can solve these two equations to find the two unknowns ρ and r. With r, we can find the position vector r along the direction of the vector L.

### Finding the Velocity Vector

We can define the derivative of the L vectors using the Lagrange interpolation formula:

$\bold{L}' = \frac{2t-t_2-t_3}{(t_1-t_2)(t_1-t_3)}\bold{L}_1 + \frac{2t - t_1 - t_3}{(t_2-t_1)(t_2-t_3)}\bold{L}_2 + \frac{2t-t_1-t_3}{(t_3-t_1)(t_3 - t_2)}\bold{L}_3$
$\bold{L}'' = \frac{2}{(t_1-t_2)(t_1-t_3)}\bold{L}_1 + \frac{2}{(t_2-t_1)(t_2-t_3)}\bold{L}_2 + \frac{2}{(t_3-t_1)(t_3 - t_2)}\bold{L}_3$

We can use these values to find:

$\bold{v} = \rho' \bold{L} + \rho \bold{L}' + \bold{r}_e$

## Differential Orbit Correction

We've seen that with a small number of observations we can determine the complete orbit of a satellite. However, the initial measurements that we have taken can be prone to errors. Any errors that we have in our measurements are made all the more large and glaring because we have relied on so few observations. However, if we increase our number of observations, we can use new values to improve our original "estimate" of the orbit. To do this, we use a method called differential orbit correction.