 PLANET POSITIONS

When planning deep space missions it is obvious that accurate positions of the planets must be known to plot the interplanetary trajectories. The computational methods presented in this page are those described by Jean Meeus in his book Astronomical Formulae for Calculators, Fourth Edition, Willmann-Bell Inc., 1988.

Julian Date

Julian Date (JD) is the system of time measurement for scientific use by the astronomy community. It is the interval of time in days and fractions of a day since 4713 BC January-1 Greenwich noon. The Julian day begins at Greenwich mean noon, that is at 12:00 Universal Time (UT).

To convert a calendar date to Julian Date, perform the following steps:

Y = year
M = month
D = day (includes hours, minutes & seconds as fraction of day)

If M = 1 or 2, take

Y = Y – 1
M = M + 12

If the date is equal to or after 1582-Oct-15 (i.e. the date is in the Gregorian calendar), calculate

A = INT( Y / 100 )
B = 2 – A + INT( A / 4 )

If the date is before 1582-Oct-15 (i.e. the date is in the Julian calendar), it is not necessary to calculate A and B.

The Julian date is then

JD = INT(365.25 × Y) + INT(30.6001 × (M + 1)) + D + 1720994.5 + B

where B is added only if the date is in the Gregorian calendar.

EXAMPLE

Calculate the JD corresponding to 1976-July-20, 12:00 UT.

Y = 1976
M = 7
D = 20.5

A = INT(1976 / 100) = 19
B = 2 – 19 + INT( 19 / 4 ) = –13

JD = INT(365.25 × 1976) + INT(30.6001 × (7 + 1)) + 20.5 + 1720994.5 – 13 = 2442980.0

The orbital elements of the major planets can be expressed as polynomials of the form

a0 + a1T + a2T2 + a3T3

where T is the time measured in Julian centuries of 36525 ephemeris days from the epoch 1900 January 0.5 ET = JD 2415020.0. In other words,

T = ( JD – 2415020.0 ) / 36525

This quantity is negative before the beginning of the year 1900, positive afterwards. Note that T is expressed in centuries, and thus should be taken with a sufficient number of decimals (at least nine).

The orbital elements are:

L = mean longitude of the planet
a = semimajor axis of the orbit (this element is a constant for each planet)
e = eccentricity of the orbit
i = inclination on the plane of the ecliptic
w = argument of perihelion
W = longitude of ascending node

The longitude of the perihelion can be calculated from

p = w + W

and the planet's mean anomaly is

M = L – p = L – wW

The perihelion and aphelion distances are

Rp = a (1 – e)
Ra = a (1 + e)

The quantities L and p are measured in two different planes, namely from the vernal equinox along the ecliptic to the orbit's ascending node, and then from this node along the orbit.

Table 1 gives the coefficients ai for the orbital elements of the major planets Mercury to Neptune. The elements for Earth are not given in Table 1.

 TABLE 1Elements for the mean equinox of the date a0 a1 a2 a3 MERCURY L 178.179 078 +149 474.070 78 +0.000 3011 a 0.387 0986 e 0.205 614 21 +0.000 020 46 –0.000 000 030 i 7.002 881 +0.001 8608 –0.000 0183 w 28.753 753 +0.370 2806 +0.000 1208 W 47.145 944 +1.185 2083 +0.000 1739 VENUS L 342.767 053 +58 519.211 91 +0.000 3097 a 0.723 3316 e 0.006 820 69 –0.000 047 74 +0.000 000 091 i 3.393 631 +0.001 0058 –0.000 0010 w 54.384 186 +0.508 1861 –0.001 3864 W 75.779 647 +0.899 8500 +0.000 4100 MARS L 293.737 334 +19 141.695 51 +0.000 3107 a 1.523 6883 e 0.093 312 90 +0.000 092 064 –0.000 000 077 i 1.850 333 –0.000 6750 +0.000 0126 w 285.431 761 +1.069 7667 +0.000 1313 +0.000 004 14 W 48.786 442 +0.770 9917 –0.000 0014 –0.000 005 33 JUPITER L 238.049 257 +3036.301 986 +0.000 3347 –0.000 001 65 a 5.202 561 e 0.048 334 75 +0.000 164 180 –0.000 000 4676 –0.000 000 0017 i 1.308 736 –0.005 6961 +0.000 0039 w 273.277 558 +0.559 4317 +0.000 704 05 +0.000 005 08 W 99.443 414 +1.010 5300 +0.000 352 22 –0.000 008 51 SATURN L 266.564 377 +1223.509 884 +0.000 3245 –0.000 0058 a 9.554 747 e 0.055 892 32 –0.000 345 50 –0.000 000 728 +0.000 000 000 74 i 2.492 519 –0.003 9189 –0.000 015 49 +0.000 000 04 w 338.307 800 +1.085 2207 +0.000 978 54 +0.000 009 92 W 112.790 414 +0.873 1951 –0.000 152 18 –0.000 005 31 URANUS L 244.197 470 +429.863 546 +0.000 3160 –0.000 000 60 a 19.218 14 e 0.046 3444 –0.000 026 58 +0.000 000 077 i 0.772 464 +0.000 6253 +0.000 0395 w 98.071 581 +0.985 7650 –0.001 0745 –0.000 000 61 W 73.477 111 +0.498 6678 +0.001 3117 NEPTUNE L 84.457 994 +219.885 914 +0.000 3205 –0.000 000 60 a 30.109 57 e 0.008 997 04 +0.000 006 330 –0.000 000 002 i 1.779 242 –0.009 5436 –0.000 0091 w 276.045 975 +0.325 6394 +0.000 140 95 +0.000 004 113 W 130.681 389 +1.098 9350 +0.000 249 87 –0.000 004 718

