xref: /plan9/sys/src/cmd/astro/venus.c (revision 7dd7cddf99dd7472612f1413b4da293630e6b1bc)
13e12c5d1SDavid du Colombier #include "astro.h"
23e12c5d1SDavid du Colombier 
33e12c5d1SDavid du Colombier 
43e12c5d1SDavid du Colombier void
venus(void)53e12c5d1SDavid du Colombier venus(void)
63e12c5d1SDavid du Colombier {
73e12c5d1SDavid du Colombier 	double pturbl, pturbb, pturbr;
83e12c5d1SDavid du Colombier 	double lograd;
93e12c5d1SDavid du Colombier 	double dele, enom, vnom, nd, sl;
103e12c5d1SDavid du Colombier 	double v0, t0, m0, j0, s0;
113e12c5d1SDavid du Colombier 	double lsun, elong, ci, dlong;
123e12c5d1SDavid du Colombier 
133e12c5d1SDavid du Colombier /*
143e12c5d1SDavid du Colombier  *	here are the mean orbital elements
153e12c5d1SDavid du Colombier  */
163e12c5d1SDavid du Colombier 
173e12c5d1SDavid du Colombier 	ecc = .00682069 - .00004774*capt + 0.091e-6*capt2;
183e12c5d1SDavid du Colombier 	incl = 3.393631 + .0010058*capt - 0.97e-6*capt2;
193e12c5d1SDavid du Colombier 	node = 75.779647 + .89985*capt + .00041*capt2;
203e12c5d1SDavid du Colombier 	argp = 130.163833 + 1.408036*capt - .0009763*capt2;
213e12c5d1SDavid du Colombier 	mrad = .7233316;
223e12c5d1SDavid du Colombier 	anom = 212.603219 + 1.6021301540*eday + .00128605*capt2;
233e12c5d1SDavid du Colombier 	motion = 1.6021687039;
243e12c5d1SDavid du Colombier 
253e12c5d1SDavid du Colombier /*
263e12c5d1SDavid du Colombier  *	mean anomalies of perturbing planets
273e12c5d1SDavid du Colombier  */
283e12c5d1SDavid du Colombier 
293e12c5d1SDavid du Colombier 	v0 = 212.60 + 1.602130154*eday;
303e12c5d1SDavid du Colombier 	t0 = 358.63  + .985608747*eday;
313e12c5d1SDavid du Colombier 	m0 = 319.74 + 0.524032490*eday;
323e12c5d1SDavid du Colombier 	j0 = 225.43 + .083090842*eday;
333e12c5d1SDavid du Colombier 	s0 = 175.8  + .033459258*eday;
343e12c5d1SDavid du Colombier 
353e12c5d1SDavid du Colombier 	v0 *= radian;
363e12c5d1SDavid du Colombier 	t0 *= radian;
373e12c5d1SDavid du Colombier 	m0 *= radian;
383e12c5d1SDavid du Colombier 	j0 *= radian;
393e12c5d1SDavid du Colombier 	s0 *= radian;
403e12c5d1SDavid du Colombier 
413e12c5d1SDavid du Colombier 	incl *= radian;
423e12c5d1SDavid du Colombier 	node *= radian;
433e12c5d1SDavid du Colombier 	argp *= radian;
443e12c5d1SDavid du Colombier 	anom = fmod(anom, 360.)*radian;
453e12c5d1SDavid du Colombier 
463e12c5d1SDavid du Colombier /*
473e12c5d1SDavid du Colombier  *	computation of long period terms affecting the mean anomaly
483e12c5d1SDavid du Colombier  */
493e12c5d1SDavid du Colombier 
503e12c5d1SDavid du Colombier 	anom +=
513e12c5d1SDavid du Colombier 		   (2.761-0.022*capt)*radsec*sin(
523e12c5d1SDavid du Colombier 		  13.*t0 - 8.*v0 + 43.83*radian + 4.52*radian*capt)
533e12c5d1SDavid du Colombier 		 + 0.268*radsec*cos(4.*m0 - 7.*t0 + 3.*v0)
543e12c5d1SDavid du Colombier 		 + 0.019*radsec*sin(4.*m0 - 7.*t0 + 3.*v0)
553e12c5d1SDavid du Colombier 		 - 0.208*radsec*sin(s0 + 1.4*radian*capt);
563e12c5d1SDavid du Colombier 
573e12c5d1SDavid du Colombier /*
583e12c5d1SDavid du Colombier  *	computation of elliptic orbit
593e12c5d1SDavid du Colombier  */
603e12c5d1SDavid du Colombier 
613e12c5d1SDavid du Colombier 	enom = anom + ecc*sin(anom);
623e12c5d1SDavid du Colombier 	do {
633e12c5d1SDavid du Colombier 		dele = (anom - enom + ecc * sin(enom)) /
643e12c5d1SDavid du Colombier 			(1 - ecc*cos(enom));
653e12c5d1SDavid du Colombier 		enom += dele;
66*7dd7cddfSDavid du Colombier 	} while(fabs(dele) > converge);
673e12c5d1SDavid du Colombier 	vnom = 2*atan2(sqrt((1+ecc)/(1-ecc))*sin(enom/2),
683e12c5d1SDavid du Colombier 		cos(enom/2));
693e12c5d1SDavid du Colombier 	rad = mrad*(1 - ecc*cos(enom));
703e12c5d1SDavid du Colombier 
713e12c5d1SDavid du Colombier 	lambda = vnom + argp;
723e12c5d1SDavid du Colombier 
733e12c5d1SDavid du Colombier /*
743e12c5d1SDavid du Colombier  *	perturbations in longitude
753e12c5d1SDavid du Colombier  */
763e12c5d1SDavid du Colombier 
773e12c5d1SDavid du Colombier 	icosadd(venfp, vencp);
783e12c5d1SDavid du Colombier 	pturbl = cosadd(4, v0, t0, m0, j0);
793e12c5d1SDavid du Colombier 	pturbl *= radsec;
803e12c5d1SDavid du Colombier 
813e12c5d1SDavid du Colombier /*
823e12c5d1SDavid du Colombier  *	perturbations in latidude
833e12c5d1SDavid du Colombier  */
843e12c5d1SDavid du Colombier 
853e12c5d1SDavid du Colombier 	pturbb = cosadd(3, v0, t0, j0);
863e12c5d1SDavid du Colombier 	pturbb *= radsec;
873e12c5d1SDavid du Colombier 
883e12c5d1SDavid du Colombier /*
893e12c5d1SDavid du Colombier  *	perturbations in log radius vector
903e12c5d1SDavid du Colombier  */
913e12c5d1SDavid du Colombier 
923e12c5d1SDavid du Colombier 	pturbr = cosadd(4, v0, t0, m0, j0);
933e12c5d1SDavid du Colombier 
943e12c5d1SDavid du Colombier /*
953e12c5d1SDavid du Colombier  *	reduction to the ecliptic
963e12c5d1SDavid du Colombier  */
973e12c5d1SDavid du Colombier 
983e12c5d1SDavid du Colombier 	lambda += pturbl;
993e12c5d1SDavid du Colombier 	nd = lambda - node;
1003e12c5d1SDavid du Colombier 	lambda = node + atan2(sin(nd)*cos(incl),cos(nd));
1013e12c5d1SDavid du Colombier 
1023e12c5d1SDavid du Colombier 	sl = sin(incl)*sin(nd);
1033e12c5d1SDavid du Colombier 	beta = atan2(sl, pyth(sl)) + pturbb;
1043e12c5d1SDavid du Colombier 
1053e12c5d1SDavid du Colombier 	lograd = pturbr*2.30258509;
1063e12c5d1SDavid du Colombier 	rad *= 1 + lograd;
1073e12c5d1SDavid du Colombier 
1083e12c5d1SDavid du Colombier 
1093e12c5d1SDavid du Colombier 	motion *= radian*mrad*mrad/(rad*rad);
1103e12c5d1SDavid du Colombier 
1113e12c5d1SDavid du Colombier /*
1123e12c5d1SDavid du Colombier  *	computation of magnitude
1133e12c5d1SDavid du Colombier  */
1143e12c5d1SDavid du Colombier 
1153e12c5d1SDavid du Colombier 	lsun = 99.696678 + 0.9856473354*eday;
1163e12c5d1SDavid du Colombier 	lsun *= radian;
1173e12c5d1SDavid du Colombier 	elong = lambda - lsun;
1183e12c5d1SDavid du Colombier 	ci = (rad - cos(elong))/sqrt(1 + rad*rad - 2*rad*cos(elong));
1193e12c5d1SDavid du Colombier 	dlong = atan2(pyth(ci), ci)/radian;
1203e12c5d1SDavid du Colombier 	mag = -4 + .01322*dlong + .0000004247*dlong*dlong*dlong;
1213e12c5d1SDavid du Colombier 
1223e12c5d1SDavid du Colombier 	semi = 8.41;
1233e12c5d1SDavid du Colombier 
1243e12c5d1SDavid du Colombier 	helio();
1253e12c5d1SDavid du Colombier 	geo();
1263e12c5d1SDavid du Colombier }
127