From edcd42bc2d587722a6bc5b0d53eb28ace2aa8b69 Mon Sep 17 00:00:00 2001 From: Edward d'Auvergne Date: Mon, 7 Dec 2015 21:03:33 +0100 Subject: [PATCH] Optimisation of the celestialBody ephemeris code. By storing repetitive intermediate calculations, the number of mathematical operations for a single call to CelestialBody::updatePosition() has decreased by 12. This matches the changes to MoonPos::updatePosition(). --- simgear/ephemeris/celestialBody.cxx | 22 ++++++++++++++++------ 1 file changed, 16 insertions(+), 6 deletions(-) diff --git a/simgear/ephemeris/celestialBody.cxx b/simgear/ephemeris/celestialBody.cxx index 4e887392..ad4a76f2 100644 --- a/simgear/ephemeris/celestialBody.cxx +++ b/simgear/ephemeris/celestialBody.cxx @@ -52,7 +52,8 @@ void CelestialBody::updatePosition(double mjd, Star *ourSun) { double eccAnom, v, ecl, actTime, - xv, yv, xh, yh, zh, xg, yg, zg, xe, ye, ze; + xv, yv, xh, yh, zh, xg, yg, zg, xe, ye, ze, + cosN, sinN, cosvw, sinvw, sinvw_cosi, cosecl, sinecl; updateOrbElements(mjd); actTime = sgCalcActTime(mjd); @@ -66,10 +67,19 @@ void CelestialBody::updatePosition(double mjd, Star *ourSun) v = atan2(yv, xv); // the planet's true anomaly r = sqrt (xv*xv + yv*yv); // the planet's distance + // repetitive calculations, minimised for speed + cosN = cos(N); + sinN = sin(N); + cosvw = cos(v+w); + sinvw = sin(v+w); + sinvw_cosi = sinvw * cos(i); + cosecl = cos(ecl); + sinecl = sin(ecl); + // calculate the planet's position in 3D space - xh = r * (cos(N) * cos(v+w) - sin(N) * sin(v+w) * cos(i)); - yh = r * (sin(N) * cos(v+w) + cos(N) * sin(v+w) * cos(i)); - zh = r * (sin(v+w) * sin(i)); + xh = r * (cosN * cosvw - sinN * sinvw_cosi); + yh = r * (sinN * cosvw + cosN * sinvw_cosi); + zh = r * (sinvw * sin(i)); // calculate the ecliptic longitude and latitude xg = xh + ourSun->getxs(); @@ -80,8 +90,8 @@ void CelestialBody::updatePosition(double mjd, Star *ourSun) latEcl = atan2(zh, sqrt(xh*xh+yh*yh)); xe = xg; - ye = yg * cos(ecl) - zg * sin(ecl); - ze = yg * sin(ecl) + zg * cos(ecl); + ye = yg * cosecl - zg * sinecl; + ze = yg * sinecl + zg * cosecl; rightAscension = atan2(ye, xe); declination = atan2(ze, sqrt(xe*xe + ye*ye)); /* SG_LOG(SG_GENERAL, SG_INFO, "Planet found at : "