EXAMPLE

Calculate the orbital elements of Mars on 1976-July-20, 12:00 UT.

From our previous calculation we have JD 2442980.0; therefore,

T = (2442980.0 – 2415020.0) / 36525 = 0.765503080

Consequently, from Table 1, we find

L = 293.737334 + (19141.69551 × 0.765503080) + (0.0003107 × 0.7655030802) = 14946.764387o = 186.764387o

Similarly,

a = 1.5236883 AU
e = 0.093383330
i = 1.849824o
w = 286.250750o
W = 49.376635o

The longitude of perihelion and the mean anomaly are,

p = 286.250750 + 49.376635 = 335.627385o
M = 186.764387 – 335.627385 = –148.862998 = 211.137002o

Since the plane of the ecliptic is the plane of Earth's orbit around the Sun, for Earth i = 0. The angles w and W are, therefore, not determined. The semimajor axis of Earth is a = 1.0000002 AU. The remaining orbital elements are calculated as follows:

L = 99.69668 + 36000.76892 T + 0.0003025 T2

e = 0.01675104 – 0.0000418 T – 0.000000126 T2

M = 358.47583 + 35999.04975 T – 0.000150 T2 – 0.0000033 T3

p = L - M

EXAMPLE

Calculate the orbital elements of Earth on 1976-July-20, 12:00 UT.

T = 0.765503080

L = 99.69668 + (36000.76892 × 0.765503080) + (0.0003025 × 0.7655030802) = 27658.396351o = 298.396351o

Similarly,

a = 1.0000002 AU
e = 0.016718968
M = 195.859204o

p = 298.396351 – 195.859204 = 102.537147o

From the values of e and M, calculate the eccentric anomaly E from

E = M + e × sin E

The above is a transcendental equation in E that must be solved by iteration. It is also important that the angles M and E be expressed in radians. The equation can be solved in degrees if we replace e with eo = e × 180/p, giving us

E = M + eo × sin E

Then calculate the true anomaly n from

tan(n/2) = [ (1 + e) / (1 – e) ]1/2 × tan(E/2)

The radius vector of the planet can be calculated by one of the following

r = a × ( 1 – e × cos E )

r = a × ( 1 – e2 ) / ( 1 + e × cos n )

The planet's argument of latitude is

u = L + n – M – W

The ecliptical longitude l can be deduced from (lW), which is given by

tan(lW) = cos i × tan u

If i < 90o, as for the major planets, (lW) and u must lie in the same quadrant. When a programmable calculator or computer is used, in order to avoid the use of other tests, the preceding formula can be better written as

tan(lW) = cos i × sin u / cos u

and then the conversion from rectangular to polar coordinates should be applied to the numerator and the denominator of the fraction in the right-hand side. This will give (lW) directly in the correct quadrant.

The planet's ecliptic latitude b is given by

sin b = sin u × sin i

with –90o < b < 90o.

We have now obtained the heliocentric ecliptical coordinates l, b, and r of the planet for the given instant.

EXAMPLE

Calculate the heliocentric ecliptical coordinates of Mars on 1976-July-20, 12:00 UT.

As previously calculated, we have

L = 186.764387o
a = 1.5236883 AU
e = 0.093383330
i = 1.849824o
w = 286.250750o
W = 49.376635o
p = 335.627385o
M = 211.137002o

Following the steps outlined above, we have

eo = 0.093383330 × 180/p = 5.3504707

E = 211.137002 + 5.3504707 × sin E
E = 208.577611o

tan(n/2) = [ (1 + 0.093383330) / (1 – 0.093383330) ]1/2 × tan(208.577611 / 2)
n = 206.114239o

r = 1.5236883 × ( 1 – 0.093383330 × cos 208.577611 ) = 1.648641 AU

u = 186.764387 + 206.114239 – 211.137002 – 49.376635 = 132.364988o

tan(l – 49.376635) = cos 1.849824 × sin 132.364988 / cos 132.364988
l = 181.756494o

sin b = sin 132.364988 × sin 1.849824
b = 1.366666o

