1 cd5bae78 2004-04-21 devnull #include "astro.h"
 2 cd5bae78 2004-04-21 devnull 
 3 cd5bae78 2004-04-21 devnull void
 4 cd5bae78 2004-04-21 devnull mars(void)
 5 cd5bae78 2004-04-21 devnull {
 6 cd5bae78 2004-04-21 devnull 	double pturbl, pturbb, pturbr;
 7 cd5bae78 2004-04-21 devnull 	double lograd;
 8 cd5bae78 2004-04-21 devnull 	double dele, enom, vnom, nd, sl;
 9 cd5bae78 2004-04-21 devnull 	double lsun, elong, ci, dlong;
10 cd5bae78 2004-04-21 devnull 
11 cd5bae78 2004-04-21 devnull 
12 cd5bae78 2004-04-21 devnull 	ecc = .09331290 + .000092064*capt;
13 cd5bae78 2004-04-21 devnull 	incl = 1.850333 - 6.75e-4*capt;
14 cd5bae78 2004-04-21 devnull 	node = 48.786442 + .770992*capt;
15 cd5bae78 2004-04-21 devnull 	argp = 334.218203 + 1.840758*capt + 1.30e-4*capt2;
16 cd5bae78 2004-04-21 devnull 	mrad = 1.5236915;
17 cd5bae78 2004-04-21 devnull 	anom = 319.529425 + .5240207666*eday + 1.808e-4*capt2;
18 cd5bae78 2004-04-21 devnull 	motion = 0.5240711638;
19 cd5bae78 2004-04-21 devnull 
20 cd5bae78 2004-04-21 devnull 
21 cd5bae78 2004-04-21 devnull 	incl = incl*radian;
22 cd5bae78 2004-04-21 devnull 	node = node*radian;
23 cd5bae78 2004-04-21 devnull 	argp = argp*radian;
24 cd5bae78 2004-04-21 devnull 	anom = fmod(anom,360.)*radian;
25 cd5bae78 2004-04-21 devnull 
26 cd5bae78 2004-04-21 devnull 	enom = anom + ecc*sin(anom);
27 cd5bae78 2004-04-21 devnull 	do {
28 cd5bae78 2004-04-21 devnull 		dele = (anom - enom + ecc * sin(enom)) /
29 cd5bae78 2004-04-21 devnull 			(1. - ecc*cos(enom));
30 cd5bae78 2004-04-21 devnull 		enom += dele;
31 cd5bae78 2004-04-21 devnull 	} while(fabs(dele) > converge);
32 cd5bae78 2004-04-21 devnull 	vnom = 2.*atan2(sqrt((1.+ecc)/(1.-ecc))*sin(enom/2.),cos(enom/2.));
33 cd5bae78 2004-04-21 devnull 	rad = mrad*(1. - ecc*cos(enom));
34 cd5bae78 2004-04-21 devnull 
35 cd5bae78 2004-04-21 devnull 	lambda = vnom + argp;
36 cd5bae78 2004-04-21 devnull 	pturbl = 0.;
37 cd5bae78 2004-04-21 devnull 	lambda = lambda + pturbl*radsec;
38 cd5bae78 2004-04-21 devnull 	pturbb = 0.;
39 cd5bae78 2004-04-21 devnull 	pturbr = 0.;
40 cd5bae78 2004-04-21 devnull 
41 cd5bae78 2004-04-21 devnull /*
42 cd5bae78 2004-04-21 devnull  *	reduce to the ecliptic
43 cd5bae78 2004-04-21 devnull  */
44 cd5bae78 2004-04-21 devnull 
45 cd5bae78 2004-04-21 devnull 	nd = lambda - node;
46 cd5bae78 2004-04-21 devnull 	lambda = node + atan2(sin(nd)*cos(incl),cos(nd));
47 cd5bae78 2004-04-21 devnull 
48 cd5bae78 2004-04-21 devnull 	sl = sin(incl)*sin(nd) + pturbb*radsec;
49 cd5bae78 2004-04-21 devnull 	beta = atan2(sl, pyth(sl));
50 cd5bae78 2004-04-21 devnull 
51 cd5bae78 2004-04-21 devnull 	lograd = pturbr*2.30258509;
52 cd5bae78 2004-04-21 devnull 	rad *= 1. + lograd;
53 cd5bae78 2004-04-21 devnull 
54 cd5bae78 2004-04-21 devnull 
55 cd5bae78 2004-04-21 devnull 	motion *= radian*mrad*mrad/(rad*rad);
56 cd5bae78 2004-04-21 devnull 	semi = 4.68;
57 cd5bae78 2004-04-21 devnull 
58 cd5bae78 2004-04-21 devnull 	lsun = 99.696678 + 0.9856473354*eday;
59 cd5bae78 2004-04-21 devnull 	lsun *= radian;
60 cd5bae78 2004-04-21 devnull 	elong = lambda - lsun;
61 cd5bae78 2004-04-21 devnull 	ci = (rad - cos(elong))/sqrt(1. + rad*rad - 2.*rad*cos(elong));
62 cd5bae78 2004-04-21 devnull 	dlong = atan2(pyth(ci), ci)/radian;
63 cd5bae78 2004-04-21 devnull 	mag = -1.30 + .01486*dlong;
64 cd5bae78 2004-04-21 devnull 
65 cd5bae78 2004-04-21 devnull 	helio();
66 cd5bae78 2004-04-21 devnull 	geo();
67 cd5bae78 2004-04-21 devnull }