Finding the position of the Moon is in many ways similar to finding the position of the Sun.

The Moon moves along a path that is near the ecliptic but is slightly offset – by 5.14 degrees. If this were not so, then we would get eclipses every new and full moon.

The other information we need to calculate the Moon’s position, is the ascending node, which is the point where the orbit of the moon crosses the ecliptic, which is 318.351648 degrees, and the eccentricity of the orbit, which is 0.054900. The perigee of the Moon is also needed which is 36.340410 degrees.

We then need to know the position of the Sun too, so we also need the solar ecliptic longitude, perigee, and eccentricity.

The last thing we need now is a starting point for our calculation, so for the epoch given, we need to supply the ecliptic longitude for that time, from which the current position can be calculated.

		public static void CalcMoonPos(DateTime dDate, DateTime dEpoch, double fMEpochLong, double fMPeriLong, double fMAscNode, double fMIncl, double fMEcc, double fSEpochEclLong, double fSPeriEclLong, double fSEcc, ref double fMRA, ref double fMDecl)
		{
			double fN, fSM, fSE, fSLambda;
			double fL, fMM, fMN, fME, fAE, fMEC, fA3, fA4, fMV, fMM1, fL1, fL2, fN1, fX, fY;
			double fT, fMLambda, fMBeta, fJD1, fJD2, fDays;

			fJD1 = UraniaTime.GetJulianDay(dDate, 0);
			fJD2 = UraniaTime.GetJulianDay(dEpoch, 0);
			fDays = (fJD1 - fJD2);
			fDays += 1;

			fN = (360.0/365.242191) * fDays;
			fN = Trig.PutIn360Deg(fN);
			fSM = fN + fSEpochEclLong - fSPeriEclLong;
			fSM = Trig.PutIn360Deg(fSM);

			fSE = (360.0 / Math.PI) * fSEcc * Math.Sin(Trig.DegToRad(fSM));
			fSLambda = fN + fSE + fSEpochEclLong;

			fL = (13.176396 * fDays) + fMEpochLong;
			fL = Trig.PutIn360Deg(fL);
			
			fMM = fL - (0.111404 * fDays) - fMPeriLong;
			fMM = Trig.PutIn360Deg(fMM);

			fMN = fMAscNode - (0.0529539 * fDays);
			fMN = Trig.PutIn360Deg(fMN);

			fME = 1.2739 * Trig.Sin((2.0 * (fL - fSLambda)) - fMM);
			fAE = 0.1858 * Trig.Sin(fSM);
			fA3 = 0.37 * Trig.Sin(fSM);

			fMM1 = fMM + fME - fAE + fA3;
			
			fMEC = 6.2886 * Trig.Sin(fMM1);
			fA4 = 0.214 * Trig.Sin(2.0 * fMM1);
			fL1 = fL + fME + fMEC - fAE + fA4;

			fMV = 0.6583 * Trig.Sin(2.0 * (fL1 - fSLambda));
			fL2 = fL1 + fMV;

			fN1 = fMN - (0.16 * Trig.Sin(fSM));
			fY = Trig.Sin(fL2 - fN1) * Trig.Cos(fMIncl);
			fX = Trig.Cos(fL2 - fN1);

			fT = Trig.Atan(fY / fX);

			fT = Trig.TanQuadrant(fX, fY, fT);
			
			fMLambda = fT + fN1;
			fMBeta = Trig.Asin(Trig.Sin(fL2 - fN1) * Trig.Sin(fMIncl));
			UraniaCoord.ConvEclToEqu(23.441884,fMLambda, fMBeta, ref fMRA, ref fMDecl);
		}

Share