1 |
// Upgraded to Delphi 2009: Sebastian Zierer
|
2 |
|
3 |
(* ***** BEGIN LICENSE BLOCK *****
|
4 |
* Version: MPL 1.1
|
5 |
*
|
6 |
* The contents of this file are subject to the Mozilla Public License Version
|
7 |
* 1.1 (the "License"); you may not use this file except in compliance with
|
8 |
* the License. You may obtain a copy of the License at
|
9 |
* http://www.mozilla.org/MPL/
|
10 |
*
|
11 |
* Software distributed under the License is distributed on an "AS IS" basis,
|
12 |
* WITHOUT WARRANTY OF ANY KIND, either express or implied. See the License
|
13 |
* for the specific language governing rights and limitations under the
|
14 |
* License.
|
15 |
*
|
16 |
* The Original Code is TurboPower SysTools
|
17 |
*
|
18 |
* The Initial Developer of the Original Code is
|
19 |
* TurboPower Software
|
20 |
*
|
21 |
* Portions created by the Initial Developer are Copyright (C) 1996-2002
|
22 |
* the Initial Developer. All Rights Reserved.
|
23 |
*
|
24 |
* Contributor(s):
|
25 |
*
|
26 |
* ***** END LICENSE BLOCK ***** *)
|
27 |
|
28 |
{*********************************************************}
|
29 |
{* SysTools: StJupsat.pas 4.04 *}
|
30 |
{*********************************************************}
|
31 |
{* SysTools: Astronomical Routines *}
|
32 |
{* (for the four "Gallilean" moons of Jupiter *}
|
33 |
{* Callisto, Europa, Ganymede, and Io) *}
|
34 |
{*********************************************************}
|
35 |
|
36 |
{$I StDefine.inc}
|
37 |
|
38 |
{ ************************************************************** }
|
39 |
{ Sources: }
|
40 |
{ 1. Astronomical Algorithms, Jean Meeus, Willmann-Bell, 1991. }
|
41 |
{ }
|
42 |
{ 2. Planetary and Lunar Coordinates (1984-2000), U.S. Govt, }
|
43 |
{ 1983. }
|
44 |
{ }
|
45 |
{ 3. Supplement to the American Ephemeris and Nautical Almanac,}
|
46 |
{ U.S. Govt, 1964. }
|
47 |
{ }
|
48 |
{ 4. MPO96-98 source files, Brian D. Warner, 1995-98. }
|
49 |
{ }
|
50 |
{ ************************************************************** }
|
51 |
|
52 |
(* **************************************************************
|
53 |
|
54 |
The formulae in this unit are based on DYNAMICAL time, which is based
|
55 |
on the actual rotation of the Earth and is gradually slowing. Universal
|
56 |
Time is based on a fixed rotation rate. To directly compare results in
|
57 |
the astronomical literature (and Meeus' examples), you must use a Universal
|
58 |
Time that is adjusted by the value Delta T. This value is approximately 1 min
|
59 |
in the latter part of the 20th century and will be about 80 seconds in the
|
60 |
year 2020. As an example, to compare the High Precision positions for
|
61 |
1992 December 16 (Meeus' example), you must use a Universal Time of
|
62 |
1992 December 16 at 00:00:59 which equals 00:00:00 Dynamical Time
|
63 |
|
64 |
|
65 |
The Shadows parameter is used for high precision calculations only. If True,
|
66 |
the positions are calculated as seen from the SUN, not the Earth. For eclipses
|
67 |
of the satellites by Jupiter, in effect, the position is in reference to the
|
68 |
SHADOW of Jupiter in space, not the planet itself. For shadow transits, where
|
69 |
the shadow of the satellite is projected onto the surface of the planet, the
|
70 |
position is that of the satellite's shadow in reference to Jupiter and not
|
71 |
the satellite itself.
|
72 |
|
73 |
The purpose of the Shadows parameter is to complete the prediction for
|
74 |
satellite phenomenon. For example, using Shadow := False, the result may
|
75 |
indicate that a satellite goes behind Jupiter at a given instant but not
|
76 |
if the satellite is visible because it is in Jupiter's shadow. Setting
|
77 |
Shadow := True for the same time will indicate if the planet is in or out
|
78 |
of Jupiter's shadow.
|
79 |
|
80 |
(Shadow := FALSE) and (abs(satellite X-coordinate) = 1)
|
81 |
-------------------------------------------------------
|
82 |
If the X-value is negative and heading towards 0, the satellite is entering
|
83 |
the front of the planet.
|
84 |
If the X-value is negative and increasing, the satellite is coming from
|
85 |
behind the planet.
|
86 |
If the X-value is positive and heading towards 0, the satellite is going
|
87 |
behind the planet.
|
88 |
If the X-value is positive and increasing, the satellite is leaving
|
89 |
the front of the planet.
|
90 |
|
91 |
(Shadow := TRUE) and (abs(satellite X-coordinate) = 1)
|
92 |
-------------------------------------------------------
|
93 |
If the X-value is negative and heading towards 0, the satellite's shadow is
|
94 |
entering the planet's disc.
|
95 |
If the X-value is negative and increasing, the satellite is leaving the
|
96 |
planet's shadow.
|
97 |
If the X-value is positive and heading towards 0, the satellite entering
|
98 |
the planet's shadow.
|
99 |
If the X-value is positive and increasing, the satellite's shadow is
|
100 |
leaving the planet.
|
101 |
|
102 |
The X and Y coordinates are based on the equatorial radius of Jupiter. Because
|
103 |
the planet is considerably flattened by its rapid rotation, the polar diameter
|
104 |
is less than 1. To avoid dealing with an elliptical disc for Jupiter, multiply
|
105 |
the Y values only by 1.071374. This creates a "circular Jupiter" and so makes
|
106 |
determining if the satellite is above or below Jupiter easier
|
107 |
(abs(Y-coordinate) > 1).
|
108 |
|
109 |
****************************************************************** *)
|
110 |
|
111 |
unit StJupsat;
|
112 |
|
113 |
interface
|
114 |
|
115 |
type
|
116 |
TStJupSatPos = packed record
|
117 |
X : Double;
|
118 |
Y : Double;
|
119 |
end;
|
120 |
|
121 |
TStJupSats = packed record
|
122 |
Io : TStJupSatPos;
|
123 |
Europa : TStJupSatPos;
|
124 |
Ganymede : TStJupSatPos;
|
125 |
Callisto : TStJupSatPos;
|
126 |
end;
|
127 |
|
128 |
function GetJupSats(JD : TDateTime; HighPrecision, Shadows : Boolean) : TStJupSats;
|
129 |
|
130 |
implementation
|
131 |
uses
|
132 |
StDate, StAstro, StAstroP, StJup, StMath;
|
133 |
|
134 |
type
|
135 |
SunCoordsRec = packed record
|
136 |
X, Y, Z : Double;
|
137 |
L, B, R : Double;
|
138 |
end;
|
139 |
|
140 |
TranformRec = packed record
|
141 |
A, B, C : array[1..6] of Double;
|
142 |
end;
|
143 |
|
144 |
|
145 |
function SunCoords(JD : Double) : SunCoordsRec;
|
146 |
var
|
147 |
L, B, R,
|
148 |
T0, TM,
|
149 |
RS,
|
150 |
OB, A : Double;
|
151 |
begin
|
152 |
T0 := (JD - StdDate) / 365250;
|
153 |
TM := T0/100;
|
154 |
RS := radcor * 3600;
|
155 |
OB := 0.4090928042223
|
156 |
- 4680.93/RS * TM
|
157 |
- 1.55/RS * sqr(TM)
|
158 |
+ 1999.25/RS * sqr(TM) * TM
|
159 |
- 51.38/RS * sqr(sqr(TM))
|
160 |
- 249.67/RS * sqr(sqr(TM)) * TM
|
161 |
- 39.05/RS * sqr(sqr(TM)) * sqr(TM)
|
162 |
+ 7.12/RS * sqr(sqr(TM)) * sqr(TM) * TM
|
163 |
+ 27.87/RS * sqr(sqr(sqr(TM)));
|
164 |
|
165 |
L := 175347046
|
166 |
+ 3341656 * cos(4.6692568 + 6283.07585*T0)
|
167 |
+ 34894 * cos(4.6261000 + 12566.1517*T0)
|
168 |
+ 3497 * cos(2.7441000 + 5753.3849*T0)
|
169 |
+ 3418 * cos(2.8289000 + 3.5231*T0)
|
170 |
+ 3136 * cos(3.6277000 + 77713.7715*T0)
|
171 |
+ 2676 * cos(4.4181000 + 7860.4194*T0)
|
172 |
+ 2343 * cos(6.1352000 + 3930.2097*T0)
|
173 |
+ 1324 * cos(0.7425000 + 11506.7698*T0)
|
174 |
+ 1273 * cos(2.0371000 + 529.6910*T0)
|
175 |
+ 1199 * cos(1.1096000 + 1577.3435*T0)
|
176 |
+ 990 * cos(5.2330000 + 5884.9270*T0)
|
177 |
+ 902 * cos(2.0450000 + 26.1490*T0)
|
178 |
+ 857 * cos(3.5080000 + 398.149*T0)
|
179 |
+ 780 * cos(1.1790000 + 5223.694*T0)
|
180 |
+ 753 * cos(2.5330000 + 5507.553*T0)
|
181 |
+ 505 * cos(4.5830000 + 18849.228*T0)
|
182 |
+ 492 * cos(4.2050000 + 775.523*T0)
|
183 |
+ 357 * cos(2.9200000 + 0.067*T0)
|
184 |
+ 317 * cos(5.8490000 + 11790.626*T0)
|
185 |
+ 284 * cos(1.8990000 + 796.298*T0)
|
186 |
+ 271 * cos(0.3150000 + 10977.079*T0)
|
187 |
+ 243 * cos(0.3450000 + 5486.778*T0)
|
188 |
+ 206 * cos(4.8060000 + 2544.314*T0)
|
189 |
+ 205 * cos(1.8690000 + 5573.143*T0)
|
190 |
+ 202 * cos(2.4580000 + 6069.777*T0)
|
191 |
+ 156 * cos(0.8330000 + 213.299*T0)
|
192 |
+ 132 * cos(3.4110000 + 2942.463*T0)
|
193 |
+ 126 * cos(1.0830000 + 20.775*T0)
|
194 |
+ 115 * cos(0.6450000 + 0.980*T0)
|
195 |
+ 103 * cos(0.6360000 + 4694.003*T0)
|
196 |
+ 102 * cos(0.9760000 + 15720.839*T0)
|
197 |
+ 102 * cos(4.2670000 + 7.114*T0)
|
198 |
+ 99 * cos(6.2100000 + 2146.170*T0)
|
199 |
+ 98 * cos(0.6800000 + 155.420*T0)
|
200 |
+ 86 * cos(5.9800000 +161000.690*T0)
|
201 |
+ 85 * cos(1.3000000 + 6275.960*T0)
|
202 |
+ 85 * cos(3.6700000 + 71430.700*T0)
|
203 |
+ 80 * cos(1.8100000 + 17260.150*T0);
|
204 |
|
205 |
A := 628331966747.0
|
206 |
+ 206059 * cos(2.678235 + 6283.07585*T0)
|
207 |
+ 4303 * cos(2.635100 + 12566.1517*T0)
|
208 |
+ 425 * cos(1.590000 + 3.523*T0)
|
209 |
+ 119 * cos(5.796000 + 26.298*T0)
|
210 |
+ 109 * cos(2.966000 + 1577.344*T0)
|
211 |
+ 93 * cos(2.590000 + 18849.23*T0)
|
212 |
+ 72 * cos(1.140000 + 529.69*T0)
|
213 |
+ 68 * cos(1.870000 + 398.15*T0)
|
214 |
+ 67 * cos(4.410000 + 5507.55*T0)
|
215 |
+ 59 * cos(2.890000 + 5223.69*T0)
|
216 |
+ 56 * cos(2.170000 + 155.42*T0)
|
217 |
+ 45 * cos(0.400000 + 796.30*T0)
|
218 |
+ 36 * cos(0.470000 + 775.52*T0)
|
219 |
+ 29 * cos(2.650000 + 7.11*T0)
|
220 |
+ 21 * cos(5.340000 + 0.98*T0)
|
221 |
+ 19 * cos(1.850000 + 5486.78*T0)
|
222 |
+ 19 * cos(4.970000 + 213.30*T0)
|
223 |
+ 17 * cos(2.990000 + 6275.96*T0)
|
224 |
+ 16 * cos(0.030000 + 2544.31*T0);
|
225 |
L := L + (A * T0);
|
226 |
|
227 |
A := 52919
|
228 |
+ 8720 * cos(1.0721 + 6283.0758*T0)
|
229 |
+ 309 * cos(0.867 + 12566.152*T0)
|
230 |
+ 27 * cos(0.050 + 3.52*T0)
|
231 |
+ 16 * cos(5.190 + 26.30*T0)
|
232 |
+ 16 * cos(3.68 + 155.42*T0)
|
233 |
+ 10 * cos(0.76 + 18849.23*T0)
|
234 |
+ 9 * cos(2.06 + 77713.77*T0)
|
235 |
+ 7 * cos(0.83 + 775.52*T0)
|
236 |
+ 5 * cos(4.66 + 1577.34*T0);
|
237 |
L := L + (A * sqr(T0));
|
238 |
|
239 |
A := 289 * cos(5.844 + 6283.076*T0)
|
240 |
+ 35
|
241 |
+ 17 * cos(5.49 + 12566.15*T0)
|
242 |
+ 3 * cos(5.20 + 155.42*T0)
|
243 |
+ 1 * cos(4.72 + 3.52*T0);
|
244 |
L := L + (A * sqr(T0) * T0);
|
245 |
|
246 |
A := 114 * cos(3.142);
|
247 |
L := L + (A * sqr(sqr(T0)));
|
248 |
L := L / 1.0E+8;
|
249 |
|
250 |
{solar latitude}
|
251 |
B := 280 * cos(3.199 + 84334.662*T0)
|
252 |
+ 102 * cos(5.422 + 5507.553*T0)
|
253 |
+ 80 * cos(3.88 + 5223.69*T0)
|
254 |
+ 44 * cos(3.70 + 2352.87*T0)
|
255 |
+ 32 * cos(4.00 + 1577.34*T0);
|
256 |
|
257 |
A := 9 * cos(3.90 + 5507.550*T0)
|
258 |
+ 6 * cos(1.73 + 5223.690*T0);
|
259 |
B := B + (A * T0);
|
260 |
B := B / 1.0E+8;
|
261 |
|
262 |
|
263 |
{solar radius vector (astronomical units)}
|
264 |
R := 100013989
|
265 |
+ 1670700 * cos(3.0984635 + 6283.07585*T0)
|
266 |
+ 13956 * cos(3.05525 + 12566.15170*T0)
|
267 |
+ 3084 * cos(5.1985 + 77713.7715*T0)
|
268 |
+ 1628 * cos(1.1739 + 5753.3849*T0)
|
269 |
+ 1576 * cos(2.8649 + 7860.4194*T0)
|
270 |
+ 925 * cos(5.453 + 11506.770*T0)
|
271 |
+ 542 * cos(4.564 + 3930.210*T0)
|
272 |
+ 472 * cos(3.661 + 5884.927*T0)
|
273 |
+ 346 * cos(0.964 + 5507.553*T0)
|
274 |
+ 329 * cos(5.900 + 5223.694*T0)
|
275 |
+ 307 * cos(0.299 + 5573.143*T0)
|
276 |
+ 243 * cos(4.273 + 11790.629*T0)
|
277 |
+ 212 * cos(5.847 + 1577.344*T0)
|
278 |
+ 186 * cos(5.022 + 10977.079*T0)
|
279 |
+ 175 * cos(3.012 + 18849.228*T0)
|
280 |
+ 110 * cos(5.055 + 5486.778*T0)
|
281 |
+ 98 * cos(0.89 + 6069.78*T0)
|
282 |
+ 86 * cos(5.69 + 15720.84*T0)
|
283 |
+ 86 * cos(1.27 +161000.69*T0)
|
284 |
+ 65 * cos(0.27 + 17260.15*T0)
|
285 |
+ 63 * cos(0.92 + 529.69*T0)
|
286 |
+ 57 * cos(2.01 + 83996.85*T0)
|
287 |
+ 56 * cos(5.24 + 71430.70*T0)
|
288 |
+ 49 * cos(3.25 + 2544.31*T0)
|
289 |
+ 47 * cos(2.58 + 775.52*T0)
|
290 |
+ 45 * cos(5.54 + 9437.76*T0)
|
291 |
+ 43 * cos(6.01 + 6275.96*T0)
|
292 |
+ 39 * cos(5.36 + 4694.00*T0)
|
293 |
+ 38 * cos(2.39 + 8827.39*T0)
|
294 |
+ 37 * cos(0.83 + 19651.05*T0)
|
295 |
+ 37 * cos(4.90 + 12139.55*T0)
|
296 |
+ 36 * cos(1.67 + 12036.46*T0)
|
297 |
+ 35 * cos(1.84 + 2942.46*T0)
|
298 |
+ 33 * cos(0.24 + 7084.90*T0)
|
299 |
+ 32 * cos(0.18 + 5088.63*T0)
|
300 |
+ 32 * cos(1.78 + 398.15*T0)
|
301 |
+ 28 * cos(1.21 + 6286.60*T0)
|
302 |
+ 28 * cos(1.90 + 6279.55*T0)
|
303 |
+ 26 * cos(4.59 + 10447.39*T0);
|
304 |
R := R / 1.0E+8;
|
305 |
|
306 |
A := 103019 * cos(1.107490 + 6283.075850*T0)
|
307 |
+ 1721 * cos(1.0644 + 12566.1517*T0)
|
308 |
+ 702 * cos(3.142)
|
309 |
+ 32 * cos(1.02 + 18849.23*T0)
|
310 |
+ 31 * cos(2.84 + 5507.55*T0)
|
311 |
+ 25 * cos(1.32 + 5223.69*T0)
|
312 |
+ 18 * cos(1.42 + 1577.34*T0)
|
313 |
+ 10 * cos(5.91 + 10977.08*T0)
|
314 |
+ 9 * cos(1.42 + 6275.96*T0)
|
315 |
+ 9 * cos(0.27 + 5486.78*T0);
|
316 |
R := R + (A * T0 / 1.0E+8);
|
317 |
|
318 |
A := 4359 * cos(5.7846 + 6283.0758*T0)
|
319 |
+ 124 * cos(5.579 + 12566.152*T0)
|
320 |
+ 12 * cos(3.14)
|
321 |
+ 9 * cos(3.63 + 77713.77*T0)
|
322 |
+ 6 * cos(1.87 + 5573.14*T0)
|
323 |
+ 3 * cos(5.47 + 18849.23*T0);
|
324 |
R := R + (A * sqr(T0) / 1.0E+8);
|
325 |
|
326 |
L := (L + PI);
|
327 |
L := Frac(L / 2.0 / PI) * 2.0 * Pi;
|
328 |
if L < 0 then
|
329 |
L := L + (2.0*PI);
|
330 |
B := -B;
|
331 |
|
332 |
Result.L := L;
|
333 |
Result.B := B;
|
334 |
Result.R := R;
|
335 |
Result.X := R * cos(B) * cos(L);
|
336 |
Result.Y := R * (cos(B) * sin(L) * cos(OB) - sin(B) * sin(OB));
|
337 |
Result.Z := R * (cos(B) * sin(L) * sin(OB) + sin(B) * cos(OB));
|
338 |
end;
|
339 |
|
340 |
|
341 |
{-------------------------------------------------------------------------}
|
342 |
|
343 |
function JupSatsLo(AJD : Double) : TStJupSats;
|
344 |
var
|
345 |
DateDif, {d}
|
346 |
ArgJup, {V}
|
347 |
AnomE, {M}
|
348 |
AnomJ, {N}
|
349 |
DeltaLong, {J}
|
350 |
ECenterE, {A}
|
351 |
ECenterJ, {B}
|
352 |
K,
|
353 |
RVE, {R}
|
354 |
RVJ, {r}
|
355 |
EJDist, {Delta}
|
356 |
Phase, {Psi}
|
357 |
Lambda, {Lambda}
|
358 |
DS, DE, {DS, DE}
|
359 |
Mu1, Mu2,
|
360 |
Mu3, Mu4, {Mu1 - Mu4}
|
361 |
G, H, {G, H}
|
362 |
TmpDbl1,
|
363 |
TmpDbl2,
|
364 |
R1, R2,
|
365 |
R3, R4 {R1 - R4}
|
366 |
: Double;
|
367 |
|
368 |
begin
|
369 |
AJD := DateTimeToAJD(AJD);
|
370 |
DateDif := AJD - 2451545.0;
|
371 |
ArgJup := 172.74 + (0.00111588 * DateDif);
|
372 |
ArgJup := Frac(ArgJup/360.0) * 360.0;
|
373 |
if (ArgJup < 0) then
|
374 |
ArgJup := 360.0 + ArgJup;
|
375 |
ArgJup := ArgJup / radcor;
|
376 |
|
377 |
AnomE := 357.529 + (0.9856003 * DateDif);
|
378 |
AnomE := Frac(AnomE/360.0) * 360.0;
|
379 |
if (AnomE < 0) then
|
380 |
AnomE := 360.0 + AnomE;
|
381 |
AnomE := AnomE / radcor;
|
382 |
|
383 |
AnomJ := 20.020 + (0.0830853 * DateDif + (0.329 * sin(ArgJup)));
|
384 |
AnomJ := Frac(AnomJ/360.0) * 360.0;
|
385 |
if (AnomJ < 0) then
|
386 |
AnomJ := 360.0 + AnomJ;
|
387 |
AnomJ := AnomJ / radcor;
|
388 |
|
389 |
DeltaLong := 66.115 + (0.9025179 * DateDif - (0.329 * sin(ArgJup)));
|
390 |
DeltaLong := Frac(DeltaLong/360.0) * 360.0;
|
391 |
if (DeltaLong < 0) then
|
392 |
DeltaLong := 360.0 + DeltaLong;
|
393 |
DeltaLong := DeltaLong / radcor;
|
394 |
|
395 |
ECenterE := 1.915 * sin(AnomE) + 0.020 * sin(2*AnomE);
|
396 |
ECenterE := ECenterE / radcor;
|
397 |
|
398 |
ECenterJ := 5.555 * sin(AnomJ) + 0.168 * sin(2*AnomJ);
|
399 |
ECenterJ := ECenterJ / radcor;
|
400 |
|
401 |
K := (DeltaLong + ECenterE - ECenterJ);
|
402 |
|
403 |
RVE := 1.00014 - (0.01671 * cos(AnomE)) - (0.00014 * cos(2*AnomE));
|
404 |
RVJ := 5.20872 - (0.25208 * cos(AnomJ)) - (0.00611 * cos(2*AnomJ));
|
405 |
|
406 |
EJDist := sqrt(sqr(RVJ) + sqr(RVE) - (2 * RVJ * RVE * cos(K)));
|
407 |
|
408 |
Phase := RVE/EJDist * sin(K);
|
409 |
Phase := StInvSin(Phase);
|
410 |
|
411 |
if ((sin(K) < 0) and (Phase > 0)) or
|
412 |
((sin(K) > 0) and (Phase < 0)) then
|
413 |
Phase := -Phase;
|
414 |
|
415 |
Lambda := 34.35 + (0.083091 * DateDif) + (0.329 * sin(ArgJup));
|
416 |
Lambda := Lambda / radcor + ECenterJ;
|
417 |
|
418 |
DS := 3.12 * sin(Lambda + 42.8 / radcor);
|
419 |
DE := DS - 2.22 * sin(Phase) * cos(Lambda + 22/radcor)
|
420 |
- 1.30 * ((RVJ - EJDist) / EJDist) * sin(Lambda - 100.5/radcor);
|
421 |
DE := DE / radcor;
|
422 |
|
423 |
Mu1 := 163.8067 + 203.4058643 * (DateDif - (EJDist / 173));
|
424 |
Mu1 := Frac(Mu1/360.0) * 360.0;
|
425 |
if (Mu1 < 0) then
|
426 |
Mu1 := 360.0 + Mu1;
|
427 |
Mu1 := Mu1 / radcor + Phase - ECenterJ;
|
428 |
|
429 |
Mu2 := 358.4108 + 101.2916334 * (DateDif - (EJDist / 173));
|
430 |
Mu2 := Frac(Mu2/360.0) * 360.0;
|
431 |
if (Mu2 < 0) then
|
432 |
Mu2 := 360.0 + Mu2;
|
433 |
Mu2 := Mu2 / radcor + Phase - ECenterJ;
|
434 |
|
435 |
Mu3 := 5.7129 + 50.2345179 * (DateDif - (EJDist / 173));
|
436 |
Mu3 := Frac(Mu3/360.0) * 360.0;
|
437 |
if (Mu3 < 0) then
|
438 |
Mu3 := 360.0 + Mu3;
|
439 |
Mu3 := Mu3 / radcor + Phase - ECenterJ;
|
440 |
|
441 |
Mu4 := 224.8151 + 21.4879801 * (DateDif - (EJDist / 173));
|
442 |
Mu4 := Frac(Mu4/360.0) * 360.0;
|
443 |
if (Mu4 < 0) then
|
444 |
Mu4 := 360.0 + Mu4;
|
445 |
Mu4 := Mu4 / radcor + Phase - ECenterJ;
|
446 |
|
447 |
G := 331.18 + 50.310482 * (DateDif - (EJDist / 173));
|
448 |
G := Frac(G/360.0) * 360.0;
|
449 |
if (G < 0) then
|
450 |
G := 360.0 + G;
|
451 |
G := G / radcor;
|
452 |
H := 87.40 + 21.569231 * (DateDif - (EJDist / 173));
|
453 |
H := Frac(H/360.0) * 360.0;
|
454 |
if (H < 0) then
|
455 |
H := 360.0 + H;
|
456 |
H := H / radcor;
|
457 |
|
458 |
TmpDbl1 := 0.473 * sin(2 * (Mu1 - Mu2)) / radcor;
|
459 |
TmpDbl2 := 1.065 * sin(2 * (Mu2 - Mu3)) / radcor;
|
460 |
|
461 |
R1 := 5.9073 - 0.0244 * cos(2 * (Mu1 - Mu2));
|
462 |
R2 := 9.3991 - 0.0882 * cos(2 * (Mu2 - Mu3));
|
463 |
R3 := 14.9924 - 0.0216 * cos(G);
|
464 |
R4 := 26.3699 - 0.1935 * cos(H);
|
465 |
|
466 |
Mu1 := Mu1 + TmpDbl1;
|
467 |
Mu2 := Mu2 + TmpDbl2;
|
468 |
Mu3 := Mu3 + (0.165 * sin(G)) / radcor;
|
469 |
Mu4 := Mu4 + (0.841 * sin(H)) / radcor;
|
470 |
|
471 |
Result.Io.X := R1 * sin(Mu1);
|
472 |
Result.Io.Y := -R1 * cos(Mu1) * sin(DE);
|
473 |
|
474 |
Result.Europa.X := R2 * sin(Mu2);
|
475 |
Result.Europa.Y := -R2 * cos(Mu2) * sin(DE);
|
476 |
|
477 |
Result.Ganymede.X := R3 * sin(Mu3);
|
478 |
Result.Ganymede.Y := -R3 * cos(Mu3) * sin(DE);
|
479 |
|
480 |
Result.Callisto.X := R4 * sin(Mu4);
|
481 |
Result.Callisto.Y := -R4 * cos(Mu4) * sin(DE);
|
482 |
end;
|
483 |
|
484 |
{-------------------------------------------------------------------------}
|
485 |
|
486 |
function JupSatsHi(AJD : Double; Shadows : Boolean) : TStJupSats;
|
487 |
var
|
488 |
I : longint;
|
489 |
SunPos : SunCoordsRec;
|
490 |
STUT : TStDateTimeRec;
|
491 |
JupPos : TStEclipticalCord;
|
492 |
|
493 |
SatX : array[1..5] of Double;
|
494 |
SatY : array[1..5] of Double;
|
495 |
SatZ : array[1..5] of Double;
|
496 |
|
497 |
TD1,
|
498 |
TD2,
|
499 |
Angle, {Temporary Double values}
|
500 |
LTime, {Tau}
|
501 |
AJDT, {AJD adjusted for light time (Tau)}
|
502 |
JupX,
|
503 |
JupY,
|
504 |
JupZ, {Jupiter's geocentric rectangular coordinates}
|
505 |
EJDist, {Delta}
|
506 |
Jup1,
|
507 |
Jup2, {/\, Alpha}
|
508 |
DateDif, {t}
|
509 |
L1, L2,
|
510 |
L3, L4, {script L1-4}
|
511 |
Pi1, Pi2,
|
512 |
Pi3, Pi4, {Pi1-4}
|
513 |
W1, W2,
|
514 |
W3, W4, {Omega1-4}
|
515 |
Inequality, {upside down L}
|
516 |
PhiLambda,
|
517 |
NodeJup, {Psi}
|
518 |
AnomJup, {G}
|
519 |
AnomSat, {G'}
|
520 |
LongPerJ,
|
521 |
S1, S2,
|
522 |
S3, S4, {Sum1-4}
|
523 |
TL1, TL2,
|
524 |
TL3, TL4, {Capital L1-4}
|
525 |
B1, B2,
|
526 |
B3, B4, {tangent of latitude}
|
527 |
R1, R2,
|
528 |
R3, R4, {radius vector}
|
529 |
T0, {Julian Centuries}
|
530 |
Precession, {P}
|
531 |
Inclination {I}
|
532 |
|
533 |
: Double;
|
534 |
Transforms : array[1..5] of TranformRec;
|
535 |
|
536 |
begin
|
537 |
FillChar(Result, SizeOf(TStJupSats), #0);
|
538 |
AJD := DateTimeToAJD(AJD);
|
539 |
SunPos := SunCoords(AJD);
|
540 |
|
541 |
if not Shadows then begin
|
542 |
TD1 := 5;
|
543 |
AJDT := AJD - 0.0057755183 * TD1; {first guess}
|
544 |
repeat
|
545 |
JupPos := ComputeJupiter(AJDT);
|
546 |
|
547 |
JupX := JupPos.R0 * cos(JupPos.B0) * cos(JupPos.L0)
|
548 |
+ SunPos.R * cos(SunPos.L);
|
549 |
JupY := JupPos.R0 * cos(JupPos.B0) * sin(JupPos.L0)
|
550 |
+ SunPos.R * sin(SunPos.L);
|
551 |
JupZ := JupPos.R0 * sin(JupPos.B0);
|
552 |
|
553 |
EJDist := sqrt(sqr(JupX) + sqr(JupY) + sqr(JupZ));
|
554 |
TD2 := abs(EJDist - TD1);
|
555 |
if abs(TD2) > 0.0005 then begin
|
556 |
AJDT := AJD - 0.0057755183 * ((EJDist + TD1) / 2);
|
557 |
TD1 := EJDist;
|
558 |
end;
|
559 |
until (TD2 <= 0.0005);
|
560 |
end else begin
|
561 |
JupPos := ComputeJupiter(AJD);
|
562 |
|
563 |
JupX := JupPos.R0 * cos(JupPos.B0) * cos(JupPos.L0);
|
564 |
JupY := JupPos.R0 * cos(JupPos.B0) * sin(JupPos.L0);
|
565 |
JupZ := JupPos.R0 * sin(JupPos.B0);
|
566 |
EJDist := sqrt(sqr(JupX+SunPos.X) +
|
567 |
sqr(JupY+SunPos.Y) + sqr(JupZ+SunPos.Z));
|
568 |
end;
|
569 |
|
570 |
Jup1 := StInvTan2(JupX, JupY);
|
571 |
Jup2 := ArcTan(JupZ / sqrt(sqr(JupX) + sqr(JupY)));
|
572 |
|
573 |
DateDif := AJD - 2443000.5 - (0.0057755183 * EJDist);
|
574 |
|
575 |
L1 := 106.07947 + 203.488955432 * DateDif;
|
576 |
L1 := Frac(L1/360.0) * 360.0;
|
577 |
if (L1 < 0) then
|
578 |
L1 := 360.0 + L1;
|
579 |
L1 := L1 / radcor;
|
580 |
|
581 |
L2 := 175.72938 + 101.374724550 * DateDif;
|
582 |
L2 := Frac(L2/360.0) * 360.0;
|
583 |
if (L2 < 0) then
|
584 |
L2 := 360.0 + L2;
|
585 |
L2 := L2 / radcor;
|
586 |
|
587 |
L3 := 120.55434 + 50.317609110 * DateDif;
|
588 |
L3 := Frac(L3/360.0) * 360.0;
|
589 |
if (L3 < 0) then
|
590 |
L3 := 360.0 + L3;
|
591 |
L3 := L3 / radcor;
|
592 |
|
593 |
L4 := 84.44868 + 21.571071314 * DateDif;
|
594 |
L4 := Frac(L4/360.0) * 360.0;
|
595 |
if (L4 < 0) then
|
596 |
L4 := 360.0 + L4;
|
597 |
L4 := L4 / radcor;
|
598 |
|
599 |
Pi1 := 58.3329 + 0.16103936 * DateDif;
|
600 |
Pi1 := Frac(Pi1/360.0) * 360.0;
|
601 |
if (Pi1 < 0) then
|
602 |
Pi1 := 360.0 + Pi1;
|
603 |
Pi1 := Pi1 / radcor;
|
604 |
|
605 |
Pi2 := 132.8959 + 0.04647985 * DateDif;
|
606 |
Pi2 := Frac(Pi2/360.0) * 360.0;
|
607 |
if (Pi2 < 0) then
|
608 |
Pi2 := 360.0 + Pi2;
|
609 |
Pi2 := Pi2 / radcor;
|
610 |
|
611 |
Pi3 := 187.2887 + 0.00712740 * DateDif;
|
612 |
Pi3 := Frac(Pi3/360.0) * 360.0;
|
613 |
if (Pi3 < 0) then
|
614 |
Pi3 := 360.0 + Pi3;
|
615 |
Pi3 := Pi3 / radcor;
|
616 |
|
617 |
Pi4 := 335.3418 + 0.00183998 * DateDif;
|
618 |
Pi4 := Frac(Pi4/360.0) * 360.0;
|
619 |
if (Pi4 < 0) then
|
620 |
Pi4 := 360.0 + Pi4;
|
621 |
Pi4 := Pi4 / radcor;
|
622 |
|
623 |
W1 := 311.0793 - 0.13279430 * DateDif;
|
624 |
W1 := Frac(W1/360.0) * 360.0;
|
625 |
if (W1 < 0) then
|
626 |
W1 := 360.0 + W1;
|
627 |
W1 := W1 / radcor;
|
628 |
|
629 |
W2 := 100.5099 - 0.03263047 * DateDif;
|
630 |
W2 := Frac(W2/360.0) * 360.0;
|
631 |
if (W2 < 0) then
|
632 |
W2 := 360.0 + W2;
|
633 |
W2 := W2 / radcor;
|
634 |
|
635 |
W3 := 119.1688 - 0.00717704 * DateDif;
|
636 |
W3 := Frac(W3/360.0) * 360.0;
|
637 |
if (W3 < 0) then
|
638 |
W3 := 360.0 + W3;
|
639 |
W3 := W3 / radcor;
|
640 |
|
641 |
W4 := 322.5729 - 0.00175934 * DateDif;
|
642 |
W4 := Frac(W4/360.0) * 360.0;
|
643 |
if (W4 < 0) then
|
644 |
W4 := 360.0 + W4;
|
645 |
W4 := W4 / radcor;
|
646 |
|
647 |
Inequality := 0.33033 * sin((163.679 + 0.0010512*DateDif) / radcor)
|
648 |
+ 0.03439 * sin((34.486 - 0.0161731*DateDif) / radcor);
|
649 |
Inequality := Inequality / radcor;
|
650 |
|
651 |
PhiLambda := 191.8132 + 0.17390023 * DateDif;
|
652 |
PhiLambda := Frac(PhiLambda / 360.0) * 360.0;
|
653 |
if (PhiLambda < 0) then
|
654 |
PhiLambda := 360.0 + PhiLambda;
|
655 |
PhiLambda := PhiLambda / radcor;
|
656 |
|
657 |
NodeJup := 316.5182 - 0.00000208 * DateDif;
|
658 |
NodeJup := Frac(NodeJup / 360.0) * 360.0;
|
659 |
if (NodeJup < 0) then
|
660 |
NodeJup := 360.0 + NodeJup;
|
661 |
NodeJup := NodeJup / radcor;
|
662 |
|
663 |
AnomJup := 30.23756 + 0.0830925701 * DateDif;
|
664 |
AnomJup := Frac(AnomJup / 360.0) * 360.0;
|
665 |
if (AnomJup < 0) then
|
666 |
AnomJup := 360.0 + AnomJup;
|
667 |
AnomJup := AnomJup / radcor + Inequality;
|
668 |
|
669 |
AnomSat := 31.97853 + 0.0334597339 * DateDif;
|
670 |
AnomSat := Frac(AnomSat / 360.0) * 360.0;
|
671 |
if (AnomSat < 0) then
|
672 |
AnomSat := 360.0 + AnomSat;
|
673 |
AnomSat := AnomSat / radcor;
|
674 |
|
675 |
LongPerJ := 13.469942 / radcor;
|
676 |
|
677 |
S1 := 0.47259 * sin(2*(L1-L2))
|
678 |
- 0.03480 * sin(Pi3-Pi4)
|
679 |
- 0.01756 * sin(Pi1 + Pi3 - 2*LongPerJ - 2*AnomJup)
|
680 |
+ 0.01080 * sin(L2 - 2*L3 + Pi3)
|
681 |
+ 0.00757 * sin(PhiLambda)
|
682 |
+ 0.00663 * sin(L2 - 2*L3 + Pi4)
|
683 |
+ 0.00453 * sin(L1 - Pi3)
|
684 |
+ 0.00453 * sin(L2 - 2*L3 + Pi2)
|
685 |
- 0.00354 * sin(L1-L2)
|
686 |
- 0.00317 * sin(2*NodeJup - 2*LongPerJ)
|
687 |
- 0.00269 * sin(L2 - 2*L3 + Pi1)
|
688 |
+ 0.00263 * sin(L1 - Pi4)
|
689 |
+ 0.00186 * sin(L1 - Pi1)
|
690 |
- 0.00186 * sin(AnomJup)
|
691 |
+ 0.00167 * sin(Pi2 - Pi3)
|
692 |
+ 0.00158 * sin(4*(L1-L2))
|
693 |
- 0.00155 * sin(L1 - L3)
|
694 |
- 0.00142 * sin(NodeJup + W3 - 2*LongPerJ - 2*AnomJup)
|
695 |
- 0.00115 * sin(2*(L1 - 2*L2 + W2))
|
696 |
+ 0.00089 * sin(Pi2 - Pi4)
|
697 |
+ 0.00084 * sin(W2 - W3)
|
698 |
+ 0.00084 * sin(L1 + Pi3 - 2*LongPerJ - 2*AnomJup)
|
699 |
+ 0.00053 * sin(NodeJup - W2);
|
700 |
|
701 |
S2 := 1.06476 * sin(2*(L2-L3))
|
702 |
+ 0.04253 * sin(L1 - 2*L2 + Pi3)
|
703 |
+ 0.03579 * sin(L2 - Pi3)
|
704 |
+ 0.02383 * sin(L1 - 2*L2 + Pi4)
|
705 |
+ 0.01977 * sin(L2 - Pi4)
|
706 |
- 0.01843 * sin(PhiLambda)
|
707 |
+ 0.01299 * sin(Pi3 - Pi4)
|
708 |
- 0.01142 * sin(L2 - L3)
|
709 |
+ 0.01078 * sin(L2 - Pi2)
|
710 |
- 0.01058 * sin(AnomJup)
|
711 |
+ 0.00870 * sin(L2 - 2*L3 + Pi2)
|
712 |
- 0.00775 * sin(2*(NodeJup - LongPerJ))
|
713 |
+ 0.00524 * sin(2*(L1-L2))
|
714 |
- 0.00460 * sin(L1-L3)
|
715 |
+ 0.00450 * sin(L2 - 2*L3 + Pi1)
|
716 |
+ 0.00327 * sin(NodeJup - 2*AnomJup + W3 - 2*LongPerJ)
|
717 |
- 0.00296 * sin(Pi1 + Pi3 - 2*LongPerJ - 2*AnomJup)
|
718 |
- 0.00151 * sin(2*AnomJup)
|
719 |
+ 0.00146 * sin(NodeJup - W3)
|
720 |
+ 0.00125 * sin(NodeJup - W4)
|
721 |
- 0.00117 * sin(L1 -2*L3 + Pi3)
|
722 |
- 0.00095 * sin(2*(L2-W2))
|
723 |
+ 0.00086 * sin(2*(L1-2*L2 +W2))
|
724 |
- 0.00086 * sin(5*AnomSat - 2*AnomJup + 52.225/radcor)
|
725 |
- 0.00078 * sin(L2-L4)
|
726 |
- 0.00064 * sin(L1 - 2*L3 + Pi4)
|
727 |
- 0.00063 * sin(3*L3 - 7*L4 + 4*Pi4)
|
728 |
+ 0.00061 * sin(Pi1 - Pi4)
|
729 |
+ 0.00058 * sin(2*(NodeJup - LongPerJ - AnomJup))
|
730 |
+ 0.00058 * sin(W3 - W4)
|
731 |
+ 0.00056 * sin(2*(L2-L4))
|
732 |
+ 0.00055 * sin(2*(L1-L3))
|
733 |
+ 0.00052 * sin(3*L3 - 7*L4 + Pi3 + 3*Pi4)
|
734 |
- 0.00043 * sin(L1 - Pi3)
|
735 |
+ 0.00042 * sin(Pi3 - Pi2)
|
736 |
+ 0.00041 * sin(5*(L2-L3))
|
737 |
+ 0.00041 * sin(Pi4 - LongPerJ)
|
738 |
+ 0.00038 * sin(L2 - Pi1)
|
739 |
+ 0.00032 * sin(W2 - W3)
|
740 |
+ 0.00032 * sin(2*(L3 - AnomJup - LongPerJ))
|
741 |
+ 0.00029 * sin(Pi1 - Pi3);
|
742 |
|
743 |
S3 := 0.16477 * sin(L3 - Pi3)
|
744 |
+ 0.09062 * sin(L3 - Pi4)
|
745 |
- 0.06907 * sin(L2 - L3)
|
746 |
+ 0.03786 * sin(Pi3 - Pi4)
|
747 |
+ 0.01844 * sin(2*(L3-L4))
|
748 |
- 0.01340 * sin(AnomJup)
|
749 |
+ 0.00703 * sin(L2 - 2*L3 + Pi3)
|
750 |
- 0.00670 * sin(2*(NodeJup - LongPerJ))
|
751 |
- 0.00540 * sin(L3-L4)
|
752 |
+ 0.00481 * sin(Pi1 + Pi3 -2*LongPerJ - 2*AnomJup)
|
753 |
- 0.00409 * sin(L2 - 2*L3 + Pi2)
|
754 |
+ 0.00379 * sin(L2 - 2*L3 + Pi4)
|
755 |
+ 0.00235 * sin(NodeJup - W3)
|
756 |
+ 0.00198 * sin(NodeJup - W4)
|
757 |
+ 0.00180 * sin(PhiLambda)
|
758 |
+ 0.00129 * sin(3*(L3-L4))
|
759 |
+ 0.00124 * sin(L1-L3)
|
760 |
- 0.00119 * sin(5*AnomSat - 2*AnomJup + 52.225/radcor)
|
761 |
+ 0.00109 * sin(L1-L2)
|
762 |
- 0.00099 * sin(3*L3 - 7*L4 + 4*Pi4)
|
763 |
+ 0.00091 * sin(W3 - W4)
|
764 |
+ 0.00081 * sin(3*L3 - 7*L4 + Pi3 + 3*Pi4)
|
765 |
- 0.00076 * sin(2*L2 - 3*L3 + Pi3)
|
766 |
+ 0.00069 * sin(Pi4 - LongPerJ)
|
767 |
- 0.00058 * sin(2*L3 - 3*L4 + Pi4)
|
768 |
+ 0.00057 * sin(L3 + Pi3 - 2*LongPerJ - 2*AnomJup)
|
769 |
- 0.00057 * sin(L3 - 2*L4 + Pi4)
|
770 |
- 0.00052 * sin(Pi2 - Pi3)
|
771 |
- 0.00052 * sin(L2 - 2*L3 + Pi1)
|
772 |
+ 0.00048 * sin(L3 - 2*L4 + Pi3)
|
773 |
- 0.00045 * sin(2*L2 - 3*L3 + Pi4)
|
774 |
- 0.00041 * sin(Pi2 - Pi4)
|
775 |
- 0.00038 * sin(2*AnomJup)
|
776 |
- 0.00033 * sin(Pi3 - Pi4 + W3 - W4)
|
777 |
- 0.00032 * sin(3*L3 - 7*L4 + 2*Pi3 + 2*Pi4)
|
778 |
+ 0.00030 * sin(4*(L3-L4))
|
779 |
- 0.00029 * sin(W3 + NodeJup - 2*LongPerJ - 2*AnomJup)
|
780 |
+ 0.00029 * sin(L3 + Pi4 - 2*LongPerJ - 2*AnomJup)
|
781 |
+ 0.00026 * sin(L3 - LongPerJ - AnomJup)
|
782 |
+ 0.00024 * sin(L2 - 3*L3 + 2*L4)
|
783 |
+ 0.00021 * sin(2*(L3 - LongPerJ - AnomJup))
|
784 |
- 0.00021 * sin(L3 - Pi2)
|
785 |
+ 0.00017 * sin(2*(L3 - Pi3));
|
786 |
|
787 |
S4 := 0.84109 * sin(L4 - Pi4)
|
788 |
+ 0.03429 * sin(Pi4 - Pi3)
|
789 |
- 0.03305 * sin(2*(NodeJup - LongPerJ))
|
790 |
- 0.03211 * sin(AnomJup)
|
791 |
- 0.01860 * sin(L4 - Pi3)
|
792 |
+ 0.01182 * sin(NodeJup - W4)
|
793 |
+ 0.00622 * sin(L4 + Pi4 - 2*AnomJup - 2*LongPerJ)
|
794 |
+ 0.00385 * sin(2*(L4 - Pi4))
|
795 |
- 0.00284 * sin(5*AnomSat - 2*AnomJup + 52.225/radcor)
|
796 |
- 0.00233 * sin(2*(NodeJup - Pi4))
|
797 |
- 0.00223 * sin(L3 - L4)
|
798 |
- 0.00208 * sin(L4 - LongPerJ)
|
799 |
+ 0.00177 * sin(NodeJup + W4 - 2*Pi4)
|
800 |
+ 0.00134 * sin(Pi4 - LongPerJ)
|
801 |
+ 0.00125 * sin(2*(L4 - AnomJup - LongPerJ))
|
802 |
- 0.00117 * sin(2*AnomJup)
|
803 |
- 0.00112 * sin(2*(L3 - L4))
|
804 |
+ 0.00106 * sin(3*L3 - 7*L4 + 4*Pi4)
|
805 |
+ 0.00102 * sin(L4 - AnomJup - LongPerJ)
|
806 |
+ 0.00096 * sin(2*L4 - NodeJup - W4)
|
807 |
+ 0.00087 * sin(2*(NodeJup - W4))
|
808 |
- 0.00087 * sin(3*L3 - 7*L4 + Pi3 + 3*Pi4)
|
809 |
+ 0.00085 * sin(L3 - 2*L4 + Pi4)
|
810 |
- 0.00081 * sin(2*(L4 - NodeJup))
|
811 |
+ 0.00071 * sin(L4 + Pi4 -2*LongPerJ - 3*AnomJup)
|
812 |
+ 0.00060 * sin(L1 - L4)
|
813 |
- 0.00056 * sin(NodeJup - W3)
|
814 |
- 0.00055 * sin(L3 - 2*L4 + Pi3)
|
815 |
+ 0.00051 * sin(L2 - L4)
|
816 |
+ 0.00042 * sin(2*(NodeJup - AnomJup - LongPerJ))
|
817 |
+ 0.00039 * sin(2*(Pi4 - W4))
|
818 |
+ 0.00036 * sin(NodeJup + LongPerJ - Pi4 - W4)
|
819 |
+ 0.00035 * sin(2*AnomSat - AnomJup + 188.37/radcor)
|
820 |
- 0.00035 * sin(L4 - Pi4 + 2*LongPerJ - 2*NodeJup)
|
821 |
- 0.00032 * sin(L4 + Pi4 - 2*LongPerJ - AnomJup)
|
822 |
+ 0.00030 * sin(3*L3 - 7*L4 + 2*Pi3 + 2*Pi4)
|
823 |
+ 0.00030 * sin(2*AnomSat - 2*AnomJup + 149.15/radcor)
|
824 |
+ 0.00028 * sin(L4 - Pi4 + 2*NodeJup - 2*LongPerJ)
|
825 |
- 0.00028 * sin(2*(L4 - W4))
|
826 |
- 0.00027 * sin(Pi3 - Pi4 + W3 - W4)
|
827 |
- 0.00026 * sin(5*AnomSat - 3*AnomJup + 188.37/radcor)
|
828 |
+ 0.00025 * sin(W4 - W3)
|
829 |
- 0.00025 * sin(L2 - 3*L3 + 2*L4)
|
830 |
- 0.00023 * sin(3*(L3 - L4))
|
831 |
+ 0.00021 * sin(2*L4 - 2*LongPerJ - 3*AnomJup)
|
832 |
- 0.00021 * sin(2*L3 - 3*L4 + Pi4)
|
833 |
+ 0.00019 * sin(L4 - Pi4 - AnomJup)
|
834 |
- 0.00019 * sin(2*L4 - Pi3 - Pi4)
|
835 |
- 0.00018 * sin(L4 - Pi4 + AnomJup)
|
836 |
- 0.00016 * sin(L4 + Pi3 -2*LongPerJ - 2*AnomJup);
|
837 |
|
838 |
S1 := S1/radcor;
|
839 |
S2 := S2/radcor;
|
840 |
S3 := S3/radcor;
|
841 |
S4 := S4/radcor;
|
842 |
|
843 |
TL1 := L1 + S1;
|
844 |
TL2 := L2 + S2;
|
845 |
TL3 := L3 + S3;
|
846 |
TL4 := L4 + S4;
|
847 |
|
848 |
B1 := 0.0006502 * sin(TL1 - W1)
|
849 |
+ 0.0001835 * sin(TL1 - W2)
|
850 |
+ 0.0000329 * sin(TL1 - W3)
|
851 |
- 0.0000311 * sin(TL1 - NodeJup)
|
852 |
+ 0.0000093 * sin(TL1 - W4)
|
853 |
+ 0.0000075 * sin(3*TL1 - 4*L2 - 1.9927/radcor * S1 + W2)
|
854 |
+ 0.0000046 * sin(TL1 + NodeJup - 2*LongPerJ - 2*AnomJup);
|
855 |
B1 := ArcTan(B1);
|
856 |
|
857 |
B2 := 0.0081275 * sin(TL2 - W2)
|
858 |
+ 0.0004512 * sin(TL2 - W3)
|
859 |
- 0.0003286 * sin(TL2 - NodeJup)
|
860 |
+ 0.0001164 * sin(TL2 - W4)
|
861 |
+ 0.0000273 * sin(L1 - 2*L3 + 1.0146/radcor * S2 + W2)
|
862 |
+ 0.0000143 * sin(TL2 + NodeJup - 2*LongPerJ - 2*AnomJup)
|
863 |
- 0.0000143 * sin(TL2 - W1)
|
864 |
+ 0.0000035 * sin(TL2 - NodeJup + AnomJup)
|
865 |
- 0.0000028 * sin(L1 - 2*L3 + 1.0146/radcor * S2 + W3);
|
866 |
B2 := ArcTan(B2);
|
867 |
|
868 |
B3 := 0.0032364 * sin(TL3 - W3)
|
869 |
- 0.0016911 * sin(TL3 - NodeJup)
|
870 |
+ 0.0006849 * sin(TL3 - W4)
|
871 |
- 0.0002806 * sin(TL3 - W2)
|
872 |
+ 0.0000321 * sin(TL3 + NodeJup - 2*LongPerJ - 2*AnomJup)
|
873 |
+ 0.0000051 * sin(TL3 - NodeJup + AnomJup)
|
874 |
- 0.0000045 * sin(TL3 - NodeJup - AnomJup)
|
875 |
- 0.0000045 * sin(TL3 + NodeJup - 2*LongPerJ)
|
876 |
+ 0.0000037 * sin(TL3 + NodeJup - 2*LongPerJ - 3*AnomJup)
|
877 |
+ 0.0000030 * sin(2*L2 - 3*TL3 + 4.03/radcor * S3 + W2)
|
878 |
- 0.0000021 * sin(2*L2 - 3*TL3 + 4.03/radcor * S3 + W3);
|
879 |
B3 := ArcTan(B3);
|
880 |
|
881 |
B4 := -0.0076579 * sin(TL4 - NodeJup)
|
882 |
+ 0.0044148 * sin(TL4 - W4)
|
883 |
- 0.0005106 * sin(TL4 - W3)
|
884 |
+ 0.0000773 * sin(TL4 + NodeJup - 2*LongPerJ - 2*AnomJup)
|
885 |
+ 0.0000104 * sin(TL4 - NodeJup + AnomJup)
|
886 |
- 0.0000102 * sin(TL4 - NodeJup - AnomJup)
|
887 |
+ 0.0000088 * sin(TL4 + NodeJup - 2*LongPerJ - 3*AnomJup)
|
888 |
- 0.0000038 * sin(TL4 + NodeJup - 2*LongPerJ - AnomJup);
|
889 |
B4 := ArcTan(B4);
|
890 |
|
891 |
R1 := -0.0041339 * cos(2*(L1-L2))
|
892 |
- 0.0000395 * cos(L1 - Pi3)
|
893 |
- 0.0000214 * cos(L1 - Pi4)
|
894 |
+ 0.0000170 * cos(L1 - L2)
|
895 |
- 0.0000162 * cos(L1 - Pi1)
|
896 |
- 0.0000130 * cos(4*(L1-L2))
|
897 |
+ 0.0000106 * cos(L1 - L3)
|
898 |
- 0.0000063 * cos(L1 + Pi3 - 2*LongPerJ - 2*AnomJup);
|
899 |
|
900 |
R2 := 0.0093847 * cos(L1-L2)
|
901 |
- 0.0003114 * cos(L2 - Pi3)
|
902 |
- 0.0001738 * cos(L2 - Pi4)
|
903 |
- 0.0000941 * cos(L2 - Pi2)
|
904 |
+ 0.0000553 * cos(L2 - L3)
|
905 |
+ 0.0000523 * cos(L1 - L3)
|
906 |
- 0.0000290 * cos(2*(L1-L2))
|
907 |
+ 0.0000166 * cos(2*(L2-W2))
|
908 |
+ 0.0000107 * cos(L1 - 2*L3 + Pi3)
|
909 |
- 0.0000102 * cos(L2 - Pi1)
|
910 |
- 0.0000091 * cos(2*(L1-L3));
|
911 |
|
912 |
R3 := -0.0014377 * cos(L3 - Pi3)
|
913 |
- 0.0007904 * cos(L3 - Pi4)
|
914 |
+ 0.0006342 * cos(L2 - L3)
|
915 |
- 0.0001758 * cos(2*(L3-L4))
|
916 |
+ 0.0000294 * cos(L3 - L4)
|
917 |
- 0.0000156 * cos(3*(L3-L4))
|
918 |
+ 0.0000155 * cos(L1 - L3)
|
919 |
- 0.0000153 * cos(L1 - L2)
|
920 |
+ 0.0000070 * cos(2*L2 - 3*L3 + Pi3)
|
921 |
- 0.0000051 * cos(L3 + Pi3 - 2*LongPerJ - 2*AnomJup);
|
922 |
|
923 |
R4 := -0.0073391 * cos(L4 - Pi4)
|
924 |
+ 0.0001620 * cos(L4 - Pi3)
|
925 |
+ 0.0000974 * cos(L3 - L4)
|
926 |
- 0.0000541 * cos(L4 + Pi4 - 2*LongPerJ - 2*AnomJup)
|
927 |
- 0.0000269 * cos(2*(L4-Pi4))
|
928 |
+ 0.0000182 * cos(L4- LongPerJ)
|
929 |
+ 0.0000177 * cos(2*(L3-L4))
|
930 |
- 0.0000167 * cos(2*L4 - NodeJup - W4)
|
931 |
+ 0.0000167 * cos(NodeJup - W4)
|
932 |
- 0.0000155 * cos(2*(L4-LongPerj-AnomJup))
|
933 |
+ 0.0000142 * cos(2*(L4-NodeJup))
|
934 |
+ 0.0000104 * cos(L1 - L4)
|
935 |
+ 0.0000092 * cos(L2 - L4)
|
936 |
- 0.0000089 * cos(L4 - LongPerJ - AnomJup)
|
937 |
- 0.0000062 * cos(L4 + Pi4 - 2*LongPerJ - 3*AnomJup)
|
938 |
+ 0.0000048 * cos(2*(L4-W4));
|
939 |
|
940 |
R1 := 5.90730 * (1 + R1);
|
941 |
R2 := 9.39912 * (1 + R2);
|
942 |
R3 := 14.99240 * (1 + R3);
|
943 |
R4 := 26.36990 * (1 + R4);
|
944 |
|
945 |
T0 := (AJD - 2433282.423) / 36525;
|
946 |
Precession := (1.3966626*T0 + 0.0003088*sqr(T0)) / radcor;
|
947 |
|
948 |
TL1 := TL1 + Precession;
|
949 |
TL2 := TL2 + Precession;
|
950 |
TL3 := TL3 + Precession;
|
951 |
TL4 := TL4 + Precession;
|
952 |
NodeJup := NodeJup + Precession;
|
953 |
|
954 |
T0 := (AJD - AstJulianDatePrim(1900, 1, 1, 0)) / 36525;
|
955 |
Inclination := (3.120262 + 0.0006*T0) / radcor;
|
956 |
|
957 |
SatX[1] := R1 * cos(TL1 - NodeJup) * cos(B1);
|
958 |
SatY[1] := R1 * sin(TL1 - NodeJup) * cos(B1);
|
959 |
SatZ[1] := R1 * sin(B1);
|
960 |
|
961 |
SatX[2] := R2 * cos(TL2 - NodeJup) * cos(B2);
|
962 |
SatY[2] := R2 * sin(TL2 - NodeJup) * cos(B2);
|
963 |
SatZ[2] := R2 * sin(B2);
|
964 |
|
965 |
SatX[3] := R3 * cos(TL3 - NodeJup) * cos(B3);
|
966 |
SatY[3] := R3 * sin(TL3 - NodeJup) * cos(B3);
|
967 |
SatZ[3] := R3 * sin(B3);
|
968 |
|
969 |
SatX[4] := R4 * cos(TL4 - NodeJup) * cos(B4);
|
970 |
SatY[4] := R4 * sin(TL4 - NodeJup) * cos(B4);
|
971 |
SatZ[4] := R4 * sin(B4);
|
972 |
|
973 |
SatX[5] := 0;
|
974 |
SatY[5] := 0;
|
975 |
SatZ[5] := 1;
|
976 |
|
977 |
T0 := (AJD - 2451545.0) / 36525.0;
|
978 |
TD1 := 100.464441
|
979 |
+ 1.0209550 * T0
|
980 |
+ 0.00040117 * sqr(T0)
|
981 |
+ 0.000000569 * sqr(T0) * T0;
|
982 |
TD1 := TD1 / radcor;
|
983 |
|
984 |
TD2 := 1.303270
|
985 |
- 0.0054966 * T0
|
986 |
+ 0.00000465 * sqr(T0)
|
987 |
- 0.000000004 * sqr(T0) * T0;
|
988 |
TD2 := TD2 / radcor;
|
989 |
|
990 |
for I := 1 to 5 do begin
|
991 |
Transforms[I].A[1] := SatX[I];
|
992 |
Transforms[I].B[1] := SatY[I] * cos(Inclination)
|
993 |
- SatZ[I] * sin(Inclination);
|
994 |
Transforms[I].C[1] := SatY[I] * sin(Inclination)
|
995 |
+ SatZ[I] * cos(Inclination);
|
996 |
|
997 |
Transforms[I].A[2] := Transforms[I].A[1] * cos(NodeJup - TD1)
|
998 |
- Transforms[I].B[1] * sin(NodeJup - TD1);
|
999 |
Transforms[I].B[2] := Transforms[I].A[1] * sin(NodeJup - TD1)
|
1000 |
+ Transforms[I].B[1] * cos(NodeJup - TD1);
|
1001 |
Transforms[I].C[2] := Transforms[I].C[1];
|
1002 |
|
1003 |
Transforms[I].A[3] := Transforms[I].A[2];
|
1004 |
Transforms[I].B[3] := Transforms[I].B[2] * cos(TD2)
|
1005 |
- Transforms[I].C[2] * sin(TD2);
|
1006 |
Transforms[I].C[3] := Transforms[I].B[2] * sin(TD2)
|
1007 |
+ Transforms[I].C[2] * cos(TD2);
|
1008 |
|
1009 |
Transforms[I].A[4] := Transforms[I].A[3] * cos(TD1)
|
1010 |
- Transforms[I].B[3] * sin(TD1);
|
1011 |
Transforms[I].B[4] := Transforms[I].A[3] * sin(TD1)
|
1012 |
+ Transforms[I].B[3] * cos(TD1);
|
1013 |
Transforms[I].C[4] := Transforms[I].C[3];
|
1014 |
|
1015 |
Transforms[I].A[5] := Transforms[I].A[4] * sin(Jup1)
|
1016 |
- Transforms[I].B[4] * cos(Jup1);
|
1017 |
Transforms[I].B[5] := Transforms[I].A[4] * cos(Jup1)
|
1018 |
+ Transforms[I].B[4] * sin(Jup1);
|
1019 |
Transforms[I].C[5] := Transforms[I].C[4];
|
1020 |
|
1021 |
Transforms[I].A[6] := Transforms[I].A[5];
|
1022 |
Transforms[I].B[6] := Transforms[I].C[5] * sin(Jup2)
|
1023 |
+ Transforms[I].B[5] * cos(Jup2);
|
1024 |
Transforms[I].C[6] := Transforms[I].C[5] * cos(Jup2)
|
1025 |
- Transforms[I].B[5] * sin(Jup2);
|
1026 |
end;
|
1027 |
|
1028 |
Angle := StInvTan2(Transforms[5].C[6], Transforms[5].A[6]);
|
1029 |
|
1030 |
{Io calculations}
|
1031 |
Result.Io.X := Transforms[1].A[6] * cos(Angle)
|
1032 |
- Transforms[1].C[6] * sin(Angle);
|
1033 |
Result.Io.Y := Transforms[1].A[6] * sin(Angle)
|
1034 |
+ Transforms[1].C[6] * cos(Angle);
|
1035 |
TD1 := Transforms[1].B[6];
|
1036 |
|
1037 |
{correct for light time}
|
1038 |
TD2 := abs(TD1) / 17295 * sqrt(1 - sqr(Result.Io.X/R1));
|
1039 |
Result.Io.X := Result.Io.X + TD2;
|
1040 |
|
1041 |
{correct for perspective}
|
1042 |
TD2 := EJDist / (EJDist + TD1/2095);
|
1043 |
Result.Io.X := Result.Io.X * TD2;
|
1044 |
Result.Io.Y := Result.Io.Y * TD2;
|
1045 |
|
1046 |
{Europa calculations}
|
1047 |
Result.Europa.X := Transforms[2].A[6] * cos(Angle)
|
1048 |
- Transforms[2].C[6] * sin(Angle);
|
1049 |
Result.Europa.Y := Transforms[2].A[6] * sin(Angle)
|
1050 |
+ Transforms[2].C[6] * cos(Angle);
|
1051 |
TD1 := Transforms[2].B[6];
|
1052 |
|
1053 |
{correct for light time}
|
1054 |
TD2 := abs(TD1) / 21819 * sqrt(1 - sqr(Result.Europa.X/R2));
|
1055 |
Result.Europa.X := Result.Europa.X + TD2;
|
1056 |
|
1057 |
{correct for perspective}
|
1058 |
TD2 := EJDist / (EJDist + TD1/2095);
|
1059 |
Result.Europa.X := Result.Europa.X * TD2;
|
1060 |
Result.Europa.Y := Result.Europa.Y * TD2;
|
1061 |
|
1062 |
{Ganymede calculations}
|
1063 |
Result.Ganymede.X := Transforms[3].A[6] * cos(Angle)
|
1064 |
- Transforms[3].C[6] * sin(Angle);
|
1065 |
Result.Ganymede.Y := Transforms[3].A[6] * sin(Angle)
|
1066 |
+ Transforms[3].C[6] * cos(Angle);
|
1067 |
TD1 := Transforms[3].B[6];
|
1068 |
|
1069 |
{correct for light time}
|
1070 |
TD2 := abs(TD1) / 27558 * sqrt(1 - sqr(Result.Ganymede.X/R3));
|
1071 |
Result.Ganymede.X := Result.Ganymede.X + TD2;
|
1072 |
|
1073 |
{correct for perspective}
|
1074 |
TD2 := EJDist / (EJDist + TD1/2095);
|
1075 |
Result.Ganymede.X := Result.Ganymede.X * TD2;
|
1076 |
Result.Ganymede.Y := Result.Ganymede.Y * TD2;
|
1077 |
|
1078 |
{Callisto calculations}
|
1079 |
Result.Callisto.X := Transforms[4].A[6] * cos(Angle)
|
1080 |
- Transforms[4].C[6] * sin(Angle);
|
1081 |
Result.Callisto.Y := Transforms[4].A[6] * sin(Angle)
|
1082 |
+ Transforms[4].C[6] * cos(Angle);
|
1083 |
TD1 := Transforms[4].B[6];
|
1084 |
|
1085 |
{correct for light time}
|
1086 |
TD2 := abs(TD1) / 36548 * sqrt(1 - sqr(Result.Callisto.X/R4));
|
1087 |
Result.Callisto.X := Result.Callisto.X + TD2;
|
1088 |
|
1089 |
{correct for perspective}
|
1090 |
TD2 := EJDist / (EJDist + TD1/2095);
|
1091 |
Result.Callisto.X := Result.Callisto.X * TD2;
|
1092 |
Result.Callisto.Y := Result.Callisto.Y * TD2;
|
1093 |
end;
|
1094 |
|
1095 |
{-------------------------------------------------------------------------}
|
1096 |
|
1097 |
function GetJupSats(JD : TDateTime; HighPrecision, Shadows : Boolean) : TStJupSats;
|
1098 |
begin
|
1099 |
if not HighPrecision then
|
1100 |
Result := JupSatsLo(JD)
|
1101 |
else
|
1102 |
Result := JupSatsHi(JD, Shadows);
|
1103 |
end;
|
1104 |
|
1105 |
end.
|