Earth's orbit

With my solar meridian data and a decent clock (my computer's clock), it was possible to tease out information about the Earth's orbit. The meridian altitudes gave the declination of the Sun. With the Sun's declination, the Sun's ecliptic longitude was given succinctly by the formula:

E_long (radians) = inverse sin ( sin(dec)/sin(EPS)), where EPS is the tilt of the Earth's axis to its orbit, as found earlier.

The solar ecliptic longitude starts at zero when the Sun's declination is zero and increasing (the spring equinox). If the Earth's orbit were circular, then the E_long would be simply: E_long = (360 deg/365.242 days) * (days since the spring equinox).

Here is a plot of my determination of the solar ecliptic longitude from October 2020 to October 2021:

Here is the linear fit to the ecliptic longitude. If the Earth's orbit were circular, then the residuals from the fit would be random and small:

It is clear that there are deviations from linearity. Near the center of the plot, the data values tend to be slightly above the fit, and at the ends of the plot the data values fall slightly below the fit. The slope is 0.98492 degrees/day, the average motion of the Sun along the ecliptic. Here are the residuals:

There are several things to note here. The Sun gets almost two degrees ahead or behind the longitude it should have if the Earth's orbit is circular. Following Kepler, we assume that the orbit is slightly elliptical, so let's try to find the eccentricity and the time of perihelion and aphelion. The time derivative of this curve tells us when the longitude is increasing or decreasing the fastest. This ecliptic velocity will be a maximum at perihelion and a minimum at aphelion. I fitted a curve to each side of this data plot and then took the derivative of that fit to find the velocities.

Here is the cubic fit to the ascending branch of the data plot. y = M0 + M1*days + M2*days^2 + M3*days^3 .

The excess velocity is the first derivative: y_dot = M1 + 2*M2*days + 3*M3*days^2 .

The excess velocity reaches a maximum when the second derivative equals zero: 0 = 2*M2 + 6*M3*days

Insert the values for M2 and M3 and solve for days: MJD day of perihelion = 0.0022207/(3*1.9928E-6) = 371.45

( JD= 2458849.5 + 371.45 = UT 11am, Jan 6, 2021. Modern computations give UT 13:51 on Jan 2, 2021)

Substitute this MJD day back into the equation for y_dot to get the maximum ecliptic excess velocity = 0.03777 , then add in the linear fit velocity 0.98492 to get = 1.0227 degrees/day for the perihelion angular velocity.

Now the same for the descending branch:

Here is the cubic fit to the descending branch of the data plot. y = M0 + M1*days + M2*days^2 + M3*days^3 .

The excess velocity is the first derivative: y_dot = M1 + 2*M2*days + 3*M3*days^2 .

The excess velocity reaches a minimum when the second derivative equals zero: 0 = 2*M2 + 6*M3*days

Insert the values for M2 and M3 and solve for days: MJD day of aphelion = 0.0019116/(3*1.1583E-6) = 550.01

(JD 2458849.5 + 550.01 = UT 0:15am, July 4, 2021. Modern computations give UT 22:27 July 5, 2021)

Substitute this MJD day back into the equation for y_dot to get the minimum ecliptic velocity = 0.98492 - 0.0301 = 0.9548 degrees/day


Now let's go back to the solar diameter page and get the max and min solar diameters of 0.539 and 0.525 degrees. Note that the difference between our two angular diameters is 0.539 - 0.525 = 0.014 degrees, with two significant figures at best. The angular diameter of an object that has linear size S and is at a linear distance of L will be S/L radians or 57.29*S/L degrees. If the angular diameter is 0.539 degrees, then L1 = (57.29/0.539)*S or 106.3 linear solar diameters. If the angular diameter is 0.525 degrees, then L2 = (57.29/0.525)*S or 109.12 linear solar diameters.

This means that at perihelion the Earth is 106.3 solar diameters from the Sun's center (about 92 million miles), and at aphelion it is 109.1 solar diameters (about 94 million miles). The simple formula for eccentricity e = (L2-L1)/(L2+L1) gives e = 0.013 which is about 20% smaller than the currently known value of 0.0167.


Let's see if we can use the higher accuracy numbers for the perihelion/aphelion angular velocities to get a better value for the eccentricity. I will use Va and Vp for the linear velocities of the Earth at aphelion and perihelion; Ra and Rp for the aphelion and perihelion distances from the Sun; and Wa and Wp as the angular velocities at aphelion and perihelion. I use W to substitute for the Greek lowercase omega that typically stands for angular velocity. From above, Wa = 0.9548 deg/day and Wp = 1.0227 deg/day.

Conservation of angular momentum gives us Kepler's law of equal areas in equal times.

This can be stated as Va*Ra = Vp*Rp.

Using the definition of angular velocity, Va = Wa*Ra and Vp=Wp*Rp.

Combining these two equations gives us Wa*Ra^2 = Wp*Rp^2.

Taking the square root of both sides of this equation gives sqrt(Wa)*Ra = sqrt(Wp)*Rp.

This equation can be re-expressed as Ra/Rp = sqrt(Wp/Wa).

The eccentricity of an orbit can be written as e = 1 - (2/((Ra/Rp)+1)).

Substituting in our expression for angular velocities, we get e = 1 - (2/(sqrt(Wp/Wa)+1)).

Further, substituting in my measured values I get e = 0.0172, within about 3% of the current value of 0.0167.


Tycho was smart enough to know that a good clock would ease his work immensely. Unfortunately for him, the clock technology in his day was so poor that he could not trust his clocks to be accurate. He kept several clocks of different types, only to find that they didn't agree with each other even over the course of a single day. See Rosa's presentation for details.

By finding only solar meridian altitudes, and timed with just my computer's clock, I was able to determine the Earth's perihelion and aphelion dates to within a few days, and calculated the eccentricity of the Earth's orbit to within a few percent.


Table of Contents Previous Page Next Page