Rearranged FinalPosition

feature/19
Daniel Warner 2011-03-27 14:17:31 +01:00
parent 36712dd6c3
commit 7b51b2256b
1 changed files with 35 additions and 57 deletions

View File

@ -180,6 +180,11 @@ void SGDP4::Initialize(const double& theta2, const double& betao2, const double&
void SGDP4::FindPosition(double tsince) { void SGDP4::FindPosition(double tsince) {
double temp;
double temp1;
double temp2;
double temp3;
/* /*
* local copies which we can safely modify * local copies which we can safely modify
*/ */
@ -280,21 +285,13 @@ void SGDP4::FindPosition(double tsince) {
/* /*
* long period periodics * long period periodics
*/ */
double axn; double axn = tsince_eccentricity * cos(tsince_arg_perigee);
double xll; temp = 1.0 / (a * beta * beta);
double aynl; double xll = temp * xlcof_ * axn;
double xlt; double aynl = temp * aycof_;
double ayn; double xlt = xl + xll;
double elsq; double ayn = tsince_eccentricity * sin(tsince_arg_perigee) + aynl;
{ double elsq = axn * axn + ayn * ayn;
axn = tsince_eccentricity * cos(tsince_arg_perigee);
double temp1 = 1.0 / (a * beta * beta);
xll = temp1 * xlcof_ * axn;
aynl = temp1 * aycof_;
xlt = xl + xll;
ayn = tsince_eccentricity * sin(tsince_arg_perigee) + aynl;
elsq = axn * axn + ayn * ayn;
}
if (elsq >= 1.0) { if (elsq >= 1.0) {
throw new SatelliteException("Error: sqrt(e) >= 1 (elsq >= 1.0)"); throw new SatelliteException("Error: sqrt(e) >= 1 (elsq >= 1.0)");
@ -358,52 +355,33 @@ void SGDP4::FindPosition(double tsince) {
/* /*
* short period preliminary quantities * short period preliminary quantities
*/ */
double pl; temp = 1.0 - elsq;
double r; double pl = a * temp;
double betal; double r = a * (1.0 - ecose);
double rdot; temp1 = 1.0 / r;
double rfdot; double rdot = Globals::XKE() * sqrt(a) * esine * temp1;
double u; double rfdot = Globals::XKE() * sqrt(pl) * temp1;
double sin2u; temp2 = a * temp1;
double cos2u; double betal = sqrt(temp1);
{ temp3 = 1.0 / (1.0 + betal);
double elsq = axn * axn + ayn * ayn; double cosu = temp2 * (cosepw - axn + ayn * esine * temp3);
double temp1 = 1.0 - elsq; double sinu = temp2 * (sinepw - ayn - axn * esine * temp3);
pl = a * temp1; double u = atan2(sinu, cosu);
r = a * (1.0 - ecose); double sin2u = 2.0 * sinu * cosu;
betal = sqrt(temp1); double cos2u = 2.0 * cosu * cosu - 1.0;
double temp2 = 1.0 / r;
double temp3 = a * temp2;
double temp4 = 1.0 / (1.0 + betal);
rdot = Globals::XKE() * sqrt(a) * esine * temp2;
rfdot = Globals::XKE() * sqrt(pl) * temp2;
double cosu = temp3 * (cosepw - axn + ayn * esine * temp4);
double sinu = temp3 * (sinepw - ayn - axn * esine * temp4);
u = atan2(sinu, cosu);
sin2u = 2.0 * sinu * cosu;
cos2u = 2.0 * cosu * cosu - 1.0;
}
/* /*
* update for short periodics * update for short periodics
*/ */
double rk; temp = 1.0 / pl;
double uk; temp1 = Globals::CK2() * temp;
double xnodek; temp2 = temp1 * temp;
double xinck; double rk = r * (1.0 - 1.5 * temp2 * betal * x3thm1_) + 0.5 * temp1 * x1mth2_ * cos2u;
double rdotk; double uk = u - 0.25 * temp2 * x7thm1_ * sin2u;
double rfdotk; double xnodek = xnode + 1.5 * temp2 * cosio_ * sin2u;
{ double xinck = tsince_inclination + 1.5 * temp2 * cosio_ * sinio_ * cos2u;
double temp1 = 1.0 / pl; double rdotk = rdot - xn * temp1 * x1mth2_ * sin2u;
double temp2 = Globals::CK2() * temp1; double rfdotk = rfdot + xn * temp1 * (x1mth2_ * cos2u + 1.5 * x3thm1_);
double temp3 = temp2 * temp1;
rk = r * (1.0 - 1.5 * temp3 * betal * x3thm1_) + 0.5 * temp2 * x1mth2_ * cos2u;
uk = u - 0.25 * temp3 * x7thm1_ * sin2u;
xnodek = xnode + 1.5 * temp3 * cosio_ * sin2u;
xinck = tsince_inclination + 1.5 * temp3 * cosio_ * sinio_ * cos2u;
rdotk = rdot - xn * temp2 * x1mth2_ * sin2u;
rfdotk = rfdot + xn * temp2 * (x1mth2_ * cos2u + 1.5 * x3thm1_);
}
if (rk < 0.0) { if (rk < 0.0) {
throw new SatelliteException("Error: satellite decayed (rk < 0.0)"); throw new SatelliteException("Error: satellite decayed (rk < 0.0)");