Computing planet positions in Planetoweb

Below you can see, how does Planetoweb computes planets' positions

Meaning of symbols

Source: wikipedia.org
$M$ - mean anomaly
$t$ - time since J2000.0
$P$ - orbital period
$E$ - eccentric anomaly
$e$ - orbital eccentricity
$v$ - real anomaly
$a$ - semi-major axis
$N$ - longitude of ascending node
$w$ - argument of periapsis
$i$ - inclination

Mean anomaly

„Mean anomaly - the angular distance from the periapsis to an imaginary object, which moves witch constant angular velocity.” (source: pl.wikipedia.org (pl))

Planetoweb computes it using this equation:

$M = M_0 + nt$

where $M_0$ means mean anomaly at J2000.0 epoch and $n$ is a change of angle in one day, which is:

$n = {2{\pi} \over P}$

Eccentric anomaly

We can user eccentric anomaly to compute real anomaly:

$v = 2 \arctan \left (\sqrt{{e+1 \over e-1}} \tan {E \over 2} \right )$

Having mean anomaly we can use Kepler's equation to compute eccentric anomaly:

$M = E - e \sin E$

Unfortunately the Kepler's equation is a transcendental equation and we would have to user numerical method to solve it.

„Equation of the centre”

Luckily we can skip computing eccentric anomaly and use „Equation of the center”. It's equal to the difference between the real and mean anomaly. Expanding it to terms of the order of $e^3$ gives us enough accuracy.

Real anomaly

Using „equation of the centre” we get equation on the approximated real anomaly:

$v = M + \left (2e - {e^3 \over 4} \right ) \sin M + {5 \over 4} e^2 \sin 2M + {13 \over 12} e^3 \sin 3M$

Distance to the Sun and coordinates

The centre of Planetoweb's coordinate system is in the centre of Sun. The plane of reference for inclination ($i$) is the plane of Earth's orbit. In our model it's plane:

$\begin{cases}\begin{matrix} x \in \mathbb R \\ y \in \mathbb R \\ z = 0 \end{matrix}\end{cases}$

Having semi-major axis ($a$) we can compute a distance from a planet to the Sun ($r$):

$r = {a(1 - e^2) \over 1 + e \cos v}$

and when we know the longitude of the ascending node ($N$) and the argument of periapsis ($w$) we can compute final coordinates:

$x = r \left( \cos (N) \cos (v+w) - \sin (N) \sin (v+w) \cos (i) \right)$ $y = r \left( \sin (N) \cos (v+w) + \cos (N) \sin (v+w) \cos (i) \right)$ $z = r \sin (v+w) \sin (i)$