| [2551] | 1 | #include <stdio.h> | 
|---|
|  | 2 | #include <math.h> | 
|---|
|  | 3 |  | 
|---|
|  | 4 | #include "astro.h" | 
|---|
|  | 5 |  | 
|---|
|  | 6 | static void m (double t, double k, double *mj); | 
|---|
|  | 7 |  | 
|---|
|  | 8 | #define unw(w,z)        ((w)-floor((w)/(z))*(z)) | 
|---|
|  | 9 |  | 
|---|
|  | 10 | /* given a modified Julian date, mj, return the mjd of the new | 
|---|
|  | 11 | * and full moons about then, mjn and mjf. | 
|---|
|  | 12 | * TODO: exactly which ones does it find? eg: | 
|---|
|  | 13 | *   5/28/1988 yields 5/15 and 5/31 | 
|---|
|  | 14 | *   5/29             6/14     6/29 | 
|---|
|  | 15 | */ | 
|---|
|  | 16 | void | 
|---|
|  | 17 | moonnf (double mj, double *mjn, double *mjf) | 
|---|
|  | 18 | { | 
|---|
|  | 19 | int mo, yr; | 
|---|
|  | 20 | double dy; | 
|---|
|  | 21 | double mj0; | 
|---|
|  | 22 | double k, tn, tf, t; | 
|---|
|  | 23 |  | 
|---|
|  | 24 | mjd_cal (mj, &mo, &dy, &yr); | 
|---|
|  | 25 | cal_mjd (1, 0., yr, &mj0); | 
|---|
|  | 26 | k = (yr-1900+((mj-mj0)/365))*12.3685; | 
|---|
|  | 27 | k = floor(k+0.5); | 
|---|
|  | 28 | tn = k/1236.85; | 
|---|
|  | 29 | tf = (k+0.5)/1236.85; | 
|---|
|  | 30 | t = tn; | 
|---|
|  | 31 | m (t, k, mjn); | 
|---|
|  | 32 | t = tf; | 
|---|
|  | 33 | k += 0.5; | 
|---|
|  | 34 | m (t, k, mjf); | 
|---|
|  | 35 | } | 
|---|
|  | 36 |  | 
|---|
|  | 37 | static void | 
|---|
|  | 38 | m (double t, double k, double *mj) | 
|---|
|  | 39 | { | 
|---|
|  | 40 | double t2, a, a1, b, b1, c, ms, mm, f, ddjd; | 
|---|
|  | 41 |  | 
|---|
|  | 42 | t2 = t*t; | 
|---|
|  | 43 | a = 29.53*k; | 
|---|
|  | 44 | c = degrad(166.56+(132.87-9.173e-3*t)*t); | 
|---|
|  | 45 | b = 5.8868e-4*k+(1.178e-4-1.55e-7*t)*t2+3.3e-4*sin(c)+7.5933E-1; | 
|---|
|  | 46 | ms = 359.2242+360*unw(k/1.236886e1,1)-(3.33e-5+3.47e-6*t)*t2; | 
|---|
|  | 47 | mm = 306.0253+360*unw(k/9.330851e-1,1)+(1.07306e-2+1.236e-5*t)*t2; | 
|---|
|  | 48 | f = 21.2964+360*unw(k/9.214926e-1,1)-(1.6528e-3+2.39e-6*t)*t2; | 
|---|
|  | 49 | ms = unw(ms,360); | 
|---|
|  | 50 | mm = unw(mm,360); | 
|---|
|  | 51 | f = unw(f,360); | 
|---|
|  | 52 | ms = degrad(ms); | 
|---|
|  | 53 | mm = degrad(mm); | 
|---|
|  | 54 | f = degrad(f); | 
|---|
|  | 55 | ddjd = (1.734e-1-3.93e-4*t)*sin(ms)+2.1e-3*sin(2*ms) | 
|---|
|  | 56 | -4.068e-1*sin(mm)+1.61e-2*sin(2*mm)-4e-4*sin(3*mm) | 
|---|
|  | 57 | +1.04e-2*sin(2*f)-5.1e-3*sin(ms+mm)-7.4e-3*sin(ms-mm) | 
|---|
|  | 58 | +4e-4*sin(2*f+ms)-4e-4*sin(2*f-ms)-6e-4*sin(2*f+mm) | 
|---|
|  | 59 | +1e-3*sin(2*f-mm)+5e-4*sin(ms+2*mm); | 
|---|
|  | 60 | a1 = (long)a; | 
|---|
|  | 61 | b = b+ddjd+(a-a1); | 
|---|
|  | 62 | b1 = (long)b; | 
|---|
|  | 63 | a = a1+b1; | 
|---|
|  | 64 | b = b-b1; | 
|---|
|  | 65 | *mj = a + b; | 
|---|
|  | 66 | } | 
|---|
|  | 67 |  | 
|---|
|  | 68 | /* For RCS Only -- Do Not Edit */ | 
|---|
| [2818] | 69 | static char *rcsid[2] = {(char *)rcsid, "@(#) $RCSfile: moonnf.c,v $ $Date: 2005-08-21 10:02:38 $ $Revision: 1.3 $ $Name: not supported by cvs2svn $"}; | 
|---|