1 /*
2 * - - - - - - - - - -
3 * g a l _ d s c o m
4 * - - - - - - - - - -
5 *
6 * This routine is part of the General Astrodynamics Library
7 *
8 * Description:
9 *
10 * This routine provides deep space common items used by both the secular
11 * and periodics routines.
12 *
13 * Status:
14 *
15 * SGP4 support routine.
16 *
17 * Given:
18 *
19 * epoch d
20 * ep d Eccentricity
21 * argpp d Argument of Perigee
22 * tc d
23 * inclp d Inclination
24 * nodep d Right Ascension of Ascending Node
25 * np d Mean Motion
26 *
27 * Returned:
28 *
29 * *sinim d
30 * *cosim d
31 * *sinomm d
32 * *cosomm d
33 * *snodm d
34 * *cnodm d
35 * *day d
36 * *e3 d
37 * *ee2 d
38 * *em d Eccentricity
39 * *emsq d Eccentricity squared
40 * *gam d
41 * *peo d
42 * *pgho d
43 * *pho d
44 * *pinco d
45 * *plo d
46 * *rtemsq d
47 * *se2 d
48 * *se3 d
49 * *sgh2 d
50 * *sgh3 d
51 * *sgh4 d
52 * *sh2 d
53 * *sh3 d
54 * *si2 d
55 * *si3 d
56 * *sl2 d
57 * *sl3 d
58 * *sl4 d
59 * *s1 d
60 * *s2 d
61 * *s3 d
62 * *s4 d
63 * *s5 d
64 * *s6 d
65 * *s7 d
66 * *ss1 d
67 * *ss2 d
68 * *ss3 d
69 * *ss4 d
70 * *ss5 d
71 * *ss6 d
72 * *ss7 d
73 * *sz1 d
74 * *sz2 d
75 * *sz3 d
76 * *sz11 d
77 * *sz12 d
78 * *sz13 d
79 * *sz21 d
80 * *sz22 d
81 * *sz23 d
82 * *sz31 d
83 * *sz32 d
84 * *sz33 d
85 * *xgh2 d
86 * *xgh3 d
87 * *xgh4 d
88 * *xh2 d
89 * *xh3 d
90 * *xi2 d
91 * *xi3 d
92 * *xl2 d
93 * *xl3 d
94 * *xl4 d
95 * *nm d Mean Motion
96 * *z1 d
97 * *z2 d
98 * *z3 d
99 * *z11 d
100 * *z12 d
101 * *z13 d
102 * *z21 d
103 * *z22 d
104 * *z23 d
105 * *z31 d
106 * *z32 d
107 * *z33 d
108 * *zmol d
109 * *zmos d
110 *
111 * Notes:
112 *
113 * 1) This routine is a translation from c++ to c of David Vallado's SGP4UNIT.dscom
114 * routine ( 2007 November 16 ). This routine used to be called dpper, but the
115 * functions inside weren't well organized.
116 *
117 * References:
118 *
119 * NORAD Spacetrack Report #3 1980
120 * Hoots, Roehrich
121 *
122 * NORAD Spacetrack Report #6 1986
123 * Hoots
124 *
125 * Hoots, Schumacher and Glover 2004
126 *
127 * Revisting Spacetrack Report #3
128 * Vallado, David, Crawford, Paul, Hujsak, Richard, Kelso, T.S.
129 * AIAA 2006-6753
130 *
131 * This revision:
132 *
133 * 2008 April 6
134 *
135 * Copyright (C) 2008 Paul C. L. Willmott. See notes at end.
136 *
137 *-----------------------------------------------------------------------
138 */
139
140 #include <math.h>
141 #include "gal_const.h"
142 #include "gal_dscom.h"
143
144 void
gal_dscom(double epoch,double ep,double argpp,double tc,double inclp,double nodep,double np,double * snodm,double * cnodm,double * sinim,double * cosim,double * sinomm,double * cosomm,double * day,double * e3,double * ee2,double * em,double * emsq,double * gam,double * peo,double * pgho,double * pho,double * pinco,double * plo,double * rtemsq,double * se2,double * se3,double * sgh2,double * sgh3,double * sgh4,double * sh2,double * sh3,double * si2,double * si3,double * sl2,double * sl3,double * sl4,double * s1,double * s2,double * s3,double * s4,double * s5,double * s6,double * s7,double * ss1,double * ss2,double * ss3,double * ss4,double * ss5,double * ss6,double * ss7,double * sz1,double * sz2,double * sz3,double * sz11,double * sz12,double * sz13,double * sz21,double * sz22,double * sz23,double * sz31,double * sz32,double * sz33,double * xgh2,double * xgh3,double * xgh4,double * xh2,double * xh3,double * xi2,double * xi3,double * xl2,double * xl3,double * xl4,double * nm,double * z1,double * z2,double * z3,double * z11,double * z12,double * z13,double * z21,double * z22,double * z23,double * z31,double * z32,double * z33,double * zmol,double * zmos)145 gal_dscom
146 (
147 double epoch,
148 double ep,
149 double argpp,
150 double tc,
151 double inclp,
152 double nodep,
153 double np,
154 double *snodm,
155 double *cnodm,
156 double *sinim,
157 double *cosim,
158 double *sinomm,
159 double *cosomm,
160 double *day,
161 double *e3,
162 double *ee2,
163 double *em,
164 double *emsq,
165 double *gam,
166 double *peo,
167 double *pgho,
168 double *pho,
169 double *pinco,
170 double *plo,
171 double *rtemsq,
172 double *se2,
173 double *se3,
174 double *sgh2,
175 double *sgh3,
176 double *sgh4,
177 double *sh2,
178 double *sh3,
179 double *si2,
180 double *si3,
181 double *sl2,
182 double *sl3,
183 double *sl4,
184 double *s1,
185 double *s2,
186 double *s3,
187 double *s4,
188 double *s5,
189 double *s6,
190 double *s7,
191 double *ss1,
192 double *ss2,
193 double *ss3,
194 double *ss4,
195 double *ss5,
196 double *ss6,
197 double *ss7,
198 double *sz1,
199 double *sz2,
200 double *sz3,
201 double *sz11,
202 double *sz12,
203 double *sz13,
204 double *sz21,
205 double *sz22,
206 double *sz23,
207 double *sz31,
208 double *sz32,
209 double *sz33,
210 double *xgh2,
211 double *xgh3,
212 double *xgh4,
213 double *xh2,
214 double *xh3,
215 double *xi2,
216 double *xi3,
217 double *xl2,
218 double *xl3,
219 double *xl4,
220 double *nm,
221 double *z1,
222 double *z2,
223 double *z3,
224 double *z11,
225 double *z12,
226 double *z13,
227 double *z21,
228 double *z22,
229 double *z23,
230 double *z31,
231 double *z32,
232 double *z33,
233 double *zmol,
234 double *zmos
235 )
236
237 {
238
239 /*
240 * Constants
241 */
242
243 const double zes = 0.01675 ;
244 const double zel = 0.05490 ;
245 const double c1ss = 2.9864797e-6 ;
246 const double c1l = 4.7968065e-7 ;
247 const double zsinis = 0.39785416 ;
248 const double zcosis = 0.91744867 ;
249 const double zcosgs = 0.1945905 ;
250 const double zsings = -0.98088458 ;
251
252 /*
253 * Local Variables
254 */
255
256 int lsflg ;
257
258 double a1 , a2 , a3 , a4 , a5 , a6 , a7 ,
259 a8 , a9 , a10 , betasq, cc , ctem , stem ,
260 x1 , x2 , x3 , x4 , x5 , x6 , x7 ,
261 x8 , xnodce, xnoi , zcosg , zcosgl, zcosh , zcoshl,
262 zcosi , zcosil, zsing , zsingl, zsinh , zsinhl, zsini ,
263 zsinil, zx , zy ;
264
265 /*
266 * - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
267 */
268
269 *nm = np ;
270 *em = ep ;
271 *snodm = sin ( nodep ) ;
272 *cnodm = cos ( nodep ) ;
273 *sinomm = sin ( argpp ) ;
274 *cosomm = cos ( argpp ) ;
275 *sinim = sin ( inclp ) ;
276 *cosim = cos ( inclp ) ;
277 *emsq = (*em) * (*em) ;
278 betasq = 1.0 - *emsq ;
279 *rtemsq = sqrt ( betasq ) ;
280
281 /*
282 * Initialize lunar solar terms
283 */
284
285 *peo = 0.0 ;
286 *pinco = 0.0 ;
287 *plo = 0.0 ;
288 *pgho = 0.0 ;
289 *pho = 0.0 ;
290 *day = epoch + 18261.5 + tc / 1440.0 ;
291 xnodce = fmod ( 4.5236020 - 9.2422029e-4 * (*day), GAL_2PI ) ;
292 stem = sin ( xnodce ) ;
293 ctem = cos ( xnodce ) ;
294 zcosil = 0.91375164 - 0.03568096 * ctem ;
295 zsinil = sqrt ( 1.0 - zcosil * zcosil ) ;
296 zsinhl = 0.089683511 * stem / zsinil ;
297 zcoshl = sqrt ( 1.0 - zsinhl * zsinhl ) ;
298 *gam = 5.8351514 + 0.0019443680 * (*day) ;
299 zx = 0.39785416 * stem / zsinil ;
300 zy = zcoshl * ctem + 0.91744867 * zsinhl * stem ;
301 zx = atan2 ( zx, zy ) ;
302 zx = *gam + zx - xnodce ;
303 zcosgl = cos ( zx ) ;
304 zsingl = sin ( zx ) ;
305
306 /*
307 * Do solar terms
308 */
309
310 zcosg = zcosgs ;
311 zsing = zsings ;
312 zcosi = zcosis ;
313 zsini = zsinis ;
314 zcosh = *cnodm ;
315 zsinh = *snodm ;
316 cc = c1ss ;
317 xnoi = 1.0 / *nm ;
318
319 for ( lsflg = 1; lsflg <= 2; lsflg++ ) {
320 a1 = zcosg * zcosh + zsing * zcosi * zsinh ;
321 a3 = -zsing * zcosh + zcosg * zcosi * zsinh ;
322 a7 = -zcosg * zsinh + zsing * zcosi * zcosh ;
323 a8 = zsing * zsini ;
324 a9 = zsing * zsinh + zcosg * zcosi * zcosh ;
325 a10 = zcosg * zsini ;
326 a2 = *cosim * a7 + *sinim * a8 ;
327 a4 = *cosim * a9 + *sinim * a10 ;
328 a5 = -*sinim * a7 + *cosim * a8 ;
329 a6 = -*sinim * a9 + *cosim * a10 ;
330
331 x1 = a1 * (*cosomm) + a2 * (*sinomm) ;
332 x2 = a3 * (*cosomm) + a4 * (*sinomm) ;
333 x3 = -a1 * (*sinomm) + a2 * (*cosomm) ;
334 x4 = -a3 * (*sinomm) + a4 * (*cosomm) ;
335 x5 = a5 * (*sinomm) ;
336 x6 = a6 * (*sinomm) ;
337 x7 = a5 * (*cosomm) ;
338 x8 = a6 * (*cosomm) ;
339
340 *z31 = 12.0 * x1 * x1 - 3.0 * x3 * x3 ;
341 *z32 = 24.0 * x1 * x2 - 6.0 * x3 * x4 ;
342 *z33 = 12.0 * x2 * x2 - 3.0 * x4 * x4 ;
343 *z1 = 3.0 * ( a1 * a1 + a2 * a2 ) + *z31 * (*emsq) ;
344 *z2 = 6.0 * ( a1 * a3 + a2 * a4 ) + *z32 * (*emsq) ;
345 *z3 = 3.0 * ( a3 * a3 + a4 * a4 ) + *z33 * (*emsq) ;
346 *z11 = -6.0 * a1 * a5 + (*emsq) * ( -24.0 * x1 * x7-6.0 * x3 * x5 ) ;
347 *z12 = -6.0 * ( a1 * a6 + a3 * a5 ) + (*emsq) *
348 ( -24.0 * (x2 * x7 + x1 * x8 ) - 6.0 * ( x3 * x6 + x4 * x5 ) ) ;
349 *z13 = -6.0 * a3 * a6 + (*emsq) * ( -24.0 * x2 * x8 - 6.0 * x4 * x6 ) ;
350 *z21 = 6.0 * a2 * a5 + (*emsq) * ( 24.0 * x1 * x5 - 6.0 * x3 * x7 ) ;
351 *z22 = 6.0 * ( a4 * a5 + a2 * a6) + (*emsq) *
352 ( 24.0 * ( x2 * x5 + x1 * x6 ) - 6.0 * ( x4 * x7 + x3 * x8 ) ) ;
353 *z23 = 6.0 * a4 * a6 + (*emsq) * ( 24.0 * x2 * x6 - 6.0 * x4 * x8 ) ;
354 *z1 = *z1 + *z1 + betasq * (*z31) ;
355 *z2 = *z2 + *z2 + betasq * (*z32) ;
356 *z3 = *z3 + *z3 + betasq * (*z33) ;
357 *s3 = cc * xnoi ;
358 *s2 = -0.5 * (*s3) / (*rtemsq) ;
359 *s4 = (*s3) * (*rtemsq) ;
360 *s1 = -15.0 * (*em) * (*s4) ;
361 *s5 = x1 * x3 + x2 * x4 ;
362 *s6 = x2 * x3 + x1 * x4 ;
363 *s7 = x2 * x4 - x1 * x3 ;
364
365 /*
366 * Do lunar terms -
367 */
368
369 if (lsflg == 1) {
370 *ss1 = *s1 ;
371 *ss2 = *s2 ;
372 *ss3 = *s3 ;
373 *ss4 = *s4 ;
374 *ss5 = *s5 ;
375 *ss6 = *s6 ;
376 *ss7 = *s7 ;
377 *sz1 = *z1 ;
378 *sz2 = *z2 ;
379 *sz3 = *z3 ;
380 *sz11 = *z11 ;
381 *sz12 = *z12 ;
382 *sz13 = *z13 ;
383 *sz21 = *z21 ;
384 *sz22 = *z22 ;
385 *sz23 = *z23 ;
386 *sz31 = *z31 ;
387 *sz32 = *z32 ;
388 *sz33 = *z33 ;
389 zcosg = zcosgl ;
390 zsing = zsingl ;
391 zcosi = zcosil ;
392 zsini = zsinil ;
393 zcosh = zcoshl * (*cnodm) + zsinhl * (*snodm) ;
394 zsinh = (*snodm) * zcoshl - (*cnodm) * zsinhl ;
395 cc = c1l ;
396 }
397 }
398
399 *zmol = fmod ( 4.7199672 + 0.22997150 * (*day) - (*gam), GAL_2PI ) ;
400 *zmos = fmod ( 6.2565837 + 0.017201977 * (*day), GAL_2PI ) ;
401
402 /*
403 * Do solar terms
404 */
405
406 *se2 = 2.0 * (*ss1) * (*ss6) ;
407 *se3 = 2.0 * (*ss1) * (*ss7) ;
408 *si2 = 2.0 * (*ss2) * (*sz12) ;
409 *si3 = 2.0 * (*ss2) * ( *sz13 - *sz11 ) ;
410 *sl2 = -2.0 * (*ss3) * (*sz2) ;
411 *sl3 = -2.0 * (*ss3) * ( *sz3 - *sz1 ) ;
412 *sl4 = -2.0 * (*ss3) * ( -21.0 - 9.0 * (*emsq) ) * zes ;
413 *sgh2 = 2.0 * (*ss4) * (*sz32) ;
414 *sgh3 = 2.0 * (*ss4) * ( *sz33 - *sz31 ) ;
415 *sgh4 = -18.0 * (*ss4) * zes ;
416 *sh2 = -2.0 * (*ss2) * (*sz22) ;
417 *sh3 = -2.0 * (*ss2) * ( *sz23 - *sz21 ) ;
418
419 /*
420 * Do lunar terms
421 */
422
423 *ee2 = 2.0 * (*s1) * (*s6) ;
424 *e3 = 2.0 * (*s1) * (*s7) ;
425 *xi2 = 2.0 * (*s2) * (*z12) ;
426 *xi3 = 2.0 * (*s2) * ( *z13 - *z11 ) ;
427 *xl2 = -2.0 * (*s3) * (*z2) ;
428 *xl3 = -2.0 * (*s3) * ( *z3 - *z1 ) ;
429 *xl4 = -2.0 * (*s3) * ( -21.0 - 9.0 * (*emsq) ) * zel ;
430 *xgh2 = 2.0 * (*s4) * (*z32) ;
431 *xgh3 = 2.0 * (*s4) * ( *z33 - *z31 ) ;
432 *xgh4 = -18.0 * (*s4) * zel ;
433 *xh2 = -2.0 * (*s2) * (*z22) ;
434 *xh3 = -2.0 * (*s2) * ( *z23 - *z21 ) ;
435
436 /*
437 * Finished.
438 */
439
440 }
441
442 /*
443 * gal - General Astrodynamics Library
444 * Copyright (C) 2008 Paul C. L. Willmott
445 *
446 * This program is free software; you can redistribute it and/or modify
447 * it under the terms of the GNU General Public License as published by
448 * the Free Software Foundation; either version 2 of the License, or
449 * (at your option) any later version.
450 *
451 * This program is distributed in the hope that it will be useful,
452 * but WITHOUT ANY WARRANTY; without even the implied warranty of
453 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
454 * GNU General Public License for more details.
455 *
456 * You should have received a copy of the GNU General Public License along
457 * with this program; if not, write to the Free Software Foundation, Inc.,
458 * 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
459 *
460 * Contact:
461 *
462 * Paul Willmott
463 * vp9mu@amsat.org
464 */
465
466