source: trunk/FACT++/pal/palUe2el.c

Last change on this file was 18347, checked in by tbretz, 9 years ago
File size: 9.0 KB
Line 
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
212void 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}
Note: See TracBrowser for help on using the repository browser.