home *** CD-ROM | disk | FTP | other *** search
/ Usenet 1994 January / usenetsourcesnewsgroupsinfomagicjanuary1994.iso / sources / x / volume19 / xephem / part19 / anomaly.c next >
Encoding:
C/C++ Source or Header  |  1993-05-15  |  905 b   |  45 lines
  1. #include <stdio.h>
  2. #include <math.h>
  3. #include "astro.h"
  4.  
  5. #define    TWOPI    (2*PI)
  6.  
  7. /* given the mean anomaly, ma, and the eccentricity, s, of elliptical motion,
  8.  * find the true anomaly, *nu, and the eccentric anomaly, *ea.
  9.  * all angles in radians.
  10.  */
  11. void
  12. anomaly (ma, s, nu, ea)
  13. double ma, s;
  14. double *nu, *ea;
  15. {
  16.     double m, fea;
  17.  
  18.     m = ma-TWOPI*(long)(ma/TWOPI);
  19.     if (m > PI) m -= TWOPI;
  20.     if (m < -PI) m += TWOPI;
  21.     fea = m;
  22.  
  23.     if (s < 1.0) {
  24.         /* elliptical */
  25.         double dla;
  26.         for (;;) {
  27.         dla = fea-(s*sin(fea))-m;
  28.         if (fabs(dla)<1e-6)
  29.             break;
  30.         dla /= 1-(s*cos(fea));
  31.         fea -= dla;
  32.         }
  33.         *nu = 2*atan(sqrt((1+s)/(1-s))*tan(fea/2));
  34.         } else {
  35.         /* hyperbolic */
  36.         double corr = 1;
  37.         while (fabs(corr) > 0.000001) {
  38.           corr = (m - s * sinh(fea) + fea) / (s*cosh(fea) - 1);
  39.           fea += corr;
  40.         }
  41.         *nu = 2*atan(sqrt((s+1)/(s-1))*tanh(fea/2));
  42.     }
  43.     *ea = fea;
  44. }
  45.