To find the position of the Sun, the first thing we need to do is set the epoch. The epoch is very important since a lot of the values used in astronomical calculations changes over time due to things like precession. The other reason for the epoch is to know the starting point of our calculation. mybloggingplanet.com

The trick to finding the Sun’s position is to know the position a a particular epoch, and then, based on the time interval between the epoch and now, determine how far the Sun has moved and thus determine where it is.

The things we need to be able to calculate the sun’s position are the Sun’s mean ecliptic longitude, perigee (point where the earth is closest to the Sun), eccentricity of the Earth’s orbit. All these values are adjusted to be to the values at the time of the epoch set.

We then calculate the eccentric anomaly which I covered in another tutorial, and then find the true anomaly from this, and at the end we get a value for the ecliptic longitude.

Since the Sun travels along the ecliptic, we can set the ecliptic latitude to 0.

The last step in the calculation is to convert the values to equatorial coordinates.

		public static void CalcSunPos(DateTime dDate, DateTime dEpoch, ref double fRA, ref double fDecl)
		{
			double fD;
			double fN;
			double fM;
			double fE;
			double fLambda;
			double fBeta;
			double fSolarMEL;
			double fSolarPL;
			double fSunEarthEcc;
			double fOblique;

            fBeta = 0;
			fD = UraniaTime.GetDaysBetween(dDate, dEpoch);
			fSolarMEL = GetSolarMEL(dEpoch, true);
			fSolarPL = GetSolarPerigeeLong(dEpoch, true);
			fSunEarthEcc = GetSunEarthEcc(dEpoch, true);
			fOblique = GetEarthObliquity(dEpoch, true);
			fN = (360.0 / 365.242191) * fD;
			fN = Trig.PutIn360Deg(fN);
			fM = fN + fSolarMEL - fSolarPL;
			fM = Trig.PutIn360Deg(fM);
			fM = Trig.DegToRad(fM);
   			fE = CalcEccentricAnomaly(fM, fM, fSunEarthEcc, 0.0000001);
			fTanV2 = Math.Sqrt((1.0 + fSunEarthEcc) / (1.0 - fSunEarthEcc)) * Math.Tan(fE / 2.0);
			fV = Math.Atan(fTanV2) * 2.0;

            fLambda = fN + fSolarMEL + fE;
			fLambda = Trig.PutIn360Deg(fLambda);

			UraniaCoord.ConvEclToEqu(fOblique, fLambda, fBeta, ref fRA, ref fDecl);
		}
Share