line |
stmt |
bran |
cond |
sub |
pod |
time |
code |
1
|
|
|
|
|
|
|
/* |
2
|
|
|
|
|
|
|
*+ |
3
|
|
|
|
|
|
|
* Name: |
4
|
|
|
|
|
|
|
* palRdplan |
5
|
|
|
|
|
|
|
|
6
|
|
|
|
|
|
|
* Purpose: |
7
|
|
|
|
|
|
|
* Approximate topocentric apparent RA,Dec of a planet |
8
|
|
|
|
|
|
|
|
9
|
|
|
|
|
|
|
* Language: |
10
|
|
|
|
|
|
|
* Starlink ANSI C |
11
|
|
|
|
|
|
|
|
12
|
|
|
|
|
|
|
* Type of Module: |
13
|
|
|
|
|
|
|
* Library routine |
14
|
|
|
|
|
|
|
|
15
|
|
|
|
|
|
|
* Invocation: |
16
|
|
|
|
|
|
|
* void palRdplan( double date, int np, double elong, double phi, |
17
|
|
|
|
|
|
|
* double * ra, double * dec, double * diam ); |
18
|
|
|
|
|
|
|
|
19
|
|
|
|
|
|
|
* Arguments: |
20
|
|
|
|
|
|
|
* date = double (Given) |
21
|
|
|
|
|
|
|
* MJD of observation (JD-2400000.5) in TDB. For all practical |
22
|
|
|
|
|
|
|
* purposes TT can be used instead of TDB, and for many applications |
23
|
|
|
|
|
|
|
* UT will do (except for the Moon). |
24
|
|
|
|
|
|
|
* np = int (Given) |
25
|
|
|
|
|
|
|
* Planet: 1 = Mercury |
26
|
|
|
|
|
|
|
* 2 = Venus |
27
|
|
|
|
|
|
|
* 3 = Moon |
28
|
|
|
|
|
|
|
* 4 = Mars |
29
|
|
|
|
|
|
|
* 5 = Jupiter |
30
|
|
|
|
|
|
|
* 6 = Saturn |
31
|
|
|
|
|
|
|
* 7 = Uranus |
32
|
|
|
|
|
|
|
* 8 = Neptune |
33
|
|
|
|
|
|
|
* else = Sun |
34
|
|
|
|
|
|
|
* elong = double (Given) |
35
|
|
|
|
|
|
|
* Observer's east longitude (radians) |
36
|
|
|
|
|
|
|
* phi = double (Given) |
37
|
|
|
|
|
|
|
* Observer's geodetic latitude (radians) |
38
|
|
|
|
|
|
|
* ra = double * (Returned) |
39
|
|
|
|
|
|
|
* RA (topocentric apparent, radians) |
40
|
|
|
|
|
|
|
* dec = double * (Returned) |
41
|
|
|
|
|
|
|
* Dec (topocentric apparent, radians) |
42
|
|
|
|
|
|
|
* diam = double * (Returned) |
43
|
|
|
|
|
|
|
* Angular diameter (equatorial, radians) |
44
|
|
|
|
|
|
|
|
45
|
|
|
|
|
|
|
* Description: |
46
|
|
|
|
|
|
|
* Approximate topocentric apparent RA,Dec of a planet, and its |
47
|
|
|
|
|
|
|
* angular diameter. |
48
|
|
|
|
|
|
|
|
49
|
|
|
|
|
|
|
* Authors: |
50
|
|
|
|
|
|
|
* PTW: Patrick T. Wallace |
51
|
|
|
|
|
|
|
* TIMJ: Tim Jenness (JAC, Hawaii) |
52
|
|
|
|
|
|
|
* {enter_new_authors_here} |
53
|
|
|
|
|
|
|
|
54
|
|
|
|
|
|
|
* Notes: |
55
|
|
|
|
|
|
|
* - Unlike with slaRdplan, Pluto is not supported. |
56
|
|
|
|
|
|
|
* - The longitude and latitude allow correction for geocentric |
57
|
|
|
|
|
|
|
* parallax. This is a major effect for the Moon, but in the |
58
|
|
|
|
|
|
|
* context of the limited accuracy of the present routine its |
59
|
|
|
|
|
|
|
* effect on planetary positions is small (negligible for the |
60
|
|
|
|
|
|
|
* outer planets). Geocentric positions can be generated by |
61
|
|
|
|
|
|
|
* appropriate use of the routines palDmoon and eraPlan94. |
62
|
|
|
|
|
|
|
|
63
|
|
|
|
|
|
|
* History: |
64
|
|
|
|
|
|
|
* 2012-03-07 (TIMJ): |
65
|
|
|
|
|
|
|
* Initial version, with some documentation from SLA/F. |
66
|
|
|
|
|
|
|
* Adapted with permission from the Fortran SLALIB library. |
67
|
|
|
|
|
|
|
* {enter_further_changes_here} |
68
|
|
|
|
|
|
|
|
69
|
|
|
|
|
|
|
* Copyright: |
70
|
|
|
|
|
|
|
* Copyright (C) 1997 Rutherford Appleton Laboratory |
71
|
|
|
|
|
|
|
* Copyright (C) 2012 Science and Technology Facilities Council. |
72
|
|
|
|
|
|
|
* All Rights Reserved. |
73
|
|
|
|
|
|
|
|
74
|
|
|
|
|
|
|
* Licence: |
75
|
|
|
|
|
|
|
* This program is free software; you can redistribute it and/or |
76
|
|
|
|
|
|
|
* modify it under the terms of the GNU General Public License as |
77
|
|
|
|
|
|
|
* published by the Free Software Foundation; either version 3 of |
78
|
|
|
|
|
|
|
* the License, or (at your option) any later version. |
79
|
|
|
|
|
|
|
* |
80
|
|
|
|
|
|
|
* This program is distributed in the hope that it will be |
81
|
|
|
|
|
|
|
* useful, but WITHOUT ANY WARRANTY; without even the implied |
82
|
|
|
|
|
|
|
* warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR |
83
|
|
|
|
|
|
|
* PURPOSE. See the GNU General Public License for more details. |
84
|
|
|
|
|
|
|
* |
85
|
|
|
|
|
|
|
* You should have received a copy of the GNU General Public License |
86
|
|
|
|
|
|
|
* along with this program; if not, write to the Free Software |
87
|
|
|
|
|
|
|
* Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, |
88
|
|
|
|
|
|
|
* MA 02110-1301, USA. |
89
|
|
|
|
|
|
|
|
90
|
|
|
|
|
|
|
* Bugs: |
91
|
|
|
|
|
|
|
* {note_any_bugs_here} |
92
|
|
|
|
|
|
|
*- |
93
|
|
|
|
|
|
|
*/ |
94
|
|
|
|
|
|
|
|
95
|
|
|
|
|
|
|
#include |
96
|
|
|
|
|
|
|
|
97
|
|
|
|
|
|
|
#include "pal.h" |
98
|
|
|
|
|
|
|
#include "palmac.h" |
99
|
|
|
|
|
|
|
#include "pal1sofa.h" |
100
|
|
|
|
|
|
|
|
101
|
9
|
|
|
|
|
|
void palRdplan( double date, int np, double elong, double phi, |
102
|
|
|
|
|
|
|
double * ra, double * dec, double * diam ) { |
103
|
|
|
|
|
|
|
|
104
|
|
|
|
|
|
|
/* AU in km */ |
105
|
|
|
|
|
|
|
const double AUKM = 1.49597870e8; |
106
|
|
|
|
|
|
|
|
107
|
|
|
|
|
|
|
/* Equatorial radii (km) */ |
108
|
9
|
|
|
|
|
|
const double EQRAU[] = { |
109
|
|
|
|
|
|
|
696000.0, /* Sun */ |
110
|
|
|
|
|
|
|
2439.7, |
111
|
|
|
|
|
|
|
6051.9, |
112
|
|
|
|
|
|
|
1738, |
113
|
|
|
|
|
|
|
3397, |
114
|
|
|
|
|
|
|
71492, |
115
|
|
|
|
|
|
|
60268, |
116
|
|
|
|
|
|
|
25559, |
117
|
|
|
|
|
|
|
24764 |
118
|
|
|
|
|
|
|
}; |
119
|
|
|
|
|
|
|
|
120
|
|
|
|
|
|
|
/* Local variables */ |
121
|
|
|
|
|
|
|
int i, j; |
122
|
|
|
|
|
|
|
double stl; |
123
|
|
|
|
|
|
|
double vgm[6]; |
124
|
|
|
|
|
|
|
double v[6]; |
125
|
|
|
|
|
|
|
double rmat[3][3]; |
126
|
|
|
|
|
|
|
double vse[6]; |
127
|
|
|
|
|
|
|
double vsg[6]; |
128
|
|
|
|
|
|
|
double vsp[6]; |
129
|
|
|
|
|
|
|
double vgo[6]; |
130
|
|
|
|
|
|
|
double dx,dy,dz,r,tl; |
131
|
|
|
|
|
|
|
|
132
|
|
|
|
|
|
|
/* Classify np */ |
133
|
9
|
50
|
|
|
|
|
if (np < 0 || np > 8 ) np=0; /* Sun */ |
134
|
|
|
|
|
|
|
|
135
|
|
|
|
|
|
|
/* Approximate local sidereal time */ |
136
|
9
|
|
|
|
|
|
stl = palGmst( date - palDt( palEpj(date)) / 86400.0) + elong; |
137
|
|
|
|
|
|
|
|
138
|
|
|
|
|
|
|
/* Geocentre to Moon (mean of date) */ |
139
|
9
|
|
|
|
|
|
palDmoon( date, v ); |
140
|
|
|
|
|
|
|
|
141
|
|
|
|
|
|
|
/* Nutation to true of date */ |
142
|
9
|
|
|
|
|
|
palNut( date, rmat ); |
143
|
9
|
|
|
|
|
|
eraRxp( rmat, v, vgm ); |
144
|
9
|
|
|
|
|
|
eraRxp( rmat, &(v[3]), &(vgm[3]) ); |
145
|
|
|
|
|
|
|
|
146
|
|
|
|
|
|
|
/* Moon? */ |
147
|
9
|
100
|
|
|
|
|
if (np == 3) { |
148
|
|
|
|
|
|
|
|
149
|
|
|
|
|
|
|
/* geocentre to Moon (true of date) */ |
150
|
7
|
100
|
|
|
|
|
for (i=0; i<6; i++) { |
151
|
6
|
|
|
|
|
|
v[i] = vgm[i]; |
152
|
|
|
|
|
|
|
} |
153
|
|
|
|
|
|
|
|
154
|
|
|
|
|
|
|
} else { |
155
|
|
|
|
|
|
|
|
156
|
|
|
|
|
|
|
/* Not moon: precession/nutation matrix J2000 to date */ |
157
|
8
|
|
|
|
|
|
palPrenut( 2000.0, date, rmat ); |
158
|
|
|
|
|
|
|
|
159
|
|
|
|
|
|
|
/* Sun to Earth-Moon Barycentre (J2000) */ |
160
|
8
|
|
|
|
|
|
palPlanet( date, 3, v, &j ); |
161
|
|
|
|
|
|
|
|
162
|
|
|
|
|
|
|
/* Precession and nutation to date */ |
163
|
8
|
|
|
|
|
|
eraRxp( rmat, v, vse ); |
164
|
8
|
|
|
|
|
|
eraRxp( rmat, &(v[3]), &(vse[3]) ); |
165
|
|
|
|
|
|
|
|
166
|
|
|
|
|
|
|
/* Sun to geocentre (true of date) */ |
167
|
56
|
100
|
|
|
|
|
for (i=0; i<6; i++) { |
168
|
48
|
|
|
|
|
|
vsg[i] = vse[i] - 0.012150581 * vgm[i]; |
169
|
|
|
|
|
|
|
} |
170
|
|
|
|
|
|
|
|
171
|
|
|
|
|
|
|
/* Sun ? */ |
172
|
8
|
100
|
|
|
|
|
if (np == 0) { |
173
|
|
|
|
|
|
|
|
174
|
|
|
|
|
|
|
/* Geocentre to Sun */ |
175
|
7
|
100
|
|
|
|
|
for (i=0; i<6; i++) { |
176
|
6
|
|
|
|
|
|
v[i] = -vsg[i]; |
177
|
|
|
|
|
|
|
} |
178
|
|
|
|
|
|
|
|
179
|
|
|
|
|
|
|
} else { |
180
|
|
|
|
|
|
|
|
181
|
|
|
|
|
|
|
/* Sun to Planet (J2000) */ |
182
|
7
|
|
|
|
|
|
palPlanet( date, np, v, &j ); |
183
|
|
|
|
|
|
|
|
184
|
|
|
|
|
|
|
/* Precession and nutation to date */ |
185
|
7
|
|
|
|
|
|
eraRxp( rmat, v, vsp ); |
186
|
7
|
|
|
|
|
|
eraRxp( rmat, &(v[3]), &(vsp[3]) ); |
187
|
|
|
|
|
|
|
|
188
|
|
|
|
|
|
|
/* Geocentre to planet */ |
189
|
49
|
100
|
|
|
|
|
for (i=0; i<6; i++) { |
190
|
42
|
|
|
|
|
|
v[i] = vsp[i] - vsg[i]; |
191
|
|
|
|
|
|
|
} |
192
|
|
|
|
|
|
|
|
193
|
|
|
|
|
|
|
} |
194
|
|
|
|
|
|
|
|
195
|
|
|
|
|
|
|
} |
196
|
|
|
|
|
|
|
|
197
|
|
|
|
|
|
|
/* Refer to origina at the observer */ |
198
|
9
|
|
|
|
|
|
palPvobs( phi, 0.0, stl, vgo ); |
199
|
63
|
100
|
|
|
|
|
for (i=0; i<6; i++) { |
200
|
54
|
|
|
|
|
|
v[i] -= vgo[i]; |
201
|
|
|
|
|
|
|
} |
202
|
|
|
|
|
|
|
|
203
|
|
|
|
|
|
|
/* Geometric distance (AU) */ |
204
|
9
|
|
|
|
|
|
dx = v[0]; |
205
|
9
|
|
|
|
|
|
dy = v[1]; |
206
|
9
|
|
|
|
|
|
dz = v[2]; |
207
|
9
|
|
|
|
|
|
r = sqrt( dx*dx + dy*dy + dz*dz ); |
208
|
|
|
|
|
|
|
|
209
|
|
|
|
|
|
|
/* Light time */ |
210
|
9
|
|
|
|
|
|
tl = PAL__CR * r; |
211
|
|
|
|
|
|
|
|
212
|
|
|
|
|
|
|
/* Correct position for planetary aberration */ |
213
|
36
|
100
|
|
|
|
|
for (i=0; i<3; i++) { |
214
|
27
|
|
|
|
|
|
v[i] -= tl * v[i+3]; |
215
|
|
|
|
|
|
|
} |
216
|
|
|
|
|
|
|
|
217
|
|
|
|
|
|
|
/* To RA,Dec */ |
218
|
9
|
|
|
|
|
|
eraC2s( v, ra, dec ); |
219
|
9
|
|
|
|
|
|
*ra = eraAnp( *ra ); |
220
|
|
|
|
|
|
|
|
221
|
|
|
|
|
|
|
/* Angular diameter (radians) */ |
222
|
9
|
|
|
|
|
|
*diam = 2.0 * asin( EQRAU[np] / (r * AUKM ) ); |
223
|
|
|
|
|
|
|
|
224
|
9
|
|
|
|
|
|
} |