MATLAB Examples

Building Sundials

If carefully built, a sundial is perfectly capable of keeping accurate time. The problem is that you need a different sundial for every combination of latitude, longitude, and day of the year. Fortunately MATLAB® can make this a manageable problem. By turning this procedure around, you can, if you know the correct time, use MATLAB to tell you which way is north.

Contents

Acknowledgments

The algorithms used here for MATLAB sundial creation were adapted from the Dutch web site De Zonnewijzerkring

Tracking the sun

Imagine a simple sundial that consists only of a pin (called a gnomon) sticking vertically out of a disk. What is the trace made by the tip of the shadow of the gnomon during a given day? Suppose it is February 9th, or day 40 of the year, and we are at latitude 42. Here's how the shadow will change from 8:00 AM to 4:00 PM.

set(gca,'XTick',[],'YTick',[])
box on

latitude = 42;
hour = 9:16;
day = 40;
[x,y] = sunshadow(latitude,day,hour);

line(x,y)

for n = 1:length(hour)
    line(x(n),y(n),'Marker','.')
    line([0 x(n)],[0 y(n)],'Color',0.8*[1 1 1])
    if hour(n) <= 12
        h = hour(n);
        amPmStr = 'AM';
    else
        h = hour(n) - 12;
        amPmStr = 'PM';
    end
    text(x(n),y(n),sprintf('%d\n%s\n', h, amPmStr), ...
        'VerticalAlignment','bottom', ...
        'HorizontalAlignment','left', ...
        'FontSize',6)
end

line(0,0,'Marker','.')
line(0,0,'Marker','o')
text(0,-0.25,'Pin Location', ...
    'FontSize',9, ...
    'HorizontalAlignment','center')
line([-0.5 0.5],-1*[1 1],'Marker','.')
text(0,-1.25,'Pin Length', ...
    'FontSize',9, ...
    'HorizontalAlignment','center')

% Draw the arrow pointing north
arrowX=[0 0 -.1 NaN 0.1 0] + 2;
arrowY=[0 1 0.8 NaN 0.8 1] - 1;
line(arrowX, arrowY, 'LineWidth',1, ...
    'Color','black','LineStyle','-')
text(arrowX(1),arrowY(1),' North', ...
    'FontSize',9)

set(gcf,'Color','white')
axis equal
xlim([-5 5])

This is an accurate clock that has been tailor-made for one day and latitude. How does the sundial change through the year?

Same latitude, January through June

Let's build up all the traces that result for a multiple days of the year at the same latitude.

cla
latitude = 42;
hour = 9:16;
for day = 0:21:171;
    [x,y] = sunshadow(latitude,day,hour);

    line(x,y)
    dayStr = datestr(datenum('1-Jan-2003')+day,'mmm dd');
    text(x(1),y(1),sprintf('%s ', dayStr), ...
        'HorizontalAlignment','right', ...
        'FontSize',7)
    for n = 1:length(hour)
        line(x(n),y(n),'Marker','.')
    end

end

line(0,0,'Marker','.')
line(0,0,'Marker','o')
text(0,-0.25,'Pin Location', ...
    'FontSize',9, ...
    'HorizontalAlignment','center')
line([-0.5 0.5],-1*[1 1],'Marker','.')
text(0,-1.25,'Pin Length', ...
    'FontSize',9, ...
    'HorizontalAlignment','center')

% Draw the arrow pointing north
arrowX=[0 0 -.1 NaN 0.1 0] + 2;
arrowY=[0 1 0.8 NaN 0.8 1] - 1;
line(arrowX, arrowY, 'LineWidth',1, ...
    'Color','black','LineStyle','-')
text(arrowX(1),arrowY(1),' North', ...
    'FontSize',9)

axis equal
xlim([-5 5])

Same time of day across the year: the analemma

This is what you would see if you marked the tip of the shadow at exactly 10 o'clock every day for a year. Notice the interesting figure eight effect. What's going on?

cla
latitude = 42;
hour = 10;
day = 0:5:365;
[x,y] = sunshadow(latitude,day,hour);

line(x,y)

for n = 1:5:length(day)
    dayStr = datestr(datenum('1-Jan-2003')+day(n),'mmm dd');
    if day(n) < 110
        horz = 'right';
    elseif day(n) < 182
        horz = 'left';
    elseif day(n) < 243
        horz = 'right';
    else
        horz = 'left';
    end
    line(x(n),y(n),'Marker','.')
    text(x(n),y(n),sprintf('  %s  ', dayStr), ...
        'HorizontalAlignment',horz, ...
        'FontSize',7)
end

line(0,0,'Marker','.')
line(0,0,'Marker','o')
text(0,-0.25,'Pin Location', ...
    'FontSize',9, ...
    'HorizontalAlignment','center')
line([-0.5 0.5],-1*[1 1],'Marker','.')
text(0,-1.25,'Pin Length', ...
    'FontSize',9, ...
    'HorizontalAlignment','center')

% Draw the arrow pointing north
arrowX=[0 0 -.1 NaN 0.1 0] + 2;
arrowY=[0 1 0.8 NaN 0.8 1] - 1;
line(arrowX, arrowY, 'LineWidth',1, ...
    'Color','black','LineStyle','-')
text(arrowX(1),arrowY(1),' North', ...
    'FontSize',9)

set(gcf,'Color','white')
axis equal
xlim([-3 3])

The figure eight shape is known as the analemma. It results from the fact that, although our clocks pretend the sun moves at exactly the same rate through the sky, the sun actually moves more rapidly during some parts of the year than others. For more information, see http://www.analemma.com/.

Using MATLAB to find north

Suppose you were lost in the woods, equipped with only a computer, a printer, and a copy of MATLAB. How would you orient yourself? Print out this figure, put a pin of the correct length in the specified location, align the pin's shadow with the line shown, and the arrow will point north.

cla
latitude = 42;
% t = now;
t = 731850.5196;
dv = datevec(t);
title(datestr(t))
day = floor(datenum(t)-datenum(dv(1),1,1));
hour = dv(4);
[x,y] = sunshadow(latitude,day,hour);

line(x,y)
line([0 x],[0 y],'Marker','.')

line(0,0,'Marker','.')
line(0,0,'Marker','o')
text(0,-0.25,'Pin Location', ...
    'FontSize',9, ...
    'HorizontalAlignment','center')
line([-0.5 0.5],-1*[1 1],'Marker','.')
text(0,-1.25,'Pin Length', ...
    'FontSize',9, ...
    'HorizontalAlignment','center')

% Draw the arrow pointing north
arrowX=[0 0 -.1 NaN 0.1 0] + 2;
arrowY=[0 1 0.8 NaN 0.8 1] - 1;
line(arrowX, arrowY, 'LineWidth',1, ...
    'Color','black','LineStyle','-')
text(arrowX(1),arrowY(1),' North', ...
    'FontSize',9)

axis equal
xlim([-5 5])