A perturbation is a disturbance in the regular and usually elliptical course of motion of a celestial body that is produced by some force additional to that that causes its regular motion. The most important perturbations in the motions of the major planets are caused by the gravitational influence of other planets. These perturbations must be accounted for if better accuracy is needed than that attainable using the above data alone. The perturbations in the motions of the giant planets are particularly important; in longitude, they can be larger than 0.3 degree for Jupiter, and larger than 1.0 degree for Saturn.

Since the purpose of this web site is to provide only a basic understanding of interplanetary space flight, the extra accuracy attained by taking into account the principle perturbations is unnecessary. The planetary positions derived by the methods described above are adequate for illustrative and education purposes.

In order to calculate an accurate position of the Moon, it is necessary to take into account hundreds of periodic terms in the Moon's longitude, latitude and parallax. In his book, Jean Meeus presents a rigorous method for calculating the Moon's position, however he limits his solution to about 100 of the most important periodic terms, being satisfied with a small inaccuracy. However, even this less accurate solution is too cumbersome to present here. Fortunately, Meeus provides suggestions to simplify the method to produce a low accuracy solution. The accuracy appears to be about ± 0.3o in longitude, ± 0.1o in latitude, and ± 0.01o in parallax.

Using the method described below, one obtains the geocentric longitude l and the geocentric latitude b of the center of the Moon, referred to the mean equinox of the date. If necessary, l and b can be converted to right ascension a and declination d using the following formulae.

tan a = ( sin l cos e – tan b sin e ) / cos l

sin d = sin b cos e + cos b sin e sin l

where e, the obliquity of the ecliptic, is
e = 23.452294 – 0.0130125 T – 0.00000164 T2 + 0.000000503 T3

The equatorial horizontal parallax p of the Moon too is obtained. When the parallax p is known, the distance between the centers of Earth and Moon, in kilometers, can be found from

D = 6378.14 / sin p

For the given instant, calculate the JD and then T by means of the formula previous presented. Then calculate the angles L', M, M', D and F by means of the following formulae, in which the various constants are expressed in degrees and decimals.

Moon's mean longitude:
L' = 270.434164 + 481267.8831 × T

Sun's mean anomaly:
M = 358.475833 + 35999.0498 × T

Moon's mean anomaly:
M' = 296.104608 + 477198.8491 × T

Moon's mean elongation:
D = 350.737486 + 445267.1142 × T

Mean distance of Moon from its ascending node:
F = 11.250889 + 483202.0251 × T

With the values of L', M, M', D and F calculated, l, b and p can be obtained by means of the following expressions where, again, all the coefficients are given in degrees and decimals.

l = L' + 6.288750 sin M'
+ 1.274018 sin(2D–M')
+ 0.658309 sin 2D
+ 0.213616 sin 2M'
– 0.185596 sin M
– 0.114336 sin 2F

b = 5.128189 sin F
+ 0.280606 sin(M'+F)
+ 0.277693 sin(M'–F)
+ 0.173238 sin(2D–F)
+ 0.055413 sin(2D+F–M')
+ 0.046272 sin(2D–F–M')

p = 0.950724
+ 0.051818 cos M'
+ 0.009531 cos(2D–M')
+ 0.007843 cos 2D
+ 0.002824 cos 2M'
+ 0.000857 cos(2D+M')

EXAMPLE

Calculate the geocentric coordinates and parallax of the Moon on 1968-December-24, 10:00 UT.

Using the previously demonstrated method,

JD = 2440214.9167,   T = 0.689799224

Using the above equations for angles L', M, M', D and F, we have

L' = 328.646595o
M = 350.592460o
M' = 67.500542o
D = 55.647457o
F = 323.632971o

The Moon's longitude is L' plus the sum of the periodic terms

l = 328.646595
+ 5.810070
+ 0.881713
+ 0.613362
+ 0.151046
+ 0.030337
+ 0.109184

l = 336.242307o

Similarly,

b = –2.480685o
p = 0.9717311o

The distance to the Moon is,

D = 6378.14 / sin 0.9717311 = 376,090 km

Converting geocentric ecliptical coordinates to right ascension and declination, we have

e = 23.443317o
a = 331.29323o = 22h5m10.4s
d = –14.41295o

The angular separation distance d between two celestial bodies, whose longitudes and latitudes are known, is given by the formula

cos d = sin b1 sin b2 + cos b1 cos b2 cos(l1l2)

EXAMPLE

In a previous example, we found that the coordinates of Mars on 1976-July-20 12:00 UT are

l = 181.756494o,   b = 1.366666o

We're given that on the same date and time the coordinates of Earth are

l = 297.883130o,   b = 0o

The angular separation between the planets is

cos d = sin 1.366666 × sin 0 + cos 1.366666 × cos 0 × cos(181.756494 – 297.883130)
d = 116.118642o