Andy K / Mbed 2 deprecated SOWB

Satellite Observers Workbench. NOT yet complete, just published for forum posters to \"cherry pick\" pieces of code as requiered as an example.

Dependencies:   mbed

sgp4sdp4/solar.c@0:0a841b89d614, 2010-10-11 (annotated)

Committer:
AjK
Date:
Mon Oct 11 10:34:55 2010 +0000
Revision:
0:0a841b89d614
Totally Alpha quality as this project isn\t completed. Just publishing it as it answers many questions asked in the forums

Who changed what in which revision?

UserRevisionLine numberNew contents of line
AjK 0:0a841b89d614 1 /*
AjK 0:0a841b89d614 2 * Unit Solar
AjK 0:0a841b89d614 3 * Author: Dr TS Kelso
AjK 0:0a841b89d614 4 * Original Version: 1990 Jul 29
AjK 0:0a841b89d614 5 * Current Revision: 1999 Nov 27
AjK 0:0a841b89d614 6 * Version: 1.30
AjK 0:0a841b89d614 7 * Copyright: 1990-1999, All Rights Reserved
AjK 0:0a841b89d614 8 *
AjK 0:0a841b89d614 9 * Ported to C by: Neoklis Kyriazis April 1 2001
AjK 0:0a841b89d614 10 */
AjK 0:0a841b89d614 11
AjK 0:0a841b89d614 12 #include "sgp4sdp4.h"
AjK 0:0a841b89d614 13
AjK 0:0a841b89d614 14 /* Calculates solar position vector */
AjK 0:0a841b89d614 15 void
AjK 0:0a841b89d614 16 Calculate_Solar_Position(double _time, vector_t *solar_vector)
AjK 0:0a841b89d614 17 {
AjK 0:0a841b89d614 18 double mjd,year,T,M,L,e,C,O,Lsa,nu,R,eps;
AjK 0:0a841b89d614 19
AjK 0:0a841b89d614 20 mjd = _time - 2415020.0;
AjK 0:0a841b89d614 21 year = 1900 + mjd/365.25;
AjK 0:0a841b89d614 22 T = (mjd + Delta_ET(year)/secday)/36525.0;
AjK 0:0a841b89d614 23 M = Radians(Modulus(358.47583 + Modulus(35999.04975*T,360.0)
AjK 0:0a841b89d614 24 - (0.000150 + 0.0000033*T)*Sqr(T),360.0));
AjK 0:0a841b89d614 25 L = Radians(Modulus(279.69668 + Modulus(36000.76892*T,360.0)
AjK 0:0a841b89d614 26 + 0.0003025*Sqr(T),360.0));
AjK 0:0a841b89d614 27 e = 0.01675104 - (0.0000418 + 0.000000126*T)*T;
AjK 0:0a841b89d614 28 C = Radians((1.919460 - (0.004789 + 0.000014*T)*T)*sin(M)
AjK 0:0a841b89d614 29 + (0.020094 - 0.000100*T)*sin(2*M) + 0.000293*sin(3*M));
AjK 0:0a841b89d614 30 O = Radians(Modulus(259.18 - 1934.142*T,360.0));
AjK 0:0a841b89d614 31 Lsa = Modulus(L + C - Radians(0.00569 - 0.00479*sin(O)),twopi);
AjK 0:0a841b89d614 32 nu = Modulus(M + C,twopi);
AjK 0:0a841b89d614 33 R = 1.0000002*(1 - Sqr(e))/(1 + e*cos(nu));
AjK 0:0a841b89d614 34 eps = Radians(23.452294 - (0.0130125 + (0.00000164 -
AjK 0:0a841b89d614 35 0.000000503*T)*T)*T + 0.00256*cos(O));
AjK 0:0a841b89d614 36 R = SGPAU*R;
AjK 0:0a841b89d614 37 solar_vector->x = R*cos(Lsa);
AjK 0:0a841b89d614 38 solar_vector->y = R*sin(Lsa)*cos(eps);
AjK 0:0a841b89d614 39 solar_vector->z = R*sin(Lsa)*sin(eps);
AjK 0:0a841b89d614 40 solar_vector->w = R;
AjK 0:0a841b89d614 41 } /*Procedure Calculate_Solar_Position*/
AjK 0:0a841b89d614 42
AjK 0:0a841b89d614 43 /*------------------------------------------------------------------*/
AjK 0:0a841b89d614 44
AjK 0:0a841b89d614 45 /* Calculates stellite's eclipse status and depth */
AjK 0:0a841b89d614 46 int
AjK 0:0a841b89d614 47 Sat_Eclipsed(vector_t *pos, vector_t *sol, double *depth)
AjK 0:0a841b89d614 48 {
AjK 0:0a841b89d614 49 double sd_sun, sd_earth, delta;
AjK 0:0a841b89d614 50 vector_t Rho, earth;
AjK 0:0a841b89d614 51
AjK 0:0a841b89d614 52 /* Determine partial eclipse */
AjK 0:0a841b89d614 53 sd_earth = ArcSin(xkmper/pos->w);
AjK 0:0a841b89d614 54 Vec_Sub(sol,pos,&Rho);
AjK 0:0a841b89d614 55 sd_sun = ArcSin(__sr__/Rho.w);
AjK 0:0a841b89d614 56 Scalar_Multiply(-1,pos,&earth);
AjK 0:0a841b89d614 57 delta = Angle(sol,&earth);
AjK 0:0a841b89d614 58 *depth = sd_earth - sd_sun - delta;
AjK 0:0a841b89d614 59 if( sd_earth < sd_sun )
AjK 0:0a841b89d614 60 return( 0 );
AjK 0:0a841b89d614 61 else
AjK 0:0a841b89d614 62 if( *depth >= 0 )
AjK 0:0a841b89d614 63 return( 1 );
AjK 0:0a841b89d614 64 else
AjK 0:0a841b89d614 65 return( 0 );
AjK 0:0a841b89d614 66
AjK 0:0a841b89d614 67 } /*Function Sat_Eclipsed*/
AjK 0:0a841b89d614 68
AjK 0:0a841b89d614 69 /*------------------------------------------------------------------*/