sgp4unit.cpp   sgp4unit.cpp 
skipping to change at line 50 skipping to change at line 50
* what was changed * what was changed
* 10 aug 04 david vallado * 10 aug 04 david vallado
* 2nd printing baseline working * 2nd printing baseline working
* 14 may 01 david vallado * 14 may 01 david vallado
* 2nd edition baseline * 2nd edition baseline
* 80 norad * 80 norad
* original baseline * original baseline
* ---------------------------------------------------------------- */ * ---------------------------------------------------------------- */
#include "sgp4unit.h" #include "sgp4unit.h"
#include <cmath>
FILE *dbgfile; FILE *dbgfile;
/* ----------- local functions - only ever used internally by sgp4 -------- -- */ /* ----------- local functions - only ever used internally by sgp4 -------- -- */
static void dpper static void dpper
( (
double e3, double ee2, double peo, double pgho, double pho, double e3, double ee2, double peo, double pgho, double pho,
double pinco, double plo, double se2, double se3, double sgh2, double pinco, double plo, double se2, double se3, double sgh2,
double sgh3, double sgh4, double sh2, double sh3, double si2, double sgh3, double sgh4, double sh2, double sh3, double si2,
double si3, double sl2, double sl3, double sl4, double t, double si3, double sl2, double sl3, double sl4, double t,
skipping to change at line 455 skipping to change at line 456
nm = np; nm = np;
em = ep; em = ep;
snodm = sin(nodep); snodm = sin(nodep);
cnodm = cos(nodep); cnodm = cos(nodep);
sinomm = sin(argpp); sinomm = sin(argpp);
cosomm = cos(argpp); cosomm = cos(argpp);
sinim = sin(inclp); sinim = sin(inclp);
cosim = cos(inclp); cosim = cos(inclp);
emsq = em * em; emsq = em * em;
betasq = 1.0 - emsq; betasq = 1.0 - emsq;
rtemsq = sqrt(betasq); rtemsq = std::sqrt(betasq);
/* ----------------- initialize lunar solar terms --------------- */ /* ----------------- initialize lunar solar terms --------------- */
peo = 0.0; peo = 0.0;
pinco = 0.0; pinco = 0.0;
plo = 0.0; plo = 0.0;
pgho = 0.0; pgho = 0.0;
pho = 0.0; pho = 0.0;
day = epoch + 18261.5 + tc / 1440.0; day = epoch + 18261.5 + tc / 1440.0;
xnodce = fmod(4.5236020 - 9.2422029e-4 * day, twopi); xnodce = fmod(4.5236020 - 9.2422029e-4 * day, twopi);
stem = sin(xnodce); stem = sin(xnodce);
ctem = cos(xnodce); ctem = cos(xnodce);
zcosil = 0.91375164 - 0.03568096 * ctem; zcosil = 0.91375164 - 0.03568096 * ctem;
zsinil = sqrt(1.0 - zcosil * zcosil); zsinil = std::sqrt(1.0 - zcosil * zcosil);
zsinhl = 0.089683511 * stem / zsinil; zsinhl = 0.089683511 * stem / zsinil;
zcoshl = sqrt(1.0 - zsinhl * zsinhl); zcoshl = std::sqrt(1.0 - zsinhl * zsinhl);
gam = 5.8351514 + 0.0019443680 * day; gam = 5.8351514 + 0.0019443680 * day;
zx = 0.39785416 * stem / zsinil; zx = 0.39785416 * stem / zsinil;
zy = zcoshl * ctem + 0.91744867 * zsinhl * stem; zy = zcoshl * ctem + 0.91744867 * zsinhl * stem;
zx = atan2(zx, zy); zx = atan2(zx, zy);
zx = gam + zx - xnodce; zx = gam + zx - xnodce;
zcosgl = cos(zx); zcosgl = cos(zx);
zsingl = sin(zx); zsingl = sin(zx);
/* ------------------------- do solar terms --------------------- */ /* ------------------------- do solar terms --------------------- */
zcosg = zcosgs; zcosg = zcosgs;
skipping to change at line 1212 skipping to change at line 1213
const double twopi = 2.0 * M_PI; const double twopi = 2.0 * M_PI;
/* ----------------------- earth constants ---------------------- */ /* ----------------------- earth constants ---------------------- */
// sgp4fix identify constants and allow alternate values // sgp4fix identify constants and allow alternate values
getgravconst( whichconst, tumin, mu, radiusearthkm, xke, j2, j3, j4, j 3oj2 ); getgravconst( whichconst, tumin, mu, radiusearthkm, xke, j2, j3, j4, j 3oj2 );
x2o3 = 2.0 / 3.0; x2o3 = 2.0 / 3.0;
/* ------------- calculate auxillary epoch quantities ---------- */ /* ------------- calculate auxillary epoch quantities ---------- */
eccsq = ecco * ecco; eccsq = ecco * ecco;
omeosq = 1.0 - eccsq; omeosq = 1.0 - eccsq;
rteosq = sqrt(omeosq); rteosq = std::sqrt(omeosq);
cosio = cos(inclo); cosio = cos(inclo);
cosio2 = cosio * cosio; cosio2 = cosio * cosio;
/* ------------------ un-kozai the mean motion ----------------- */ /* ------------------ un-kozai the mean motion ----------------- */
ak = pow(xke / no, x2o3); ak = pow(xke / no, x2o3);
d1 = 0.75 * j2 * (3.0 * cosio2 - 1.0) / (rteosq * omeosq); d1 = 0.75 * j2 * (3.0 * cosio2 - 1.0) / (rteosq * omeosq);
del = d1 / (ak * ak); del = d1 / (ak * ak);
adel = ak * (1.0 - del * del - del * adel = ak * (1.0 - del * del - del *
(1.0 / 3.0 + 134.0 * del * del / 81.0)); (1.0 / 3.0 + 134.0 * del * del / 81.0));
del = d1/(adel * adel); del = d1/(adel * adel);
skipping to change at line 1704 skipping to change at line 1705
* vallado, crawford, hujsak, kelso 2006 * vallado, crawford, hujsak, kelso 2006
------------------------------------------------------------------------- ---*/ ------------------------------------------------------------------------- ---*/
bool sgp4 bool sgp4
( (
gravconsttype whichconst, elsetrec& satrec, double tsince, gravconsttype whichconst, elsetrec& satrec, double tsince,
double r[3], double v[3] double r[3], double v[3]
) )
{ {
double am , axnl , aynl , betal , cosim , cnod , double am , axnl , aynl , betal , cosim , cnod ,
cos2u, coseo1, cosi , cosip , cosisq, cossu , cosu, cos2u, coseo1=0.0, cosi , cosip , cosisq, cossu , cosu,
delm , delomg, em , emsq , ecose , el2 , eo1 , delm , delomg, em , emsq , ecose , el2 , eo1 ,
ep , esine , argpm, argpp , argpdf, pl, mrt = 0.0, ep , esine , argpm, argpp , argpdf, pl, mrt = 0.0,
mvt , rdotl , rl , rvdot , rvdotl, sinim , mvt , rdotl , rl , rvdot , rvdotl, sinim ,
sin2u, sineo1, sini , sinip , sinsu , sinu , sin2u, sineo1=0.0, sini , sinip , sinsu , sinu ,
snod , su , t2 , t3 , t4 , tem5 , temp, snod , su , t2 , t3 , t4 , tem5 , temp,
temp1, temp2 , tempa, tempe , templ , u , ux , temp1, temp2 , tempa, tempe , templ , u , ux ,
uy , uz , vx , vy , vz , inclm , mm , uy , uz , vx , vy , vz , inclm , mm ,
nm , nodem, xinc , xincp , xl , xlm , mp , nm , nodem, xinc , xincp , xl , xlm , mp ,
xmdf , xmx , xmy , nodedf, xnode , nodep, tc , dndt, xmdf , xmx , xmy , nodedf, xnode , nodep, tc , dndt,
twopi, x2o3 , j2=0 , j3 , tumin, j4=0 , xke=0 , j3oj2=0, r adiusearthkm=0, twopi, x2o3 , j2=0 , j3 , tumin, j4=0 , xke=0 , j3oj2=0, r adiusearthkm=0,
mu, vkmpersec; mu, vkmpersec;
int ktr; int ktr;
/* ------------------ set mathematical constants --------------- */ /* ------------------ set mathematical constants --------------- */
skipping to change at line 1915 skipping to change at line 1916
if (pl < 0.0) if (pl < 0.0)
{ {
// printf("# error pl %f\n", pl); // printf("# error pl %f\n", pl);
satrec.error = 4; satrec.error = 4;
// sgp4fix add return // sgp4fix add return
return false; return false;
} }
else else
{ {
rl = am * (1.0 - ecose); rl = am * (1.0 - ecose);
rdotl = sqrt(am) * esine/rl; rdotl = std::sqrt(am) * esine/rl;
rvdotl = sqrt(pl) / rl; rvdotl = std::sqrt(pl) / rl;
betal = sqrt(1.0 - el2); betal = std::sqrt(1.0 - el2);
temp = esine / (1.0 + betal); temp = esine / (1.0 + betal);
sinu = am / rl * (sineo1 - aynl - axnl * temp); sinu = am / rl * (sineo1 - aynl - axnl * temp);
cosu = am / rl * (coseo1 - axnl + aynl * temp); cosu = am / rl * (coseo1 - axnl + aynl * temp);
su = atan2(sinu, cosu); su = atan2(sinu, cosu);
sin2u = (cosu + cosu) * sinu; sin2u = (cosu + cosu) * sinu;
cos2u = 1.0 - 2.0 * sinu * sinu; cos2u = 1.0 - 2.0 * sinu * sinu;
temp = 1.0 / pl; temp = 1.0 / pl;
temp1 = 0.5 * j2 * temp; temp1 = 0.5 * j2 * temp;
temp2 = temp1 * temp; temp2 = temp1 * temp;
skipping to change at line 2090 skipping to change at line 2091
tumin = 1.0 / xke; tumin = 1.0 / xke;
j2 = 0.001082616; j2 = 0.001082616;
j3 = -0.00000253881; j3 = -0.00000253881;
j4 = -0.00000165597; j4 = -0.00000165597;
j3oj2 = j3 / j2; j3oj2 = j3 / j2;
break; break;
// ------------ wgs-72 constants ------------ // ------------ wgs-72 constants ------------
case wgs72: case wgs72:
mu = 398600.8; // in km3 / s2 mu = 398600.8; // in km3 / s2
radiusearthkm = 6378.135; // km radiusearthkm = 6378.135; // km
xke = 60.0 / sqrt(radiusearthkm*radiusearthkm*radiusearthkm/m u); xke = 60.0 / std::sqrt(radiusearthkm*radiusearthkm*radiusearth km/mu);
tumin = 1.0 / xke; tumin = 1.0 / xke;
j2 = 0.001082616; j2 = 0.001082616;
j3 = -0.00000253881; j3 = -0.00000253881;
j4 = -0.00000165597; j4 = -0.00000165597;
j3oj2 = j3 / j2; j3oj2 = j3 / j2;
break; break;
case wgs84: case wgs84:
// ------------ wgs-84 constants ------------ // ------------ wgs-84 constants ------------
mu = 398600.5; // in km3 / s2 mu = 398600.5; // in km3 / s2
radiusearthkm = 6378.137; // km radiusearthkm = 6378.137; // km
xke = 60.0 / sqrt(radiusearthkm*radiusearthkm*radiusearthkm/m u); xke = 60.0 / std::sqrt(radiusearthkm*radiusearthkm*radiusearth km/mu);
tumin = 1.0 / xke; tumin = 1.0 / xke;
j2 = 0.00108262998905; j2 = 0.00108262998905;
j3 = -0.00000253215306; j3 = -0.00000253215306;
j4 = -0.00000161098761; j4 = -0.00000161098761;
j3oj2 = j3 / j2; j3oj2 = j3 / j2;
break; break;
default: default:
fprintf(stderr,"unknown gravity option (%d)\n",whichconst); fprintf(stderr,"unknown gravity option (%d)\n",whichconst);
break; break;
} }
 End of changes. 10 change blocks. 
11 lines changed or deleted 12 lines changed or added

This html diff was produced by rfcdiff 1.41. The latest version is available from http://tools.ietf.org/tools/rfcdiff/