{"id":8735,"date":"2022-12-02T17:33:21","date_gmt":"2022-12-02T22:33:21","guid":{"rendered":"https:\/\/blogs.mathworks.com\/community\/?p=8735"},"modified":"2022-12-08T12:35:39","modified_gmt":"2022-12-08T17:35:39","slug":"crepuscular_tidings","status":"publish","type":"post","link":"https:\/\/blogs.mathworks.com\/community\/2022\/12\/02\/crepuscular_tidings\/","title":{"rendered":"Crepuscular Tidings"},"content":{"rendered":"
\n
As the days get shorter and drearier here in New England, I get hungry for sunlight. This makes me somewhat obsessed with when the sun sets. But unlike our ancient ancestors, I don't need to worry that the sun will never return. I have something they didn't have: MATLAB. To cope with the shortening days, I can calculate exactly when the days will start getting longer. And for me, the days will start getting longer on December 10th.<\/div>\n
You're probably aware that the shortest day is December 21st. But the earliest sunset is on December 10th (at my latitude). And since I never wake up before dawn, then the day as I experience it is already getting longer by December 11th. It's a subtle point, but it helps carry me through the twilight of the waning year (I've written about this phenomenon before<\/a>).<\/div>\n
Using Francois Beauducel's \"Sunrise\" contribution on the File Exchange<\/a>, I can easily make a plot of the sunset time where I live (north latitude 42 degrees).<\/div>\n
\n
\n
% Where: Natick, Massachusetts.<\/span><\/span><\/div>\n<\/div>\n
\n
lat = 42.3;<\/span><\/div>\n<\/div>\n
\n
lon = -71.4;<\/span><\/div>\n<\/div>\n
\n
<\/div>\n<\/div>\n
\n
% When: this year.<\/span><\/span><\/div>\n<\/div>\n
\n
% List all the days of the year<\/span><\/span><\/div>\n<\/div>\n
\n
d = datetime(2022,1,1:365);<\/span><\/div>\n<\/div>\n
\n
<\/div>\n<\/div>\n
\n
% SUNRISE uses the DATENUM format<\/span><\/span><\/div>\n<\/div>\n
\n
[~,sset_datenum] = sunrise(lat,lon,0,[],datenum(d));<\/span><\/div>\n<\/div>\n
\n
sset = datetime(sset_datenum,ConvertFrom=\"datenum\"<\/span>);<\/span><\/div>\n<\/div>\n
\n
<\/div>\n<\/div>\n
\n
t = hour(sset) + minute(sset)\/60 + second(sset)\/3600;<\/span><\/div>\n<\/div>\n
\n
<\/div>\n<\/div>\n
\n
plot(d,t,LineWidth=3)<\/span><\/div>\n<\/div>\n
\n
grid on<\/span><\/span><\/div>\n<\/div>\n
\n
box on<\/span><\/span><\/div>\n<\/div>\n
\n
ylim([15 20])<\/span><\/div>\n<\/div>\n
\n
ylabel(\"Time of Sunset\"<\/span>)<\/span><\/div>\n<\/div>\n
\n
[~,ix] = min(t);<\/span><\/div>\n<\/div>\n
\n
line(sset(ix),t(ix), ...<\/span><\/span><\/div>\n<\/div>\n
\n
Marker=\"o\"<\/span>, ...<\/span><\/span><\/div>\n<\/div>\n
\n
MarkerSize=15, ...<\/span><\/span><\/div>\n<\/div>\n
\n
LineWidth=2, ...<\/span><\/span><\/div>\n<\/div>\n
\n
Color=\"red\"<\/span>)<\/span><\/div>\n<\/div>\n
\n
text(sset(ix),t(ix), ...<\/span><\/span><\/div>\n<\/div>\n
\n
sprintf(\"Earliest Sunset \\n%s \"<\/span>,sset(ix)), ...<\/span><\/span><\/div>\n<\/div>\n
\n
Color=\"red\"<\/span>, ...<\/span><\/span><\/div>\n<\/div>\n
\n
VerticalAlignment=\"bottom\"<\/span>, ...<\/span><\/span><\/div>\n<\/div>\n
\n
HorizontalAlignment=\"right\"<\/span>)<\/span><\/div>\n<\/div>\n
\n
setTimeYTick(gca)<\/span><\/div>\n
\n
\n
<\/div>\n
<\/div>\n
<\/div>\n
<\/div>\n<\/div>\n<\/div>\n<\/div>\n<\/div>\n
Note that all times are Standard Time. We're not worrying about any Daylight Saving Time nonsense. I like to call the day with the earliest sunset Crepusculus (after crepusculum, the Latin word for dusk). It's an important milestone, so it deserves a name!<\/div>\n
Let's take that same data and use it to animate the sunset time using the face of a clock. Here we're visualizing the hour hand at the time of sunset. You can see a little asymmetry here. The days get shorter more quickly than they get longer.<\/div>\n
\"animate_clock_solstice.gif\"<\/div>\n
For this animation, I used a File Exchange contribution of my own called Animator<\/a>.<\/div>\n
Now let's plot the span from latest to earliest sunset for several north latitudes.<\/div>\n
\n
\n
clf<\/span><\/div>\n<\/div>\n
\n
hold on<\/span><\/span><\/div>\n<\/div>\n
\n
latList = 0:10:50;<\/span><\/div>\n<\/div>\n
\n
crepusculusTimes = zeros(size(latList));<\/span><\/div>\n<\/div>\n
\n
crepusculusDates = NaT(size(latList));<\/span><\/div>\n<\/div>\n
\n
colorMap = flipud(parula(10));<\/span><\/div>\n<\/div>\n
\n
for <\/span>n = 1:length(latList)<\/span><\/div>\n<\/div>\n
\n
d = datetime(2022,1,1:365);<\/span><\/div>\n<\/div>\n
\n
lat = latList(n);<\/span><\/div>\n<\/div>\n
\n
lon = -71;<\/span><\/div>\n<\/div>\n
\n
[~,sset] = sunrise_dt(lat,lon,d);<\/span><\/div>\n<\/div>\n
\n
t = hour(sset)+minute(sset)\/60+second(sset)\/3600;<\/span><\/div>\n<\/div>\n
\n
[~,ixEarliest] = min(t);<\/span><\/div>\n<\/div>\n
\n
[~,ixLatest] = max(t);<\/span><\/div>\n<\/div>\n
\n
sset = sset(ixLatest:ixEarliest);<\/span><\/div>\n<\/div>\n
\n
crepusculusTimes(n) = t(ixEarliest);<\/span><\/div>\n<\/div>\n
\n
crepusculusDates(n) = d(ixEarliest);<\/span><\/div>\n<\/div>\n
\n
plot(d,t,LineWidth=3,Color=colorMap(n,:))<\/span><\/div>\n<\/div>\n
\n
grid on<\/span><\/span><\/div>\n<\/div>\n
\n
hold on<\/span><\/span><\/div>\n<\/div>\n
\n
box on<\/span><\/span><\/div>\n<\/div>\n
\n
end<\/span><\/span><\/div>\n<\/div>\n
\n
legend(string(latList + \" deg\"<\/span>),Location=\"southeast\"<\/span>)<\/span><\/div>\n<\/div>\n
\n
title(\"Sunset Times at Various North Latitudes\"<\/span>)<\/span><\/div>\n<\/div>\n
\n
ylabel(\"Time of Sunset\"<\/span>)<\/span><\/div>\n<\/div>\n
\n
setTimeYTick(gca)<\/span><\/div>\n<\/div>\n
\n
hold off<\/span><\/span><\/div>\n
\n
\n
<\/div>\n
<\/div>\n
<\/div>\n
<\/div>\n<\/div>\n<\/div>\n<\/div>\n<\/div>\n
As expected, day length varies considerably as you get farther away from the equator.<\/div>\n
\n
\n
legend off<\/span><\/span><\/div>\n<\/div>\n
\n
set(findobj(gcf,\"Type\"<\/span>,\"line\"<\/span>),LineWidth=1)<\/span><\/div>\n<\/div>\n
\n
h = line(crepusculusDates,crepusculusTimes,LineWidth=2,Color=\"red\"<\/span>,Marker=\"o\"<\/span>);<\/span><\/div>\n<\/div>\n
\n
legend(h,\"Crepusculus\"<\/span>,Location=\"southeast\"<\/span>)<\/span><\/div>\n<\/div>\n
\n
title(\"Date of the Earliest Sunset Varies with Latitude\"<\/span>)<\/span><\/div>\n
\n
\n
<\/div>\n
<\/div>\n
<\/div>\n
<\/div>\n<\/div>\n<\/div>\n<\/div>\n<\/div>\n
Crepusculus depends on your latitude! The earliest sunset skews earlier in the year as you get closer to the equator.<\/div>\n
Once again, let's display some of this information on the face of a clock. This time, instead of an animation, we'll do a static clock face where the color of the hour hand depends on the month of the year.<\/div>\n
\n
\n
d = datetime(2022,6,1:214);<\/span><\/div>\n<\/div>\n
\n
lat = 42.3;<\/span><\/div>\n<\/div>\n
\n
lon = -71.4;<\/span><\/div>\n<\/div>\n
\n
[~,sset] = sunrise_dt(lat,lon,d);<\/span><\/div>\n<\/div>\n
\n
t = hour(sset)+minute(sset)\/60+second(sset)\/3600;<\/span><\/div>\n<\/div>\n
\n
[~,ixEarliest] = min(t);<\/span><\/div>\n<\/div>\n
\n
[~,ixLatest] = max(t);<\/span><\/div>\n<\/div>\n
\n
d = d(ixLatest:ixEarliest);<\/span><\/div>\n<\/div>\n
\n
t = t(ixLatest:ixEarliest);<\/span><\/div>\n<\/div>\n
\n
sset = sset(ixLatest:ixEarliest);<\/span><\/div>\n<\/div>\n
\n
<\/div>\n<\/div>\n
\n
clf<\/span><\/div>\n<\/div>\n
\n
drawClockFace<\/span><\/div>\n<\/div>\n
\n
gammaHr = pi\/2 - 2*pi*mod(t,12)\/12;<\/span><\/div>\n<\/div>\n
\n
<\/div>\n<\/div>\n
\n
% colorMap = prism(7);<\/span><\/span><\/div>\n<\/div>\n
\n
colorMap = get(gca,\"ColorOrder\"<\/span>);<\/span><\/div>\n<\/div>\n
\n
<\/div>\n<\/div>\n
\n
colorIxList = [1 1 1 1 1 1 2 3 4 5 6 7];<\/span><\/div>\n<\/div>\n
\n
monthList = [\"Jan\"<\/span>,\"Feb\"<\/span>,\"Mar\"<\/span>,\"Apr\"<\/span>,\"May\"<\/span>,\"Jun\"<\/span>,\"Jul\"<\/span>,\"Aug\"<\/span>,\"Sep\"<\/span>,\"Oct\"<\/span>,\"Nov\"<\/span>,\"Dec\"<\/span>];<\/span><\/div>\n<\/div>\n
\n
monthList = monthList(unique(month(d)));<\/span><\/div>\n<\/div>\n
\n
h = zeros(12,1);<\/span><\/div>\n<\/div>\n
\n
<\/div>\n<\/div>\n
\n
nMonths = length(unique(month(d)));<\/span><\/div>\n<\/div>\n
\n
r = 0.8;<\/span><\/div>\n<\/div>\n
\n
for <\/span>i = 1:length(gammaHr)<\/span><\/div>\n<\/div>\n
\n
colorIx = colorIxList(month(d(i)));<\/span><\/div>\n<\/div>\n
\n
% colorIx = colorIxList(i);<\/span><\/span><\/div>\n<\/div>\n
\n
color = colorMap(colorIx,:);<\/span><\/div>\n<\/div>\n
\n
h(colorIx) = line([0 r*cos(gammaHr(i))],[0 r*sin(gammaHr(i))], ...<\/span><\/span><\/div>\n<\/div>\n
\n
Color=color, ...<\/span><\/span><\/div>\n<\/div>\n
\n
LineWidth=2);<\/span><\/div>\n<\/div>\n
\n
end<\/span><\/span><\/div>\n<\/div>\n
\n
r = 1.1;<\/span><\/div>\n<\/div>\n
\n
x = r*cos(gammaHr(1));<\/span><\/div>\n<\/div>\n
\n
y = r*sin(gammaHr(1));<\/span><\/div>\n<\/div>\n
\n
line([0 x],[0 y],Color=0.5*[1 1 1])<\/span><\/div>\n<\/div>\n
\n
sset.Format = \"dd-MMM HH:mm\"<\/span>;<\/span><\/div>\n<\/div>\n
\n
text(x,y,sprintf(\"Latest Sunset \\n%s \"<\/span>, string(sset(1))), ...<\/span><\/span><\/div>\n<\/div>\n
\n
Color=0.5*[1 1 1], ...<\/span><\/span><\/div>\n<\/div>\n
\n
VerticalAlignment=\"top\"<\/span>, ...<\/span><\/span><\/div>\n<\/div>\n
\n
HorizontalAlignment = \"right\"<\/span>)<\/span><\/div>\n<\/div>\n
\n
x = r*cos(gammaHr(end));<\/span><\/div>\n<\/div>\n
\n
y = r*sin(gammaHr(end));<\/span><\/div>\n<\/div>\n
\n
line([0 x],[0 y],Color=0.5*[1 1 1])<\/span><\/div>\n<\/div>\n
\n
sset.Format = \"dd-MMM HH:mm\"<\/span>;<\/span><\/div>\n<\/div>\n
\n
text(x,y,sprintf(\" Earliest Sunset \\n %s \"<\/span>, string(sset(end))), ...<\/span><\/span><\/div>\n<\/div>\n
\n
Color=0.5*[1 1 1], ...<\/span><\/span><\/div>\n<\/div>\n
\n
VerticalAlignment=\"top\"<\/span>, ...<\/span><\/span><\/div>\n<\/div>\n
\n
HorizontalAlignment = \"left\"<\/span>)<\/span><\/div>\n<\/div>\n
\n
axis equal<\/span><\/span><\/div>\n<\/div>\n
\n
axis off<\/span><\/span><\/div>\n<\/div>\n
\n
<\/div>\n<\/div>\n
\n
monthList = [\"Jan\"<\/span>,\"Feb\"<\/span>,\"Mar\"<\/span>,\"Apr\"<\/span>,\"May\"<\/span>,\"Jun\"<\/span>,\"Jul\"<\/span>,\"Aug\"<\/span>,\"Sep\"<\/span>,\"Oct\"<\/span>,\"Nov\"<\/span>,\"Dec\"<\/span>];<\/span><\/div>\n<\/div>\n
\n
monthList = monthList(unique(month(d)));<\/span><\/div>\n<\/div>\n
\n
h(h==0) = [];<\/span><\/div>\n<\/div>\n
\n
legend(h,monthList, ...<\/span><\/span><\/div>\n<\/div>\n
\n
Location=\"NortheastOutside\"<\/span>)<\/span><\/div>\n
\n
\n
<\/div>\n
<\/div>\n
<\/div>\n
<\/div>\n<\/div>\n<\/div>\n<\/div>\n<\/div>\n
This graph shows roughly half the year, from the latest sunset in June to the earliest sunset in December.<\/div>\n
Finally, here's an animation of the sunset visualization shown above as it changes with latitude. Once again, we see that equatorial days don't vary much in length.<\/div>\n
\"animate_clock.gif\"<\/div>\n
Which would you rather have: the steady day length of the tropics, or gloriously long days of summer at the expense of short winter days?<\/div>\n
\n
\n
function <\/span>setTimeYTick(ax)<\/span><\/div>\n<\/div>\n
\n
h = hours(get(ax,\"YTick\"<\/span>))-hours(12);<\/span><\/div>\n<\/div>\n
\n
h.Format = \"hh:mm\"<\/span>;<\/span><\/div>\n<\/div>\n
\n
hStr = replace(string(h), textBoundary+\"0\"<\/span>, \"\"<\/span>) + \" pm\"<\/span>;<\/span><\/div>\n<\/div>\n
\n
set(gca,YTickLabels=hStr, ...<\/span><\/span><\/div>\n<\/div>\n
\n
YDir=\"reverse\"<\/span>)<\/span><\/div>\n<\/div>\n
\n
end<\/span><\/span><\/div>\n<\/div>\n<\/div>\n<\/div>\n","protected":false},"excerpt":{"rendered":"
<\/div>\n

\nAs the days get shorter and drearier here in New England, I get hungry for sunlight. This makes me somewhat obsessed with when the sun sets. But unlike our ancient ancestors, I don't need to worry... read more >><\/a><\/p>\n","protected":false},"author":69,"featured_media":0,"comment_status":"open","ping_status":"closed","sticky":false,"template":"","format":"standard","meta":[],"categories":[1],"tags":[],"_links":{"self":[{"href":"https:\/\/blogs.mathworks.com\/community\/wp-json\/wp\/v2\/posts\/8735"}],"collection":[{"href":"https:\/\/blogs.mathworks.com\/community\/wp-json\/wp\/v2\/posts"}],"about":[{"href":"https:\/\/blogs.mathworks.com\/community\/wp-json\/wp\/v2\/types\/post"}],"author":[{"embeddable":true,"href":"https:\/\/blogs.mathworks.com\/community\/wp-json\/wp\/v2\/users\/69"}],"replies":[{"embeddable":true,"href":"https:\/\/blogs.mathworks.com\/community\/wp-json\/wp\/v2\/comments?post=8735"}],"version-history":[{"count":4,"href":"https:\/\/blogs.mathworks.com\/community\/wp-json\/wp\/v2\/posts\/8735\/revisions"}],"predecessor-version":[{"id":8768,"href":"https:\/\/blogs.mathworks.com\/community\/wp-json\/wp\/v2\/posts\/8735\/revisions\/8768"}],"wp:attachment":[{"href":"https:\/\/blogs.mathworks.com\/community\/wp-json\/wp\/v2\/media?parent=8735"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/community\/wp-json\/wp\/v2\/categories?post=8735"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/community\/wp-json\/wp\/v2\/tags?post=8735"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}