home *** CD-ROM | disk | FTP | other *** search
/ Usenet 1994 January / usenetsourcesnewsgroupsinfomagicjanuary1994.iso / sources / x / volume19 / xephem / part21 / pelement.c < prev    next >
Encoding:
C/C++ Source or Header  |  1993-05-15  |  4.9 KB  |  161 lines
  1. #include <stdio.h>
  2. #include <math.h>
  3. #include "astro.h"
  4.  
  5. #if defined(__STDC__) || defined(__cplusplus)
  6. #define P_(s) s
  7. #else
  8. #define P_(s) ()
  9. #endif
  10.  
  11. extern void range P_((double *v, double r));
  12.  
  13. void pelement P_((double mjd, double plan[8][9]));
  14.  
  15. #undef P_
  16.  
  17.  
  18. /* this array contains polynomial coefficients to find the various orbital
  19.  *   elements for the mean orbit at any instant in time for each major planet.
  20.  *   the first five elements are in the form a0 + a1*t + a2*t**2 + a3*t**3,
  21.  *   where t is the number of Julian centuries of 36525 Julian days since 1900
  22.  *   Jan 0.5. the last three elements are constants.
  23.  *
  24.  * the orbital element (column) indeces are:
  25.  *   [ 0- 3]: coefficients for mean longitude, in degrees;
  26.  *   [ 4- 7]: coefficients for longitude of the perihelion, in degrees;
  27.  *   [ 8-11]: coefficients for eccentricity;
  28.  *   [12-15]: coefficients for inclination, in degrees;
  29.  *   [16-19]: coefficients for longitude of the ascending node, in degrees;
  30.  *      [20]: semi-major axis, in AU;
  31.  *      [21]: angular diameter at 1 AU, in arcsec;
  32.  *      [22]: standard visual magnitude, ie, the visual magnitude of the planet
  33.  *          when at a distance of 1 AU from both the Sun and the Earth and
  34.  *          with zero phase angle.
  35.  *
  36.  * the planent (row) indeces are:
  37.  *   [0]: Mercury; [1]: Venus;   [2]: Mars;  [3]: Jupiter; [4]: Saturn;
  38.  *   [5]: Uranus;  [6]: Neptune; [7]: Pluto.
  39.  */
  40. #define    NPELE    (5*4 + 3)    /* 4 coeffs for ea of 5 elems, + 3 constants */
  41. static double elements[8][NPELE] = {
  42.  
  43.     {   /*     mercury... */
  44.  
  45.         178.179078,    415.2057519,    3.011e-4,    0.0,
  46.         75.899697,    1.5554889,    2.947e-4,    0.0,
  47.         .20561421,    2.046e-5,    3e-8,        0.0,
  48.         7.002881,    1.8608e-3,    -1.83e-5,    0.0,
  49.         47.145944,    1.1852083,    1.739e-4,    0.0,
  50.         .3870986,    6.74,         -0.42
  51.     },
  52.  
  53.     {   /*     venus... */
  54.  
  55.         342.767053,    162.5533664,    3.097e-4,    0.0,
  56.         130.163833,    1.4080361,    -9.764e-4,    0.0,
  57.         6.82069e-3,    -4.774e-5,    9.1e-8,        0.0,
  58.         3.393631,    1.0058e-3,    -1e-6,        0.0,
  59.         75.779647,    .89985,        4.1e-4,        0.0,
  60.         .7233316,    16.92,        -4.4
  61.     },
  62.  
  63.     {   /*     mars... */
  64.  
  65.         293.737334,    53.17137642,    3.107e-4,    0.0,
  66.         3.34218203e2, 1.8407584,    1.299e-4,    -1.19e-6,
  67.         9.33129e-2,    9.2064e-5,    7.7e-8,        0.0,
  68.         1.850333,    -6.75e-4,    1.26e-5,    0.0,
  69.         48.786442,    .7709917,    -1.4e-6,    -5.33e-6,
  70.         1.5236883,    9.36,        -1.52
  71.     },
  72.  
  73.     {   /*     jupiter... */
  74.  
  75.         238.049257,    8.434172183,    3.347e-4,    -1.65e-6,
  76.         1.2720972e1, 1.6099617,    1.05627e-3,    -3.43e-6,
  77.         4.833475e-2, 1.6418e-4,    -4.676e-7,    -1.7e-9,
  78.         1.308736,    -5.6961e-3,    3.9e-6,        0.0,
  79.         99.443414,    1.01053,    3.5222e-4,    -8.51e-6,
  80.         5.202561,    196.74,        -9.4
  81.     },
  82.  
  83.     {   /*     saturn... */
  84.  
  85.         266.564377,    3.398638567,    3.245e-4,    -5.8e-6,
  86.         9.1098214e1, 1.9584158,    8.2636e-4,    4.61e-6,
  87.         5.589232e-2, -3.455e-4,    -7.28e-7,    7.4e-10,
  88.         2.492519,    -3.9189e-3,    -1.549e-5,    4e-8,
  89.         112.790414,    .8731951,    -1.5218e-4,    -5.31e-6,
  90.         9.554747,    165.6,        -8.88
  91.     },
  92.  
  93.     {   /*     uranus... */
  94.  
  95.         244.19747,    1.194065406,    3.16e-4,    -6e-7,
  96.         1.71548692e2, 1.4844328,    2.372e-4,    -6.1e-7,
  97.         4.63444e-2,    -2.658e-5,    7.7e-8,        0.0,
  98.         .772464,    6.253e-4,    3.95e-5,    0.0,
  99.         73.477111,    .4986678,    1.3117e-3,    0.0,
  100.         19.21814,    65.8,        -7.19
  101.     },
  102.  
  103.     {   /*     neptune... */
  104.  
  105.         84.457994,    .6107942056,    3.205e-4,    -6e-7,
  106.         4.6727364e1, 1.4245744,    3.9082e-4,    -6.05e-7,
  107.         8.99704e-3,    6.33e-6,    -2e-9,        0.0,
  108.         1.779242,    -9.5436e-3,    -9.1e-6,    0.0,
  109.         130.681389,    1.098935,    2.4987e-4,    -4.718e-6,
  110.         30.10957,    62.2,        -6.87
  111.     },
  112.  
  113.     {   /*     pluto...(osculating 1984 jan 21) */
  114.  
  115.         95.3113544,    .3980332167,    0.0,        0.0,
  116.         224.017,    0.0,        0.0,        0.0,
  117.         .25515,    0.0,        0.0,        0.0,
  118.         17.1329,    0.0,        0.0,        0.0,
  119.         110.191,    0.0,        0.0,        0.0,
  120.         39.8151,    8.2,        -1.0
  121.     }
  122. };
  123.  
  124. /* given a modified Julian date, mjd, return the elements for the mean orbit
  125.  *   at that instant of all the major planets, together with their
  126.  *   mean daily motions in longitude, angular diameter and standard visual
  127.  *   magnitude.
  128.  * plan[i][j] contains all the values for all the planets at mjd, such that
  129.  *   i = 0..7: mercury, venus, mars, jupiter, saturn, unranus, neptune, pluto;
  130.  *   j = 0..8: mean longitude, mean daily motion in longitude, longitude of 
  131.  *     the perihelion, eccentricity, inclination, longitude of the ascending
  132.  *     node, length of the semi-major axis, angular diameter from 1 AU, and
  133.  *     the standard visual magnitude (see elements[][] comment, above).
  134.  */
  135. void
  136. pelement (mjd, plan)
  137. double mjd;
  138. double plan[8][9];
  139. {
  140.     register double *ep, *pp;
  141.     register double t = mjd/36525.;
  142.     double aa;
  143.     int planet, i;
  144.  
  145.     for (planet = 0; planet < 8; planet++) {
  146.         ep = elements[planet];
  147.         pp = plan[planet];
  148.         aa = ep[1]*t;
  149.         pp[0] = ep[0] + 360.*(aa-(long)aa) + (ep[3]*t + ep[2])*t*t;
  150.         range (pp, 360.);
  151.         pp[1] = (ep[1]*9.856263e-3) + (ep[2] + ep[3])/36525;
  152.  
  153.         for (i = 4; i < 20; i += 4)
  154.         pp[i/4+1] = ((ep[i+3]*t + ep[i+2])*t + ep[i+1])*t + ep[i+0];
  155.  
  156.         pp[6] = ep[20];
  157.         pp[7] = ep[21];
  158.         pp[8] = ep[22];
  159.     }
  160. }
  161.