| 1 | /*
|
|---|
| 2 | *+
|
|---|
| 3 | * Name:
|
|---|
| 4 | * palUe2el
|
|---|
| 5 |
|
|---|
| 6 | * Purpose:
|
|---|
| 7 | * Universal elements to heliocentric osculating elements
|
|---|
| 8 |
|
|---|
| 9 | * Language:
|
|---|
| 10 | * Starlink ANSI C
|
|---|
| 11 |
|
|---|
| 12 | * Type of Module:
|
|---|
| 13 | * Library routine
|
|---|
| 14 |
|
|---|
| 15 | * Invocation:
|
|---|
| 16 | * void palUe2el ( const double u[13], int jformr,
|
|---|
| 17 | * int *jform, double *epoch, double *orbinc,
|
|---|
| 18 | * double *anode, double *perih, double *aorq, double *e,
|
|---|
| 19 | * double *aorl, double *dm, int *jstat );
|
|---|
| 20 |
|
|---|
| 21 | * Arguments:
|
|---|
| 22 | * u = const double [13] (Given)
|
|---|
| 23 | * Universal orbital elements (Note 1)
|
|---|
| 24 | * (0) combined mass (M+m)
|
|---|
| 25 | * (1) total energy of the orbit (alpha)
|
|---|
| 26 | * (2) reference (osculating) epoch (t0)
|
|---|
| 27 | * (3-5) position at reference epoch (r0)
|
|---|
| 28 | * (6-8) velocity at reference epoch (v0)
|
|---|
| 29 | * (9) heliocentric distance at reference epoch
|
|---|
| 30 | * (10) r0.v0
|
|---|
| 31 | * (11) date (t)
|
|---|
| 32 | * (12) universal eccentric anomaly (psi) of date, approx
|
|---|
| 33 | * jformr = int (Given)
|
|---|
| 34 | * Requested element set (1-3; Note 3)
|
|---|
| 35 | * jform = int * (Returned)
|
|---|
| 36 | * Element set actually returned (1-3; Note 4)
|
|---|
| 37 | * epoch = double * (Returned)
|
|---|
| 38 | * Epoch of elements (TT MJD)
|
|---|
| 39 | * orbinc = double * (Returned)
|
|---|
| 40 | * inclination (radians)
|
|---|
| 41 | * anode = double * (Returned)
|
|---|
| 42 | * longitude of the ascending node (radians)
|
|---|
| 43 | * perih = double * (Returned)
|
|---|
| 44 | * longitude or argument of perihelion (radians)
|
|---|
| 45 | * aorq = double * (Returned)
|
|---|
| 46 | * mean distance or perihelion distance (AU)
|
|---|
| 47 | * e = double * (Returned)
|
|---|
| 48 | * eccentricity
|
|---|
| 49 | * aorl = double * (Returned)
|
|---|
| 50 | * mean anomaly or longitude (radians, JFORM=1,2 only)
|
|---|
| 51 | * dm = double * (Returned)
|
|---|
| 52 | * daily motion (radians, JFORM=1 only)
|
|---|
| 53 | * jstat = int * (Returned)
|
|---|
| 54 | * status: 0 = OK
|
|---|
| 55 | * -1 = illegal combined mass
|
|---|
| 56 | * -2 = illegal JFORMR
|
|---|
| 57 | * -3 = position/velocity out of range
|
|---|
| 58 |
|
|---|
| 59 | * Description:
|
|---|
| 60 | * Transform universal elements into conventional heliocentric
|
|---|
| 61 | * osculating elements.
|
|---|
| 62 |
|
|---|
| 63 | * Authors:
|
|---|
| 64 | * PTW: Patrick T. Wallace
|
|---|
| 65 | * TIMJ: Tim Jenness (JAC, Hawaii)
|
|---|
| 66 | * {enter_new_authors_here}
|
|---|
| 67 |
|
|---|
| 68 | * Notes:
|
|---|
| 69 | * - The "universal" elements are those which define the orbit for the
|
|---|
| 70 | * purposes of the method of universal variables (see reference 2).
|
|---|
| 71 | * They consist of the combined mass of the two bodies, an epoch,
|
|---|
| 72 | * and the position and velocity vectors (arbitrary reference frame)
|
|---|
| 73 | * at that epoch. The parameter set used here includes also various
|
|---|
| 74 | * quantities that can, in fact, be derived from the other
|
|---|
| 75 | * information. This approach is taken to avoiding unnecessary
|
|---|
| 76 | * computation and loss of accuracy. The supplementary quantities
|
|---|
| 77 | * are (i) alpha, which is proportional to the total energy of the
|
|---|
| 78 | * orbit, (ii) the heliocentric distance at epoch, (iii) the
|
|---|
| 79 | * outwards component of the velocity at the given epoch, (iv) an
|
|---|
| 80 | * estimate of psi, the "universal eccentric anomaly" at a given
|
|---|
| 81 | * date and (v) that date.
|
|---|
| 82 | * - The universal elements are with respect to the mean equator and
|
|---|
| 83 | * equinox of epoch J2000. The orbital elements produced are with
|
|---|
| 84 | * respect to the J2000 ecliptic and mean equinox.
|
|---|
| 85 | * - Three different element-format options are supported:
|
|---|
| 86 | *
|
|---|
| 87 | * Option JFORM=1, suitable for the major planets:
|
|---|
| 88 | *
|
|---|
| 89 | * EPOCH = epoch of elements (TT MJD)
|
|---|
| 90 | * ORBINC = inclination i (radians)
|
|---|
| 91 | * ANODE = longitude of the ascending node, big omega (radians)
|
|---|
| 92 | * PERIH = longitude of perihelion, curly pi (radians)
|
|---|
| 93 | * AORQ = mean distance, a (AU)
|
|---|
| 94 | * E = eccentricity, e
|
|---|
| 95 | * AORL = mean longitude L (radians)
|
|---|
| 96 | * DM = daily motion (radians)
|
|---|
| 97 | *
|
|---|
| 98 | * Option JFORM=2, suitable for minor planets:
|
|---|
| 99 | *
|
|---|
| 100 | * EPOCH = epoch of elements (TT MJD)
|
|---|
| 101 | * ORBINC = inclination i (radians)
|
|---|
| 102 | * ANODE = longitude of the ascending node, big omega (radians)
|
|---|
| 103 | * PERIH = argument of perihelion, little omega (radians)
|
|---|
| 104 | * AORQ = mean distance, a (AU)
|
|---|
| 105 | * E = eccentricity, e
|
|---|
| 106 | * AORL = mean anomaly M (radians)
|
|---|
| 107 | *
|
|---|
| 108 | * Option JFORM=3, suitable for comets:
|
|---|
| 109 | *
|
|---|
| 110 | * EPOCH = epoch of perihelion (TT MJD)
|
|---|
| 111 | * ORBINC = inclination i (radians)
|
|---|
| 112 | * ANODE = longitude of the ascending node, big omega (radians)
|
|---|
| 113 | * PERIH = argument of perihelion, little omega (radians)
|
|---|
| 114 | * AORQ = perihelion distance, q (AU)
|
|---|
| 115 | * E = eccentricity, e
|
|---|
| 116 | *
|
|---|
| 117 | * - It may not be possible to generate elements in the form
|
|---|
| 118 | * requested through JFORMR. The caller is notified of the form
|
|---|
| 119 | * of elements actually returned by means of the JFORM argument:
|
|---|
| 120 | *
|
|---|
| 121 | * JFORMR JFORM meaning
|
|---|
| 122 | *
|
|---|
| 123 | * 1 1 OK - elements are in the requested format
|
|---|
| 124 | * 1 2 never happens
|
|---|
| 125 | * 1 3 orbit not elliptical
|
|---|
| 126 | *
|
|---|
| 127 | * 2 1 never happens
|
|---|
| 128 | * 2 2 OK - elements are in the requested format
|
|---|
| 129 | * 2 3 orbit not elliptical
|
|---|
| 130 | *
|
|---|
| 131 | * 3 1 never happens
|
|---|
| 132 | * 3 2 never happens
|
|---|
| 133 | * 3 3 OK - elements are in the requested format
|
|---|
| 134 | *
|
|---|
| 135 | * - The arguments returned for each value of JFORM (cf Note 6: JFORM
|
|---|
| 136 | * may not be the same as JFORMR) are as follows:
|
|---|
| 137 | *
|
|---|
| 138 | * JFORM 1 2 3
|
|---|
| 139 | * EPOCH t0 t0 T
|
|---|
| 140 | * ORBINC i i i
|
|---|
| 141 | * ANODE Omega Omega Omega
|
|---|
| 142 | * PERIH curly pi omega omega
|
|---|
| 143 | * AORQ a a q
|
|---|
| 144 | * E e e e
|
|---|
| 145 | * AORL L M -
|
|---|
| 146 | * DM n - -
|
|---|
| 147 | *
|
|---|
| 148 | * where:
|
|---|
| 149 | *
|
|---|
| 150 | * t0 is the epoch of the elements (MJD, TT)
|
|---|
| 151 | * T " epoch of perihelion (MJD, TT)
|
|---|
| 152 | * i " inclination (radians)
|
|---|
| 153 | * Omega " longitude of the ascending node (radians)
|
|---|
| 154 | * curly pi " longitude of perihelion (radians)
|
|---|
| 155 | * omega " argument of perihelion (radians)
|
|---|
| 156 | * a " mean distance (AU)
|
|---|
| 157 | * q " perihelion distance (AU)
|
|---|
| 158 | * e " eccentricity
|
|---|
| 159 | * L " longitude (radians, 0-2pi)
|
|---|
| 160 | * M " mean anomaly (radians, 0-2pi)
|
|---|
| 161 | * n " daily motion (radians)
|
|---|
| 162 | * - means no value is set
|
|---|
| 163 | *
|
|---|
| 164 | * - At very small inclinations, the longitude of the ascending node
|
|---|
| 165 | * ANODE becomes indeterminate and under some circumstances may be
|
|---|
| 166 | * set arbitrarily to zero. Similarly, if the orbit is close to
|
|---|
| 167 | * circular, the true anomaly becomes indeterminate and under some
|
|---|
| 168 | * circumstances may be set arbitrarily to zero. In such cases,
|
|---|
| 169 | * the other elements are automatically adjusted to compensate,
|
|---|
| 170 | * and so the elements remain a valid description of the orbit.
|
|---|
| 171 |
|
|---|
| 172 | * See Also:
|
|---|
| 173 | * - Sterne, Theodore E., "An Introduction to Celestial Mechanics",
|
|---|
| 174 | * Interscience Publishers Inc., 1960. Section 6.7, p199.
|
|---|
| 175 | * - Everhart, E. & Pitkin, E.T., Am.J.Phys. 51, 712, 1983.
|
|---|
| 176 |
|
|---|
| 177 | * History:
|
|---|
| 178 | * 2012-03-09 (TIMJ):
|
|---|
| 179 | * Initial version
|
|---|
| 180 | * Adapted with permission from the Fortran SLALIB library.
|
|---|
| 181 | * {enter_further_changes_here}
|
|---|
| 182 |
|
|---|
| 183 | * Copyright:
|
|---|
| 184 | * Copyright (C) 1999 Rutherford Appleton Laboratory
|
|---|
| 185 | * Copyright (C) 2012 Science and Technology Facilities Council.
|
|---|
| 186 | * All Rights Reserved.
|
|---|
| 187 |
|
|---|
| 188 | * Licence:
|
|---|
| 189 | * This program is free software; you can redistribute it and/or
|
|---|
| 190 | * modify it under the terms of the GNU General Public License as
|
|---|
| 191 | * published by the Free Software Foundation; either version 3 of
|
|---|
| 192 | * the License, or (at your option) any later version.
|
|---|
| 193 | *
|
|---|
| 194 | * This program is distributed in the hope that it will be
|
|---|
| 195 | * useful, but WITHOUT ANY WARRANTY; without even the implied
|
|---|
| 196 | * warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
|
|---|
| 197 | * PURPOSE. See the GNU General Public License for more details.
|
|---|
| 198 | *
|
|---|
| 199 | * You should have received a copy of the GNU General Public License
|
|---|
| 200 | * along with this program; if not, write to the Free Software
|
|---|
| 201 | * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston,
|
|---|
| 202 | * MA 02110-1301, USA.
|
|---|
| 203 |
|
|---|
| 204 | * Bugs:
|
|---|
| 205 | * {note_any_bugs_here}
|
|---|
| 206 | *-
|
|---|
| 207 | */
|
|---|
| 208 |
|
|---|
| 209 | #include "pal.h"
|
|---|
| 210 | #include "palmac.h"
|
|---|
| 211 |
|
|---|
| 212 | void palUe2el ( const double u[], int jformr,
|
|---|
| 213 | int *jform, double *epoch, double *orbinc,
|
|---|
| 214 | double *anode, double *perih, double *aorq, double *e,
|
|---|
| 215 | double *aorl, double *dm, int *jstat ) {
|
|---|
| 216 |
|
|---|
| 217 | /* Canonical days to seconds */
|
|---|
| 218 | const double CD2S = PAL__GCON / PAL__SPD;
|
|---|
| 219 |
|
|---|
| 220 | int i;
|
|---|
| 221 | double pmass, date, pv[6];
|
|---|
| 222 |
|
|---|
| 223 | /* Unpack the universal elements */
|
|---|
| 224 | pmass = u[0] - 1.0;
|
|---|
| 225 | date = u[2];
|
|---|
| 226 | for (i=0; i<3; i++) {
|
|---|
| 227 | pv[i] = u[i+3];
|
|---|
| 228 | pv[i+3] = u[i+6] * CD2S;
|
|---|
| 229 | }
|
|---|
| 230 |
|
|---|
| 231 | /* Convert the position and velocity etc into conventional elements */
|
|---|
| 232 | palPv2el( pv, date, pmass, jformr, jform, epoch, orbinc, anode,
|
|---|
| 233 | perih, aorq, e, aorl, dm, jstat );
|
|---|
| 234 | }
|
|---|
