The Low Precision Ephemeris (LPE)

Table of Contents

Overview

The Low Precision Ephemeris (LPE) is a dynamic library for the computation of solar, lunar and planetary positions within a precision of about one minute of arc (1'). Compared to state of the art ephemerides like that of the Jet Propulsion Laboratory (JPL), this is a quite moderate precision. On the other hand, the low precision comes with a high performance since the library is written purely in x86 code using Microsoft's MASM32 assembler in version 6.14.8444. On an average notebook with an Intel Celeron CPU clocked with 2.5 GHz, the computation of a geocentric planetary position costs about 2.4 microseconds which means a performance of about 420.000 geocentric planetary positions per second.

The LPE can be useful in the following areas:

The library is assembled as a DLL which can be included in high-level programming languages like C and Java. For Java, no special wrapper code is necessary, the mapping to Java data is already included in the DLL. See below for details on how to access the LPE from other programs. The DLL file itself is smaller than 10 KBytes. Additionally, ephemeris files with a size of 5 KB per century will be loaded into memory.

Idea of compression

I tried to find out how Paul Ahnert's classical tables Astronomisch-chronologische Tafeln für Sonne, Mond und Planeten (Barth, Leipzig 1969) are constructed. They are based on a "mean" elliptical orbit, most orbital elements considered as fixed or changing linearly in time. The mean anomaly is computed with tables for the change in days, months, years and centuries. By adding the corresponding terms and reducing on the circle, one obtains the value for the mean anomaly to be used.

How can this work for an arbitrary century? The planetary perturbations can give errors of several degrees - how can they be taken into respect with such an approach? Ahnert's idea obviously was to put a secular perturbation information into the century's values of his table of mean anomalies, yielding sufficiently good results for many practical purposes. Compared to some brute force approaches using Cebysheff polynomials or the like for a piecewise approximation of the orbit, Ahnert's concept seemed more natural to me, since it builds on the elliptical orbit of the planet and only modifies it slightly by "correcting" the mean anomaly in such a way that the computed position will fit better.

The following figure shows the method I therefore pursued: I projected the true position of a planet (Pw in the figure) on the plane of its mean orbit of date. The projection point has a certain anomaly which can be converted to a corresponding mean anomaly by stepping through the Keplerian procedure to determine the true anomaly in backward direction. This mean anomaly will differ from the mean anomaly that can be computed using the mean orbital elements. These difference terms which I call ΔM values can be computed and saved once in binary data files and then used as corrections in the LPE algorithm.

Similarly, there is a correction term ΔR to obtain the true distance from the Sun from the approximating point on the mean orbit. For each outer planet (Jupiter to Pluto), the values of ΔM and ΔR will be tabulated in 100 day steps. Different to ΔM, the ΔR value will be applied as a factor to the distance after the solution of the Kepler problem.

Organisation of the data

As mentioned, the LPE uses binary data files for the computation of planetary positions. These data files are computed per century, with a size of about 5 KB each, containing linearly and/or quadratically regressed mean orbital elements and additionally, for the outer planets from Jupiter on, in 100-day intervals, correction terms for mean anomaly and radius vector. The correction terms are stored in bit fields of variable size, depending on the minimum and maximum size of the correction terms in the century. The intercepts, slopes and quadratic coefficients for the mean orbital elements are stored as REAL8 in the usual IEE748 format. For compatibility with a Java class working with the same data, the first release of LPE adapts the data format formerly used in this Java class (see http://www.astrotexte.ch/sources/java/doc/AstroOL/doc-files/LowPrecCalculator.html. This includes some settings which are inconvenient for an assembly ephemeris — like the Big Endian byte order in which the real8 numbers are stored. In future releases, an optimized format for the LPE will be used.

The following assembler definitions document the organization of the first part of the data files (for developers, code will serve better as a documentation than any other text).

; Coefficients for quadratic interpolation
qCoeff struc
  c0     real8 0.  ; intercept
  c1     real8 0.  ; slope
  c2     real8 0.  ; quadratic coefficient
qCoeff ends

; Coefficients for linear interpolation
lCoeff struc
  c0     real8 0.  ; intercept
  c1     real8 0.  ; slope
lCoeff ends

; Orbital elements for linear interpolation
; Each element contains two coefficients - intercept and slope
lOrbEl struc
  a    lCoeff<>  ; semi-axis
  e    lCoeff<>  ; eccentricity
  i    lCoeff<>  ; inclination
  node lCoeff<>  ; ascending node
  peri lCoeff<>  ; perihelion
  M    lCoeff<>  ; mean anomaly
  n    lCoeff<>  ; daily motion
lOrbEl ends

; Orbital elements for quadratic interpolation
; Each elements contains the three coefficents of the parabola
; used for quadratic interpolation
qOrbEl struc
  a    qCoeff<>  ; semi-axis        
  e    qCoeff<>  ; eccentricity     
  i    qCoeff<>  ; inclination      
  node qCoeff<>  ; ascending node   
  peri qCoeff<>  ; perihelion       
  M    qCoeff<>  ; mean anomaly     
  n    qCoeff<>  ; daily motion     
qOrbEl ends

; Pair of 2 Byte Integers: Minimum and Maximum
wordMiniMax struc
  min dw ?
  max dw ?
wordMiniMax ends

; Minima and Maxima of r and m
rmMiniMax struc
  m wordMiniMax<>
  r wordMiniMax<>
rmMiniMax ends

; Minima und Maxima of the r and m values contained in the file
MiniMax struc
  jup    rmMiniMax<>
  sat    rmMiniMax<>
  ura    rmMiniMax<>
  nep    rmMiniMax<>
  plu    rmMiniMax<>
MiniMax ends

; Structure of the file
fileStruc struc
  jd0      db 3 dup(0)   ; 3-byte int: Julian start date
  jd1      db 3 dup(0)   ; 3-byte int: Julian end date
  incr     db 2 dup(0)   ; 2-byte int: increment in days
  mer      lOrbEl<>
  ven      lOrbEl<>
  mar      lOrbEl<>
  jup      qOrbEl<>
  sat      qOrbEl<>
  ura      qOrbEl<>
  nep      qOrbEl<>
  plu      qOrbEl<>
  minimax  MiniMax<>
  data     db 0 ; From here on, the bit field for delta M / delta R values start
fileStruc ends

Precision

The following diagram shows the precision of LPE for the various planets, the data of the latest Swiss Ephemeris release (1.75/2008) serving as reference. Beginning with January 1, 1950, the positions have been calculated for each of the 36525 days of the century and compared to the Swiss Ephemeris positions. The graphic shows the distribution of the deviations in minutes of arc for each planet – for the current epoch 1950-2050 All the planets are in scope, as well as Sun and Moon - only Mars gives a deviation of up to 5' of arc.

These figures being in scope for the current century, the differences deteriorate quickly when leaving this time range, as the following table shows. It gives the maximum deviations for each planet per century file. The deviations increase when going to the past - up to values of about 1 degree in 50 C.E. There seems to be a systematic secular error which should be analyzed in a future release. One source of the error clearly is the fixed orbit for the earth, copied from Montenbruck's Grundlagen der Ephemeridenrechnung, but this contributes with only 6' in 50 C.E. Also, the analytical moon theory used seems to be optimized for the current times.
Maximum deviations (in ' of arc) of computed positions
Start year Sun Moon Mercury Venus Mars Jupiter Saturn Uranus Neptune Pluto
20501.0511.8441.4182.7035.7830.7430.6870.6750.6871.130
19501.1071.7251.4692.6144.3280.7410.7180.6960.6671.418
18501.1121.6601.3922.8003.8400.7330.7150.6770.6651.288
17501.1631.9281.3932.6743.1700.7290.6840.7220.7030.936
16501.0031.8851.4002.4133.2290.7270.6860.6860.6610.915
15500.9513.0491.1881.8672.8190.6880.6880.6960.6630.836
14500.8323.8621.5001.6583.6220.6770.6390.6750.6651.783
13500.9325.8291.8262.4745.0570.6550.6430.6640.6461.714
12501.1337.8921.9343.2764.9970.5850.6150.6390.6440.944
11501.4139.9242.5694.2596.1350.5330.5920.6220.6371.004
10501.65513.0092.7615.2576.3620.5250.5880.6300.6430.785
9501.99415.7083.4426.4307.2660.4720.5390.6110.6101.952
8502.25019.0084.0197.2948.7220.4240.5120.5900.6311.932
7502.61923.1434.8208.4668.4290.5610.4920.5870.6140.932
6503.02226.5005.34110.1938.5940.6960.4560.5340.5910.950
5503.45230.8825.87011.5939.7750.7250.4210.5450.5880.689
4503.89435.7096.98613.99610.8590.8670.4090.5270.5501.837
3504.55840.3337.98515.15611.8290.9610.4320.4910.5470.660
2505.11846.1618.82617.86913.8991.1570.5110.4760.5491.016
1505.78751.3939.87719.46714.8611.3140.5830.4180.5160.825
506.42957.31411.48922.42816.1091.4880.6610.4350.4900.663

Usage in C (or other languages)

Usage of LPE in C – or in any other language supporting "C" calling conventions for external function calls – is quite simple. Dependent on your IDE, you will have to include the DLL in one or the other way. In your code, you declare the following three functions - this is the full public interface of the DLL:

int lpeGetPlanet( unsigned int pPlanet, double pJD, double* pLBR);
int lpeGetPlanets( double pJD, double* pLBR);
int lpeSetDirectory( char* dir );  

Usage in Java

LPE contains a built-in public interface to Java, using the Java Native Interface (JNI) on the Java side. Once implemented, the class name which binds the library is fixed. In the LPE case, the Java class which calls the library must have the name LpeLoader and must be located in the root of the package tree (i.e. no package name). If you want to use a different class name, you must either build an additional JNI wrapper DLL at your own, or you have to adapt the public interface of the LPE assembler code and rebuild it.

The following code contains the minimum code to make the LPE accessible for Java calls. The functions can then be called as static methods like LpeLoader.lpeGetPlanets() and so on.

class LpeLoader {
  
// at class load, locate and load the LPE library
  static {
    System.loadLibrary("lpe");
    }
    
    
// Available functions from the library
  public static native int lpeGetPlanet(int ipl, double jd, double[] LPR);
  public static native int lpeGetPlanets(double jd, double[] coor);
  public static native int lpeSetDirectory(String dir);
    
}

Further development

The LPE is made available to the public with a CVS repository. You are welcome to improve the code - may it be for your own version or contributions for forthcoming releases.

There is a GNU makefile containing the build instructions for lpe.dll and a test program lpe_test.exe which will be assembled and executed automatically with each build of lpe.dll.

The output of lpe_test.exe accords to the very simple Test Anything Protocol (TAP). If everything went fine, it looks like this:

1..23
ok 1
ok 2
[ and so on: ok 3, ok 4, ...]
ok 23

For future releases, the following items will be considered:

Download

The download section contains two zip files, binaries.zip and complete.zip. binaries.zip contains only the DLL and the data files. That's all that is necessary to use the ephemeris from other programs. Instead, complete.zip is an extended version of all the ressources used during development, including some C and Java sample programs. Also the "Ephemeris Compiler" is included, a C program for producing the data files.

For download and for all information around this project, visit the Project Page.

The development is under the GNU public license.

Winterthur, January 2009
Rüdiger Plantiko