| 1 | #include "slalib.h" | 
|---|
| 2 | #include "slamac.h" | 
|---|
| 3 | void slaPertel (int jform, double date0, double date1, | 
|---|
| 4 | double epoch0, double orbi0, double anode0, | 
|---|
| 5 | double perih0, double aorq0, double e0, double am0, | 
|---|
| 6 | double *epoch1, double *orbi1, double *anode1, | 
|---|
| 7 | double *perih1, double *aorq1, double *e1, double *am1, | 
|---|
| 8 | int *jstat ) | 
|---|
| 9 | /* | 
|---|
| 10 | **  - - - - - - - - - - | 
|---|
| 11 | **   s l a P e r t e l | 
|---|
| 12 | **  - - - - - - - - - - | 
|---|
| 13 | ** | 
|---|
| 14 | **  Update the osculating orbital elements of an asteroid or comet by | 
|---|
| 15 | **  applying planetary perturbations. | 
|---|
| 16 | ** | 
|---|
| 17 | **  Given (format and dates): | 
|---|
| 18 | **     jform   int      choice of element set (2 or 3; Note 1) | 
|---|
| 19 | **     date0   double   date of osculation (TT MJD) for the given elements | 
|---|
| 20 | **     date1   double   date of osculation (TT MJD) for the updated elements | 
|---|
| 21 | ** | 
|---|
| 22 | **  Given (the unperturbed elements): | 
|---|
| 23 | **     epoch0  double   epoch (TT MJD) of the given element set (Note 2) | 
|---|
| 24 | **     orbi0   double   inclination (radians) | 
|---|
| 25 | **     anode0  double   longitude of the ascending node (radians) | 
|---|
| 26 | **     perih0  double   argument of perihelion (radians) | 
|---|
| 27 | **     aorq0   double   mean distance or perihelion distance (AU) | 
|---|
| 28 | **     e0      double   eccentricity | 
|---|
| 29 | **     am0     double   mean anomaly (radians, jform=2 only) | 
|---|
| 30 | ** | 
|---|
| 31 | **  Returned (the updated elements): | 
|---|
| 32 | **     epoch1  double*  epoch (TT MJD) of the updated element set (Note 2) | 
|---|
| 33 | **     orbi1   double*  inclination (radians) | 
|---|
| 34 | **     anode1  double*  longitude of the ascending node (radians) | 
|---|
| 35 | **     perih1  double*  argument of perihelion (radians) | 
|---|
| 36 | **     aorq1   double*  mean distance or perihelion distance (AU) | 
|---|
| 37 | **     e1      double*  eccentricity | 
|---|
| 38 | **     am1     double*  mean anomaly (radians, jform=2 only) | 
|---|
| 39 | ** | 
|---|
| 40 | **  Returned (status flag): | 
|---|
| 41 | **     jstat   int*     status: +102 = warning, distant epoch | 
|---|
| 42 | **                              +101 = warning, large timespan ( > 100 years) | 
|---|
| 43 | **                          +1 to +8 = coincident with major planet (Note 6) | 
|---|
| 44 | **                                 0 = OK | 
|---|
| 45 | **                                -1 = illegal jform | 
|---|
| 46 | **                                -2 = illegal e0 | 
|---|
| 47 | **                                -3 = illegal aorq0 | 
|---|
| 48 | **                                -4 = internal error | 
|---|
| 49 | **                                -5 = numerical error | 
|---|
| 50 | ** | 
|---|
| 51 | **  Notes: | 
|---|
| 52 | ** | 
|---|
| 53 | **  1  Two different element-format options are available: | 
|---|
| 54 | ** | 
|---|
| 55 | **     Option jform=2, suitable for minor planets: | 
|---|
| 56 | ** | 
|---|
| 57 | **     epoch   = epoch of elements (TT MJD) | 
|---|
| 58 | **     orbi    = inclination i (radians) | 
|---|
| 59 | **     anode   = longitude of the ascending node, big omega (radians) | 
|---|
| 60 | **     perih   = argument of perihelion, little omega (radians) | 
|---|
| 61 | **     aorq    = mean distance, a (AU) | 
|---|
| 62 | **     e       = eccentricity, e | 
|---|
| 63 | **     am      = mean anomaly M (radians) | 
|---|
| 64 | ** | 
|---|
| 65 | **     Option jform=3, suitable for comets: | 
|---|
| 66 | ** | 
|---|
| 67 | **     epoch   = epoch of perihelion (TT MJD) | 
|---|
| 68 | **     orbi    = inclination i (radians) | 
|---|
| 69 | **     anode   = longitude of the ascending node, big omega (radians) | 
|---|
| 70 | **     perih   = argument of perihelion, little omega (radians) | 
|---|
| 71 | **     aorq    = perihelion distance, q (AU) | 
|---|
| 72 | **     e       = eccentricity, e | 
|---|
| 73 | ** | 
|---|
| 74 | **  2  date0, date1, epoch0 and epoch1 are all instants of time in | 
|---|
| 75 | **     the TT timescale (formerly Ephemeris Time, ET), expressed | 
|---|
| 76 | **     as Modified Julian Dates (JD-2400000.5). | 
|---|
| 77 | ** | 
|---|
| 78 | **     date0 is the instant at which the given (i.e. unperturbed) | 
|---|
| 79 | **     osculating elements are correct. | 
|---|
| 80 | ** | 
|---|
| 81 | **     date1 is the specified instant at which the updated osculating | 
|---|
| 82 | **     elements are correct. | 
|---|
| 83 | ** | 
|---|
| 84 | **     epoch0 and epoch1 will be the same as date0 and date1 | 
|---|
| 85 | **     (respectively) for the jform=2 case, normally used for minor | 
|---|
| 86 | **     planets.  For the jform=3 case, the two epochs will refer to | 
|---|
| 87 | **     perihelion passage and so will not, in general, be the same as | 
|---|
| 88 | **     date0 and/or date1 though they may be similar to one another. | 
|---|
| 89 | ** | 
|---|
| 90 | **  3  The elements are with respect to the J2000 ecliptic and equinox. | 
|---|
| 91 | ** | 
|---|
| 92 | **  4  Unused elements (am0 and am1 for jform=3) are not accessed. | 
|---|
| 93 | ** | 
|---|
| 94 | **  5  See the slaPertue routine for details of the algorithm used. | 
|---|
| 95 | ** | 
|---|
| 96 | **  6  This routine is not intended to be used for major planets, which | 
|---|
| 97 | **     is why jform=1 is not available and why there is no opportunity | 
|---|
| 98 | **     to specify either the longitude of perihelion or the daily | 
|---|
| 99 | **     motion.  However, if jform=2 elements are somehow obtained for a | 
|---|
| 100 | **     major planet and supplied to the routine, sensible results will, | 
|---|
| 101 | **     in fact, be produced.  This happens because the slaPertue  routine | 
|---|
| 102 | **     that is called to perform the calculations checks the separation | 
|---|
| 103 | **     between the body and each of the planets and interprets a | 
|---|
| 104 | **     suspiciously small value (1E-3 AU) as an attempt to apply it to | 
|---|
| 105 | **     the planet concerned.  If this condition is detected, the | 
|---|
| 106 | **     contribution from that planet is ignored, and the status is set to | 
|---|
| 107 | **     the planet number (Mercury=1,...,Neptune=8) as a warning. | 
|---|
| 108 | ** | 
|---|
| 109 | **  Reference: | 
|---|
| 110 | ** | 
|---|
| 111 | **     Sterne, Theodore E., "An Introduction to Celestial Mechanics", | 
|---|
| 112 | **     Interscience Publishers Inc., 1960.  Section 6.7, p199. | 
|---|
| 113 | ** | 
|---|
| 114 | **  Called:  slaEl2ue,  slaPertue,  slaUe2el | 
|---|
| 115 | ** | 
|---|
| 116 | **  Last revision:   14 March 1999 | 
|---|
| 117 | ** | 
|---|
| 118 | **  Copyright 1999 P.T.Wallace.  All rights reserved. | 
|---|
| 119 | */ | 
|---|
| 120 | { | 
|---|
| 121 | double u[13], dm; | 
|---|
| 122 | int j, jf; | 
|---|
| 123 |  | 
|---|
| 124 |  | 
|---|
| 125 | /* Check that the elements are either minor-planet or comet format. */ | 
|---|
| 126 | if ( jform < 2 || jform > 3 ) { | 
|---|
| 127 | *jstat = -1; | 
|---|
| 128 | return; | 
|---|
| 129 | } | 
|---|
| 130 |  | 
|---|
| 131 | /* Transform the elements from conventional to universal form. */ | 
|---|
| 132 | slaEl2ue ( date0, jform, epoch0, orbi0, anode0, perih0, | 
|---|
| 133 | aorq0, e0, am0, 0.0, u, &j ); | 
|---|
| 134 | if ( j ) { | 
|---|
| 135 | *jstat = j; | 
|---|
| 136 | return; | 
|---|
| 137 | } | 
|---|
| 138 |  | 
|---|
| 139 | /* Update the universal elements. */ | 
|---|
| 140 | slaPertue ( date1, u, &j ); | 
|---|
| 141 | if ( j > 0 ) { | 
|---|
| 142 | *jstat = j; | 
|---|
| 143 | } else if ( j < 0 ) { | 
|---|
| 144 | *jstat = -5; | 
|---|
| 145 | return; | 
|---|
| 146 | } | 
|---|
| 147 |  | 
|---|
| 148 | /* Transform from universal to conventional elements. */ | 
|---|
| 149 | slaUe2el ( u, 2, | 
|---|
| 150 | &jf, epoch1, orbi1, anode1, perih1, aorq1, e1, am1, &dm, &j ); | 
|---|
| 151 | if ( jf != jform || j ) *jstat = -5; | 
|---|
| 152 |  | 
|---|
| 153 | } | 
|---|