line |
stmt |
bran |
cond |
sub |
pod |
time |
code |
1
|
|
|
|
|
|
|
/* |
2
|
|
|
|
|
|
|
*+ |
3
|
|
|
|
|
|
|
* Name: |
4
|
|
|
|
|
|
|
* palPlanel |
5
|
|
|
|
|
|
|
|
6
|
|
|
|
|
|
|
* Purpose: |
7
|
|
|
|
|
|
|
* Transform conventional elements into position and velocity |
8
|
|
|
|
|
|
|
|
9
|
|
|
|
|
|
|
* Language: |
10
|
|
|
|
|
|
|
* Starlink ANSI C |
11
|
|
|
|
|
|
|
|
12
|
|
|
|
|
|
|
* Type of Module: |
13
|
|
|
|
|
|
|
* Library routine |
14
|
|
|
|
|
|
|
|
15
|
|
|
|
|
|
|
* Invocation: |
16
|
|
|
|
|
|
|
* void palPlanel ( double date, int jform, double epoch, double orbinc, |
17
|
|
|
|
|
|
|
* double anode, double perih, double aorq, double e, |
18
|
|
|
|
|
|
|
* double aorl, double dm, double pv[6], int *jstat ); |
19
|
|
|
|
|
|
|
|
20
|
|
|
|
|
|
|
* Arguments: |
21
|
|
|
|
|
|
|
* date = double (Given) |
22
|
|
|
|
|
|
|
* Epoch (TT MJD) of osculation (Note 1) |
23
|
|
|
|
|
|
|
* jform = int (Given) |
24
|
|
|
|
|
|
|
* Element set actually returned (1-3; Note 3) |
25
|
|
|
|
|
|
|
* epoch = double (Given) |
26
|
|
|
|
|
|
|
* Epoch of elements (TT MJD) (Note 4) |
27
|
|
|
|
|
|
|
* orbinc = double (Given) |
28
|
|
|
|
|
|
|
* inclination (radians) |
29
|
|
|
|
|
|
|
* anode = double (Given) |
30
|
|
|
|
|
|
|
* longitude of the ascending node (radians) |
31
|
|
|
|
|
|
|
* perih = double (Given) |
32
|
|
|
|
|
|
|
* longitude or argument of perihelion (radians) |
33
|
|
|
|
|
|
|
* aorq = double (Given) |
34
|
|
|
|
|
|
|
* mean distance or perihelion distance (AU) |
35
|
|
|
|
|
|
|
* e = double (Given) |
36
|
|
|
|
|
|
|
* eccentricity |
37
|
|
|
|
|
|
|
* aorl = double (Given) |
38
|
|
|
|
|
|
|
* mean anomaly or longitude (radians, JFORM=1,2 only) |
39
|
|
|
|
|
|
|
* dm = double (Given) |
40
|
|
|
|
|
|
|
* daily motion (radians, JFORM=1 only) |
41
|
|
|
|
|
|
|
* u = double [13] (Returned) |
42
|
|
|
|
|
|
|
* Universal orbital elements (Note 1) |
43
|
|
|
|
|
|
|
* (0) combined mass (M+m) |
44
|
|
|
|
|
|
|
* (1) total energy of the orbit (alpha) |
45
|
|
|
|
|
|
|
* (2) reference (osculating) epoch (t0) |
46
|
|
|
|
|
|
|
* (3-5) position at reference epoch (r0) |
47
|
|
|
|
|
|
|
* (6-8) velocity at reference epoch (v0) |
48
|
|
|
|
|
|
|
* (9) heliocentric distance at reference epoch |
49
|
|
|
|
|
|
|
* (10) r0.v0 |
50
|
|
|
|
|
|
|
* (11) date (t) |
51
|
|
|
|
|
|
|
* (12) universal eccentric anomaly (psi) of date, approx |
52
|
|
|
|
|
|
|
* jstat = int * (Returned) |
53
|
|
|
|
|
|
|
* status: 0 = OK |
54
|
|
|
|
|
|
|
* - -1 = illegal JFORM |
55
|
|
|
|
|
|
|
* - -2 = illegal E |
56
|
|
|
|
|
|
|
* - -3 = illegal AORQ |
57
|
|
|
|
|
|
|
* - -4 = illegal DM |
58
|
|
|
|
|
|
|
* - -5 = numerical error |
59
|
|
|
|
|
|
|
|
60
|
|
|
|
|
|
|
* Description: |
61
|
|
|
|
|
|
|
* Heliocentric position and velocity of a planet, asteroid or comet, |
62
|
|
|
|
|
|
|
* starting from orbital elements. |
63
|
|
|
|
|
|
|
|
64
|
|
|
|
|
|
|
* Authors: |
65
|
|
|
|
|
|
|
* PTW: Pat Wallace (STFC) |
66
|
|
|
|
|
|
|
* TIMJ: Tim Jenness (JAC, Hawaii) |
67
|
|
|
|
|
|
|
* {enter_new_authors_here} |
68
|
|
|
|
|
|
|
|
69
|
|
|
|
|
|
|
* Notes: |
70
|
|
|
|
|
|
|
* - DATE is the instant for which the prediction is required. It is |
71
|
|
|
|
|
|
|
* in the TT timescale (formerly Ephemeris Time, ET) and is a |
72
|
|
|
|
|
|
|
* Modified Julian Date (JD-2400000.5). |
73
|
|
|
|
|
|
|
* - The elements are with respect to the J2000 ecliptic and equinox. |
74
|
|
|
|
|
|
|
* - A choice of three different element-set options is available: |
75
|
|
|
|
|
|
|
* |
76
|
|
|
|
|
|
|
* Option JFORM = 1, suitable for the major planets: |
77
|
|
|
|
|
|
|
* |
78
|
|
|
|
|
|
|
* EPOCH = epoch of elements (TT MJD) |
79
|
|
|
|
|
|
|
* ORBINC = inclination i (radians) |
80
|
|
|
|
|
|
|
* ANODE = longitude of the ascending node, big omega (radians) |
81
|
|
|
|
|
|
|
* PERIH = longitude of perihelion, curly pi (radians) |
82
|
|
|
|
|
|
|
* AORQ = mean distance, a (AU) |
83
|
|
|
|
|
|
|
* E = eccentricity, e (range 0 to <1) |
84
|
|
|
|
|
|
|
* AORL = mean longitude L (radians) |
85
|
|
|
|
|
|
|
* DM = daily motion (radians) |
86
|
|
|
|
|
|
|
* |
87
|
|
|
|
|
|
|
* Option JFORM = 2, suitable for minor 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 = argument of perihelion, little omega (radians) |
93
|
|
|
|
|
|
|
* AORQ = mean distance, a (AU) |
94
|
|
|
|
|
|
|
* E = eccentricity, e (range 0 to <1) |
95
|
|
|
|
|
|
|
* AORL = mean anomaly M (radians) |
96
|
|
|
|
|
|
|
* |
97
|
|
|
|
|
|
|
* Option JFORM = 3, suitable for comets: |
98
|
|
|
|
|
|
|
* |
99
|
|
|
|
|
|
|
* EPOCH = epoch of elements and perihelion (TT MJD) |
100
|
|
|
|
|
|
|
* ORBINC = inclination i (radians) |
101
|
|
|
|
|
|
|
* ANODE = longitude of the ascending node, big omega (radians) |
102
|
|
|
|
|
|
|
* PERIH = argument of perihelion, little omega (radians) |
103
|
|
|
|
|
|
|
* AORQ = perihelion distance, q (AU) |
104
|
|
|
|
|
|
|
* E = eccentricity, e (range 0 to 10) |
105
|
|
|
|
|
|
|
* |
106
|
|
|
|
|
|
|
* Unused arguments (DM for JFORM=2, AORL and DM for JFORM=3) are not |
107
|
|
|
|
|
|
|
* accessed. |
108
|
|
|
|
|
|
|
* - Each of the three element sets defines an unperturbed heliocentric |
109
|
|
|
|
|
|
|
* orbit. For a given epoch of observation, the position of the body |
110
|
|
|
|
|
|
|
* in its orbit can be predicted from these elements, which are |
111
|
|
|
|
|
|
|
* called "osculating elements", using standard two-body analytical |
112
|
|
|
|
|
|
|
* solutions. However, due to planetary perturbations, a given set |
113
|
|
|
|
|
|
|
* of osculating elements remains usable for only as long as the |
114
|
|
|
|
|
|
|
* unperturbed orbit that it describes is an adequate approximation |
115
|
|
|
|
|
|
|
* to reality. Attached to such a set of elements is a date called |
116
|
|
|
|
|
|
|
* the "osculating epoch", at which the elements are, momentarily, |
117
|
|
|
|
|
|
|
* a perfect representation of the instantaneous position and |
118
|
|
|
|
|
|
|
* velocity of the body. |
119
|
|
|
|
|
|
|
* |
120
|
|
|
|
|
|
|
* Therefore, for any given problem there are up to three different |
121
|
|
|
|
|
|
|
* epochs in play, and it is vital to distinguish clearly between |
122
|
|
|
|
|
|
|
* them: |
123
|
|
|
|
|
|
|
* |
124
|
|
|
|
|
|
|
* . The epoch of observation: the moment in time for which the |
125
|
|
|
|
|
|
|
* position of the body is to be predicted. |
126
|
|
|
|
|
|
|
* |
127
|
|
|
|
|
|
|
* . The epoch defining the position of the body: the moment in time |
128
|
|
|
|
|
|
|
* at which, in the absence of purturbations, the specified |
129
|
|
|
|
|
|
|
* position (mean longitude, mean anomaly, or perihelion) is |
130
|
|
|
|
|
|
|
* reached. |
131
|
|
|
|
|
|
|
* |
132
|
|
|
|
|
|
|
* . The osculating epoch: the moment in time at which the given |
133
|
|
|
|
|
|
|
* elements are correct. |
134
|
|
|
|
|
|
|
* |
135
|
|
|
|
|
|
|
* For the major-planet and minor-planet cases it is usual to make |
136
|
|
|
|
|
|
|
* the epoch that defines the position of the body the same as the |
137
|
|
|
|
|
|
|
* epoch of osculation. Thus, only two different epochs are |
138
|
|
|
|
|
|
|
* involved: the epoch of the elements and the epoch of observation. |
139
|
|
|
|
|
|
|
* |
140
|
|
|
|
|
|
|
* For comets, the epoch of perihelion fixes the position in the |
141
|
|
|
|
|
|
|
* orbit and in general a different epoch of osculation will be |
142
|
|
|
|
|
|
|
* chosen. Thus, all three types of epoch are involved. |
143
|
|
|
|
|
|
|
* |
144
|
|
|
|
|
|
|
* For the present routine: |
145
|
|
|
|
|
|
|
* |
146
|
|
|
|
|
|
|
* . The epoch of observation is the argument DATE. |
147
|
|
|
|
|
|
|
* |
148
|
|
|
|
|
|
|
* . The epoch defining the position of the body is the argument |
149
|
|
|
|
|
|
|
* EPOCH. |
150
|
|
|
|
|
|
|
* |
151
|
|
|
|
|
|
|
* . The osculating epoch is not used and is assumed to be close |
152
|
|
|
|
|
|
|
* enough to the epoch of observation to deliver adequate accuracy. |
153
|
|
|
|
|
|
|
* If not, a preliminary call to palPertel may be used to update |
154
|
|
|
|
|
|
|
* the element-set (and its associated osculating epoch) by |
155
|
|
|
|
|
|
|
* applying planetary perturbations. |
156
|
|
|
|
|
|
|
* - The reference frame for the result is with respect to the mean |
157
|
|
|
|
|
|
|
* equator and equinox of epoch J2000. |
158
|
|
|
|
|
|
|
* - The algorithm was originally adapted from the EPHSLA program of |
159
|
|
|
|
|
|
|
* D.H.P.Jones (private communication, 1996). The method is based |
160
|
|
|
|
|
|
|
* on Stumpff's Universal Variables. |
161
|
|
|
|
|
|
|
|
162
|
|
|
|
|
|
|
* See Also: |
163
|
|
|
|
|
|
|
* Everhart, E. & Pitkin, E.T., Am.J.Phys. 51, 712, 1983. |
164
|
|
|
|
|
|
|
|
165
|
|
|
|
|
|
|
* History: |
166
|
|
|
|
|
|
|
* 2012-03-12 (TIMJ): |
167
|
|
|
|
|
|
|
* Initial version taken directly from SLA/F. |
168
|
|
|
|
|
|
|
* Adapted with permission from the Fortran SLALIB library. |
169
|
|
|
|
|
|
|
* {enter_further_changes_here} |
170
|
|
|
|
|
|
|
|
171
|
|
|
|
|
|
|
* Copyright: |
172
|
|
|
|
|
|
|
* Copyright (C) 2002 Rutherford Appleton Laboratory |
173
|
|
|
|
|
|
|
* Copyright (C) 2012 Science and Technology Facilities Council. |
174
|
|
|
|
|
|
|
* All Rights Reserved. |
175
|
|
|
|
|
|
|
|
176
|
|
|
|
|
|
|
* Licence: |
177
|
|
|
|
|
|
|
* This program is free software; you can redistribute it and/or |
178
|
|
|
|
|
|
|
* modify it under the terms of the GNU General Public License as |
179
|
|
|
|
|
|
|
* published by the Free Software Foundation; either version 3 of |
180
|
|
|
|
|
|
|
* the License, or (at your option) any later version. |
181
|
|
|
|
|
|
|
* |
182
|
|
|
|
|
|
|
* This program is distributed in the hope that it will be |
183
|
|
|
|
|
|
|
* useful, but WITHOUT ANY WARRANTY; without even the implied |
184
|
|
|
|
|
|
|
* warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR |
185
|
|
|
|
|
|
|
* PURPOSE. See the GNU General Public License for more details. |
186
|
|
|
|
|
|
|
* |
187
|
|
|
|
|
|
|
* You should have received a copy of the GNU General Public License |
188
|
|
|
|
|
|
|
* along with this program; if not, write to the Free Software |
189
|
|
|
|
|
|
|
* Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, |
190
|
|
|
|
|
|
|
* MA 02110-1301, USA. |
191
|
|
|
|
|
|
|
|
192
|
|
|
|
|
|
|
* Bugs: |
193
|
|
|
|
|
|
|
* {note_any_bugs_here} |
194
|
|
|
|
|
|
|
*- |
195
|
|
|
|
|
|
|
*/ |
196
|
|
|
|
|
|
|
|
197
|
|
|
|
|
|
|
#include "pal.h" |
198
|
|
|
|
|
|
|
|
199
|
|
|
|
|
|
|
|
200
|
0
|
|
|
|
|
|
void palPlanel ( double date, int jform, double epoch, double orbinc, |
201
|
|
|
|
|
|
|
double anode, double perih, double aorq, double e, |
202
|
|
|
|
|
|
|
double aorl, double dm, double pv[6], int *jstat ) { |
203
|
|
|
|
|
|
|
|
204
|
|
|
|
|
|
|
int j; |
205
|
|
|
|
|
|
|
double u[13]; |
206
|
|
|
|
|
|
|
|
207
|
|
|
|
|
|
|
/* Validate elements and convert to "universal variables" parameters. */ |
208
|
0
|
|
|
|
|
|
palEl2ue( date, jform, epoch, orbinc, anode, perih, aorq, e, aorl, |
209
|
|
|
|
|
|
|
dm, u, &j ); |
210
|
|
|
|
|
|
|
|
211
|
|
|
|
|
|
|
/* Determine the position and velocity */ |
212
|
0
|
0
|
|
|
|
|
if (j == 0) { |
213
|
0
|
|
|
|
|
|
palUe2pv( date, u, pv, &j); |
214
|
0
|
0
|
|
|
|
|
if (j != 0) j = -5; |
215
|
|
|
|
|
|
|
} |
216
|
|
|
|
|
|
|
|
217
|
|
|
|
|
|
|
/* Wrap up */ |
218
|
0
|
|
|
|
|
|
*jstat = j; |
219
|
|
|
|
|
|
|
|
220
|
0
|
|
|
|
|
|
} |