VSOP87 Theory Equations Summary

VSOP87 Theory Equations Summary Part 2
by Jay Tanner

Spherical Coordinates (L,B,R)

VSOP87 series D is used to compute the heliocentric spherical LBR-coordinates or series C is used to compute the equivalent heliocentric rectangular XYZ-coordinates (see Part 1).

In the following VSOP87 spherical equations, we will encounter the following variables:

${JDTT}$ = Julian Day number corresponding to a given date and time (TT)
${t}$ = Julian millennia, corresponding to JDTT and reckoned from J2000.0

When we want to compute the position of a planet for a given date and time, we first compute the JDTT value and then compute the corresponding (${t}$) value from it. This will be the actual numerical value used as the time variable in the equations.
$$t = \frac{JDTT - 2451545}{365250}$$ The ${A,B,C}$ values in each series term in the equations are taken from the original FORTRAN program data tables. There may be thousands of terms in a sub-series, each with its own set of ${A,B,C}$ values. These values are used by the code generator tool to generate the numerical terms in the output source code.

${A}$ = Amplitude
${B}$ = Phase (not to be confused with the spherical B-coordinate)
${C}$ = Frequency

Indices used in the summation equations:
${n}$ = Sub-series order from 0 to 5
${j}$ = Index of term within sub-series (from 1 to ${k}$)
${k}$ = Total number of terms in sub-series (can be from 1 to multiple thousands)

Let:
${L}$, ${B}$, ${R}$ = Heliocentric Spherical Coordinates in radians and distance in AU

and
${L}$ = Heliocentric ecliptical longitude in radians
${B}$ = Heliocentric ecliptical latitude in radians (not to be confused with phase B)
${R}$ = Heliocentric radius vector (distance) in AU

These spherical coordinates are computed using the VSOP87 sub-series function modules generated by the Multi-Language VSOP87 Source Code Generator Tool.

1 AU (Astronomical Unit) = 149597870.691 km = 92955807.267 mi

Speed of Light = 299792.458 km/s = 186282.397 mi/s

Computing the Planetary Spherical Heliocentric Ecliptical LBR-Coordinates

In the following case of spherical coordinates, the heliocentric longitude and latitude (L, B) will be expressed in radians and the heliocentric distance (R) in astronomical units.

VSOP87 Series B = Spherical heliocentric ecliptical LBR-coordinates - Equinox J2000.0
VSOP87 Series D = Spherical heliocentric ecliptical LBR-coordinates - Equinox of Date

The VSOP87 B and D series correspond to the A and C series. The only difference is the type of coordinates used and the equinox to which they are referred. The method of computation is essentially the same.

These spherical coordinate equations are analogous to the equations used to compute the rectangular coordinates in part 1, with the computational roles of (X, Y, Z) being replaced by (L, B, R) respectively.

For each order (n), the VSOP87 (L_n B_n R_n) sub-series summations may be expressed as

$$L = \sum_{i=0}^5 ~~ l_{i} ~~ = {l_{0}+l_{1}+l_{2}+l_{3}+l_{4}+l_{5}}$$ $$B = \sum_{i=0}^5 ~~ b_{i} ~~ = {b_{0}+b_{1}+b_{2}+b_{3}+b_{4}+b_{5}}$$ $$R = \sum_{i=0}^5 ~~ r_{i} ~~ = {r_{0}+r_{1}+r_{2}+r_{3}+r_{4}+r_{5}}$$

There are six orders of LBR sub-series summations (${n}$=0 to ${n}$=5). The above sub-series summations are computed separately for each order (${n}$=0 to ${n}$=5) for each of the ${k}$ terms of the given sub-series. The essential source code to compute the numerical values of these coordinates for any given moment over spans of thousands of years can be generated by the VSOP87 Source Code Generator tool.

The six separate partial summations returned by the modules for the respective sub-series are in turn summed to obtain the final, complete ecliptical spherical LBR-coordinates of the planet.

For each order (${n}$), from 0 to 5, the VSOP87 ($l_{n} ~ b_{n} ~ r_{n}$) sub-series summations
may be expressed as:

$$l_{n} = \left(\sum_{j=1}^k ~ A_{L,n,j} \cdot {cos(B_{L,n,j} ~ + ~ t \cdot C_{L,n,j})}\right) \cdot t^n$$ $$b_{n} = \left(\sum_{j=1}^k ~ A_{B,n,j} \cdot {cos(B_{B,n,j} ~ + ~ t \cdot C_{B,n,j})}\right) \cdot t^n$$ $$r_{n} = \left(\sum_{j=1}^k ~ A_{R,n,j} \cdot {cos(B_{R,n,j} ~ + ~ t \cdot C_{R,n,j})}\right) \cdot t^n$$

The complete VSOP87 solutions in terms of spherical variables can be expressed by
the following double summations.

$$L = \sum_{n=0}^5\left(\left(\sum_{j=1}^k ~ A_{L,n,j} \cdot {cos(B_{L,n,j} ~ + ~ t \cdot C_{L,n,j})}\right) \cdot t^n\right)$$ $$B = \sum_{n=0}^5\left(\left(\sum_{j=1}^k ~ A_{B,n,j} \cdot {cos(B_{B,n,j} ~ + ~ t \cdot C_{B,n,j})}\right) \cdot t^n\right)$$ $$R = \sum_{n=0}^5\left(\left(\sum_{j=1}^k ~ A_{R,n,j} \cdot {cos(B_{R,n,j} ~ + ~ t \cdot C_{R,n,j})}\right) \cdot t^n\right)$$

Example For The Planet Mars

As an example, for the planet Mars, the complete system of general equations could be written in the following PHP code with each of the functions representing code generated by the VSOP87 Source Code Generator tool, designed to evaluate the ordered sub-series summations.

After using the tool to generate the source code for the ordered sub-functions, the programmer then only needs to write the following simple PHP code to call and sum the six returned sub-series values to obtain the respective heliocentric LBR coordinates.


 $t = ($JDTT - 2451545) / 365250;

 $L = Mars_B_L0($t) + Mars_B_L1($t) + Mars_B_L2($t)
    + Mars_B_L3($t) + Mars_B_L4($t) + Mars_B_L5($t);

 $B = Mars_B_B0($t) + Mars_B_B1($t) + Mars_B_B2($t)
    + Mars_B_B3($t) + Mars_B_B4($t) + Mars_B_B5($t);

 $R = Mars_B_R0($t) + Mars_B_R1($t) + Mars_B_R2($t)
    + Mars_B_R3($t) + Mars_B_R4($t) + Mars_B_R5($t);

These are the raw, instantaneous, geometric, heliocentric coordinates from which the apparent geocentric ecliptical and equatorial coordinates and distance are subsequently derived. By raw geometric coordinates we mean that no corrections for the location of the observer have yet been applied, such as for light-time, aberration, horizontal parallax, etc. This is the actual location of the planet in space, not its apparent location to the eye, which will require the aforementioned corrections, primarily due to the finite speed of light.

The heliocentric spherical longitude (L) and latitude (B) will be returned in radians which will need to be converted to degrees. The longitude value will then need to be modulated to fall in the range from 0 to 360 degrees. The heliocentric distance or radius vector (R) will be returned in AUs.

Computing the heliocentric coordinates for both the earth and the planet at the same moment is the first and most laborious step in computing the VSOP87 position of a planet. Once that is accomplished, it becomes a relatively trivial matter to convert heliocentric coordinates into geocentric or other coordinates.

Jay Tanner - 2022
COBOL - Language of Space Aliens