I enjoy astronomy. In this post I’m going to show off a cool file from the File Exchange and demonstrate MATLAB’s new-ish datetime functionality (including dateshift and duration). But it’s also a fun opportunity to tell a little astronomical story. A story that starts with a simple question: Who gets to decide when it’s noon? The sun, or your wristwatch?

Noon is defined to be the moment that divides morning from afternoon. Astronomically speaking, this is when the sun stops getting higher and starts getting lower. But practically speaking, it’s when Mickey’s hands point to twelve.

Who is right?

Mickey is right. How can that be? Why should the watch on your wrist have more authority than the star that we orbit?

Consider that if your timepiece is a sundial, there’s no conflict. By definition, it’s noon when the sun reaches its highest point in the sky, or when the gnomon’s shadow is the shortest (the gnomon is the pointy thing that casts the shadow). Sundials are simple and sturdy, but they have some issues. They don’t work very well at night, but more to the point, they imply an infinite number of time zones, one for every possible line of longitude.

Even this inconvenience didn’t really matter until the arrival of the railroad and the telegraph. Only then, and for the first time in history, did people travel so quickly that anyone actually cared that noon in New York wasn’t quite the same as noon in Washington DC. In particular, constructing railroad timetables was a nightmare.

To avoid the infinite time zones problem, we bin meridians into 24 groups, one for each hour of the day. One special meridian in each time zone uses the sun to determine noon, and the rest of the zone falls into line. Problem solved, right? Not quite.

Here’s where we need to use some MATLAB. I’ll be using Darin Koblick’s excellent Vectorized Solar Azimuth and Elevation Estimation. I added the file using the Add-On Explorer. Here’s what it looks like in the Add-On Manager.

This file can tell me where the sun is right now. It returns the answer as azimuth (north = 0 deg, east = 90 deg) and elevation (horizon = 0 deg, straight up = 90 deg).

Here are the coordinates where I sit.

lat =  42.299899; % deg
lon = -71.350729; % deg
alt =   0.305;    % km
tz  =  hours(timezone(lon));  % time zone adjustment from universal time

At this instant, where is the sun?

dtNow = datetime('now');
dtNow.Format = 'yyyy/MM/dd HH:mm:ss';
[az,el] = SolarAzEl(dtNow+tz,lat,lon,alt);
fprintf('Solar Azimuth:   %5.2f deg\nSolar Elevation:  %5.2f deg\n',az,el)

Solar Azimuth:   221.82 deg
Solar Elevation:  17.26 deg

It’s a little west of south (mid-afternoon) and somewhat low in the sky, as befits late January.

Here is the elevation of the sun in the sky for every minute between 11 AM and 1 PM local time. By taking the maximum elevation and finding the corresponding time, we’ve determined local noon, which is to say, the astronomical moment at which the sun is highest.

% Make a list of all the minutes from 11 AM to 1 PM
timeList = datetime(year(dtNow),month(dtNow),day(dtNow),11,0:120,0)';
% Find the corresponding solar elevations
[az,el] = SolarAzEl(timeList+tz,lat,lon,alt);
sun = [0.96 0.74 0.24];
sky = [0.85 1 1];
plot(timeList,el, ...
    'LineWidth',1,'Color',sun)
set(gca,'Color',sky)
[~,ix] = max(el);
noon = timeList(ix);
noonEl = el(ix);
line(noon,noonEl, ...
    'Color',sun, ...
    'Marker','.','MarkerSize',40)
grid on
xlabel('Time')
ylabel('Elevation Above Horizon (deg)')
noon.Format = 'HH:mm';
title(sprintf('Local Noon Occurs at %s',noon))

It’s close, but “noon” isn’t quite noon. Part of that effect is because of where I sit in the time zone. But that’s not all there is to it.

Let’s look at what happens over a period of few days. These are the solar elevations for the first two weeks of the year. You can see noon is shifting slightly later every day, and the sun is climbing slightly higher every day.

noons = datetime(year(dtNow),1,2:14);
for n = 1:length(noons)
    timeList = datetime(year(noons(n)),month(noons(n)),day(noons(n)),11,30:90,0)';
    [~,el] = SolarAzEl(timeList+tz,lat,lon,alt);
    [~,ix] = max(el);
    plot(timeList-dateshift(timeList,'start','day'),el,'Color',sun)
    hold all
    plot(timeList(ix) - dateshift(timeList(ix),'start','day'),el(ix), ...
    'Marker','.','MarkerSize',32, ...
    'Color',sun)
end
hold off
grid on
xlabel('Time')
ylabel('Elevation Above Horizon (deg)')
set(gca,'Color',sky)

It looks a little bumpy because the times are discretized to the nearest minute.

Now let’s plot where noon occurs for the entire year.

noons = datetime(year(dtNow),1,1:7:365);
els = zeros(size(noons));
for n = 1:length(noons)
    timeList = datetime(year(noons(n)),month(noons(n)),day(noons(n)),11,0:120,0)';
    [~,el] = SolarAzEl(timeList+tz,lat,lon,alt);
    [~,ix] = max(el);
    noons(n) = timeList(ix);
    els(n) = el(ix);
end
plot(noons - dateshift(noons,'start','day'),els, ...
    'Marker','.','MarkerSize',32, ...
    'Color',sun)
set(gca,'Color',sky)
grid on
xlim(duration(hours([11 12.5])))
xlabel('Time')
ylabel('Elevation (degrees)')
title('Elevation of the Sun at Local Noon')

Yow! Where did that figure 8 shape come from?

mean(noons - dateshift(noons,'start','day'))

ans = duration
   11:45:29

On average, local noon occurred at 11:45. But sometimes it was earlier and sometimes it was later by as much as 15 minutes. There’s clearly more going on here than just time zone shifts.

There are two subtle effects here. One is due to the fact that the earth’s orbit isn’t quite circular. It’s slightly elliptical, so the earth actually moves faster around the sun at certain times of the year. The other effect is due to the fact that the earth is tilted on its axis by 23.4 deg. These two factors make a the length of a day, as measured from local noon to local noon, vary by as much as a half hour. Crazy, right?

At some point hours and minutes were redefined as fractions of a single standard “mean solar day” rather than a locally measured (but inconsistent) “sundial day”. We’ve ironed out all the wrinkles in time. Put another way, first the sun defined minutes, then minutes defined the sun. The timekeeping tail wags the dog.

This is pretty basic stuff, as far as astronomy goes, but it was so much fun to find a nice file on the File Exchange and then be able to bang out this code so quickly. It’s amazing how much good stuff is waiting for you there. Go look!