spacecruft
/
celestia
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity

Corrected error in astro::anomaly() (and renamed from Anomaly in order be consistent with method name capitalization.)

ver1_5_1
Chris Laurel 2005-07-13 10:09:58 +00:00
parent 814497c85a
commit 2a0aa9cd3c
3 changed files with 23 additions and 23 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

12
src/celengine/astro.cpp
View File

@ -279,25 +279,25 @@ Point3d astro::equatorialToCelestialCart(double ra, double dec, double distance)
}
void astro::Anomaly(double meanAnomaly, double eccentricity,
void astro::anomaly(double meanAnomaly, double eccentricity,
double& trueAnomaly, double& eccentricAnomaly)
{
double e, delta, err;
double tol = 0.00000001745;
int iterations = 20; //limit while() to maximum of 20 iterations.
int iterations = 20; // limit while() to maximum of 20 iterations.
e = meanAnomaly - 2*PI * (int)(meanAnomaly/(2*PI));
e = meanAnomaly - 2*PI * (int) (meanAnomaly / (2*PI));
err = 1;
while(abs(err) > tol && iterations > 0)
{
err = e - eccentricity*sin(e) - meanAnomaly;
delta = err / (1 - eccentricity * sin(e));
delta = err / (1 - eccentricity * cos(e));
e -= delta;
iterations--;
}
trueAnomaly = 2*atan(sqrt((1+eccentricity)/(1-eccentricity))*tan(e/2));
eccentricAnomaly = e;
trueAnomaly = 2*atan(sqrt((1+eccentricity)/(1-eccentricity))*tan(e/2));
eccentricAnomaly = e;
}

2
src/celengine/astro.h
View File

@ -110,7 +110,7 @@ namespace astro
Point3f equatorialToCelestialCart(float ra, float dec, float distance);
Point3d equatorialToCelestialCart(double ra, double dec, double distance);
void Anomaly(double meanAnomaly, double eccentricity,
void anomaly(double meanAnomaly, double eccentricity,
double& trueAnomaly, double& eccentricAnomaly);
double meanEclipticObliquity(double jd);

32
src/celengine/customorbit.cpp
View File

@ -324,7 +324,7 @@ void computePlanetCoords(int p, double map, double da, double dhl, double dl,
s = gPlanetElements[p][3] + ds;
ma = map + dm;
astro::Anomaly(ma, s, nu, ea);
astro::anomaly(ma, s, nu, ea);
distance = (gPlanetElements[p][6] + da)*(1 - s*s)/(1 + s*cos(nu));
lp = radToDeg(nu) + gPlanetElements[p][2] + radToDeg(dml - dm);
lp = degToRad(lp);
@ -516,21 +516,21 @@ class EarthOrbit : public CachingOrbit
Point3d computePosition(double jd) const
{
double t, t2;
double ls, ms; /* mean longitude and mean anomoay */
double s, nu, ea; /* eccentricity, true anomaly, eccentric anomaly */
double ls, ms; // mean longitude and mean anomaly
double s, nu, ea; // eccentricity, true anomaly, eccentric anomaly
double a, b, a1, b1, c1, d1, e1, h1, dl, dr;
double eclLong, distance;
double eclLong, distance;
//Calculate the Julian centuries elapsed since 1900
// Calculate the Julian centuries elapsed since 1900
t = (jd - 2415020.0)/36525.0;
t2 = t*t;
a = 100.0021359*t;
b = 360.*(a-(int)a);
ls = 279.69668+.0003025*t2+b;
ms = meanAnomalySun(t);
ms = meanAnomalySun(t);
s = .016751-.0000418*t-1.26e-07*t2;
astro::Anomaly(degToRad(ms), s, nu, ea);
astro::anomaly(degToRad(ms), s, nu, ea);
a = 62.55209472000015*t;
b = 360*(a-(int)a);
a1 = degToRad(153.23+b);
@ -552,16 +552,16 @@ class EarthOrbit : public CachingOrbit
dr = 5.43e-06*sin(a1)+1.575e-05*sin(b1)+1.627e-05*sin(c1)+
3.076e-05*cos(d1)+9.27e-06*sin(h1);
eclLong = nu+degToRad(ls-ms+dl) + PI;
eclLong = pfmod(eclLong, TWOPI);
distance = KM_PER_AU * (1.0000002*(1-s*cos(ea))+dr);
eclLong = nu+degToRad(ls-ms+dl) + PI;
eclLong = pfmod(eclLong, TWOPI);
distance = KM_PER_AU * (1.0000002*(1-s*cos(ea))+dr);
//Correction for internal coordinate system
eclLong += PI;
return Point3d(-cos(eclLong) * distance,
0,
sin(eclLong) * distance);
// Correction for internal coordinate system
eclLong += PI;
return Point3d(-cos(eclLong) * distance,
0,
sin(eclLong) * distance);
};
double getPeriod() const