During the course of one year, I photographed 95 sunrises along the east bank of the Columbia River from my home on the west bank in Richland, WA (46.3° N latitude). I then calculated the seasonal phenomena listed in the title with the intention of explaining the variation of the photographed sunrises. The calculations use a simplified model of the Sun-Earth system and employ rotation matrices to predict the path of the Sun, as observed at any location in the northern hemisphere, throughout the year. These predictions are in good agreement with those listed by NOAA and also with the photographic data. The analysis presented here provides a novel way to calculate and understand the seasonal variations of visible sunlight.

Several years ago, I began a retirement project by taking daily photographs of sunrises visible from my home on the Columbia River in Richland, WA, during an entire year. Figure 1 shows a sample of those photographs taken in September and October 2018. The photos and their analysis were published in Ref. 1. More recently,2 the changing number of daylight hours and the elevation of the Sun at noon were studied using geometric diagrams.

Fig. 1.

(Color online) Sunrise photographs from September 15, 2018 to October 30, 2018. The times listed in the figure are sunrise times in PDT. Because this is an agricultural region, few structures appear in these photos, which span 0.44 miles (0.71 km).

Fig. 1.

(Color online) Sunrise photographs from September 15, 2018 to October 30, 2018. The times listed in the figure are sunrise times in PDT. Because this is an agricultural region, few structures appear in these photos, which span 0.44 miles (0.71 km).

Close modal

Other authors treat the same phenomena. The text by Meeus3 discusses corrections to the basic theory due to small effects, such as the refraction of sunlight in the Earth's atmosphere, that affect the number of daylight hours. References 4 and 5 are popular texts for amateur astronomers. Levine6 used vector analysis to understand the duration of daylight and the position of the rising and setting Sun.

In contrast to the above treatments, this paper studies the apparent motion of the Sun by employing rotation matrices. Using a simplified model of the solar system, discussed in Sec. II, this approach sacrifices some accuracy to provide physics students and instructors with an easier analysis and deeper understanding of the phenomena. Reference 2 provides a good background to this approach.

Our simplified model starts with a heliocentric coordinate system: The Sun is located at the origin of the orbital plane (OP) coordinate system, in which a spherical Earth circles the Sun in the xOP-yOP plane. One hemisphere of the Earth is illuminated by sunlight, represented by a unit vector S that is parallel to the line joining the centers of the Sun and Earth; the other hemisphere is in darkness. Most of the calculations in this paper will use the latitude of the author's home, which is 46.3° N.

Rather than using this heliocentric system, it is easier to proceed using the geocentric coordinate system shown in Fig. 2 with the Earth at its origin and the Sun revolving around it in a circular path in the OP plane. The vector S, when directed toward the Earth's center, makes an angle $κ$ with the negative xOP axis, where $κ$ ranges from 0° to 360° during one year. The asterisks in Fig. 2 indicate the circular boundary between the light and dark hemispheres, called the light-dark circle (LDC). As S changes the direction each day, the plane of the LDC changes its orientation, causing the number of daylight hours at a particular location on Earth to change as well.

Fig. 2.

(Color online) Geocentric coordinate system used for calculations. The figure shows a top view of the Earth obtained by looking down the +zOP axis, which is perpendicular to the plane of the paper pointing outward. The Earth's north pole is indicated by a solid black dot, and the south pole, on the under-side of the Earth, by a prominent x. The LDC is marked by asterisks. Note that the LDC plane is always perpendicular to the sunlight vector S.

Fig. 2.

(Color online) Geocentric coordinate system used for calculations. The figure shows a top view of the Earth obtained by looking down the +zOP axis, which is perpendicular to the plane of the paper pointing outward. The Earth's north pole is indicated by a solid black dot, and the south pole, on the under-side of the Earth, by a prominent x. The LDC is marked by asterisks. Note that the LDC plane is always perpendicular to the sunlight vector S.

Close modal
Instead of using the calendar date, we shall specify the day number (DN) or the angle $κ$. For December 21 (winter solstice), κ = 0° and DN = 1. While the angle κ continually increases, the simplified model assumes that it remains constant during each day and increases abruptly by $360 ° / 365 = 0.986 °$ at the beginning of each new day
(1)
On December 21, $κ=0,$ S lies along the xOP axis, and the LDC lies in the yOPzOP plane. Points on the LDC (asterisks in Fig. 2) are given in terms of OP coordinates by
(2)
(3)
(4)
where η ranges from 0° to 360° relative to the +yOP axis.

When $κ≠0$, it is helpful to construct the vector F that points from the Earth's center to a specific location on the LDC. It has a horizontal component $F h= R E cos η$ and a vertical component $F v= R E sin η$ that lie in the plane of the LDC. Let us imagine that the LDC is rotated from its initial orientation $κ=0$ to the orientation shown in Fig. 2. In that case, the vertical component $F v$ remains unchanged but the horizontal component $F h$ now has components along the –xOP axis and the yOP axis. The equations for the LDC in the OP-coordinate system are now given by

(5)
(5a)
(5b)
(5c)

We now introduce the E-coordinate system shown in Fig. 3, where the zE and xE axes are rotated by 23.5° relative to the OP-coordinate axes, and zE is aligned with the Earth's N-S axis. The xE and yE axes lie in the planet's equatorial plane. Note that the Earth rotates with respect to this coordinate system; we will later introduce a coordinate system (called Q) that is fixed with respect to points on the Earth's surface.

Fig. 3.

(Color online) The E-coordinate system is obtained by rotating the OP system by 23.5° about the yOP axis, which is perpendicular to the paper, pointing inward. A parallel beam of sunlight, indicated by the unit vector S, illuminates one hemisphere of Earth while the other is in darkness. On December 21, the north pole is pointed away from the Sun and is in complete darkness, while the south pole has 24 h of daylight, and points on the equator (EQ) have 12 h of daylight. Over the course of the year, the vector S rotates out of the plane of the page as the Sun circles the Earth in the xOP–yOP plane.

Fig. 3.

(Color online) The E-coordinate system is obtained by rotating the OP system by 23.5° about the yOP axis, which is perpendicular to the paper, pointing inward. A parallel beam of sunlight, indicated by the unit vector S, illuminates one hemisphere of Earth while the other is in darkness. On December 21, the north pole is pointed away from the Sun and is in complete darkness, while the south pole has 24 h of daylight, and points on the equator (EQ) have 12 h of daylight. Over the course of the year, the vector S rotates out of the plane of the page as the Sun circles the Earth in the xOP–yOP plane.

Close modal

The locations along the LDC that are visible in this figure are all experiencing sunset at the moment shown, since an observer rotating counter-clockwise with the Earth is travelling from sunlight into darkness. On the opposite side of the planet, hidden from view in Fig. 3, locations along the LDC are experiencing sunrise since they are traveling from darkness into sunlight. To find the number of hours of daylight and to understand the apparent shift in location and direction of sunrise shown in Fig. 1, we need to determine the xE and yE coordinates of the sunrise and sunset points for any day of the year and for any latitude L.

The E-coordinates of the LDC can be found from the OP-coordinates using the rotation matrices discussed in the  Appendix. For a rotation of θ = 23.5° about the $y OP$ axis, Eq. (A4) implies
(6)
The  Appendix also discusses the appropriate matrices for rotations about the x (Eq. (A6)) and z axes (Eq. (A8)).
In Fig. 2, the sunlight vector S  always lies in the xOPyOP plane and has OP coordinates given by $cos κ , sin κ , 0$. Using Eq. (6) with $θ=23.5°$ (the Earth's tilt relative to its orbital plane), the components of S can be transformed from the OP- to the E-coordinate system, with the following result:
(7)
This will be important when we discuss how an observer on Earth views the Sun. A helpful definition, which we shall use shortly, is given by the vector addition of SxE and SyE,
(8)

Sxy is the component of S in the xEyE plane. In a three-dimensional diagram, vector S is parallel to the line joining the center of the Sun to the center of the planet.

In Fig. 4, the 46.3° N latitude circle and the LDC are plotted in E-coordinates for February 20 (DN = 62, κ = 60°). The intersections between the LDC and the latitude circle correspond to sunrise and sunset locations. For κ = 60°, the OP coordinates of the LDC are found from Eqs. (5): $x O P , y OP , z OP = R E( − 3 / 2 cos η , 1 2 cos η , sin η)$. These are then substituted into Eq. (6) to find their E components

(9)
(9a)
(9b)
(9c)

Fig. 4.

Three-dimensional plots of the LDC (asterisks) and 46.3° N latitude circle for February 20, when $κ=60°$. Part of the latitude circle is hidden from view by the darkly shaded LDC plane.

Fig. 4.

Three-dimensional plots of the LDC (asterisks) and 46.3° N latitude circle for February 20, when $κ=60°$. Part of the latitude circle is hidden from view by the darkly shaded LDC plane.

Close modal

Intersections with the 46.3° N latitude circle are easily found, because for that circle zE/RE = sin(46.3°) = 0.7230. A simple way to proceed is to search for solutions to Eqs. (9) using a spreadsheet. For $κ=60°$, the results are shown in Table I. Clearly, $z E / R E≅0.7230$ when η ≅ 68° and 153°. If more accuracy is desired, these two angles can be used as initial estimates to locate the roots of Eq. (9c) numerically, with the results η = 68.2° and 153.1°. Evaluating Eqs. (9a) and (9b), the E-components are found to be (–0.6652, 0.1857, 0.7230)RE and (0.5278,–0.4463, 0.7230)RE for sunrise and sunset, respectively. The sunrise and sunset locations for any other latitude L and angle κ can be found in the same way. A second method using algebra to find the E coordinates of sunrise and sunset is described in the supplementary material.9

Table I.

Angle η in column 1 defines the location of a point on the LDC. Using Eqs. (5) and (9), the components of the LDC for each value of η on February 20 ($κ=60°)$ are given in the OP and E systems. All coordinates are shown as fractions of $R E$, the Earth's radius.

LDC Point on LDC in OP coordinates Point on LDC in E coordinates
η (°) xOP yOP zOP xE yE zE
67  –0.3384  0.1954  0.9205  –0.6774  0.1954  0.7092
68  –0.3244  0.1873  0.9272  –0.6672  0.1873  0.7209
69  –0.3104  0.1792  0.9336  –0.6569  0.1792  0.7324
⋯  ⋯  ⋯  ⋯  ⋯  ⋯  ⋯
152  0.7647  –0.4415  0.4695  0.5140  –0.4415  0.7354
153  0.7716  –0.4455  0.4540  0.5266  –0.4455  0.7240
154  0.7784  –0.4494  0.4384  0.5390  –0.4494  0.7124
LDC Point on LDC in OP coordinates Point on LDC in E coordinates
η (°) xOP yOP zOP xE yE zE
67  –0.3384  0.1954  0.9205  –0.6774  0.1954  0.7092
68  –0.3244  0.1873  0.9272  –0.6672  0.1873  0.7209
69  –0.3104  0.1792  0.9336  –0.6569  0.1792  0.7324
⋯  ⋯  ⋯  ⋯  ⋯  ⋯  ⋯
152  0.7647  –0.4415  0.4695  0.5140  –0.4415  0.7354
153  0.7716  –0.4455  0.4540  0.5266  –0.4455  0.7240
154  0.7784  –0.4494  0.4384  0.5390  –0.4494  0.7124
Using the information above, the hours of daylight can be easily calculated. Let vector $S R →$ ($S S →$) be the location of the sunrise (sunset) point relative to the center of the 46.3° latitude circle, i.e., $S R →= − 0.6652 , 0.1857 R E$ and $S S →$ = (0.5279, −0.4459)RE. Note that $S S →= S R →=0.6909 R E= R E cos 46.3 ° .$ The angle β between these two vectors is found using the dot product
(10)
yielding $β=155.4°$. The number of hours of daylight is then found from
(11)

Figures 5(a)–5(h) show the distribution of sunlight and darkness in the northern hemisphere for eight dates spaced throughout the year. In each figure, the $x E , y E$ coordinates of the LDC (asterisks) are plotted for $0≤η≤180°$, i.e., points in the northern hemisphere, using Eqs. (6). In the upper left corner for each figure, the angle κ is given, followed by the date and the day number (DN). In the upper right corner, the angle φ (defined above) and the number of daylight hours are shown. The vector Sxy, defined in Eq. (8), is drawn for each figure as the perpendicular bisector of the line connecting the sunrise and sunset points.

Fig. 5.

(Color online) Top view of the Earth looking down on the North Pole (N), for eight days throughout the year. The circles, in decreasing diameters, represent the Equator, 46.3° N latitude, and the Arctic Circle (66.5° N latitude). Asterisks denote the LDC, which is the boundary between sunlight and darkness. Sxy, defined in Eq. (8), points to the location of noon on the latitude circle. The upper left of each figure lists the angle κ, the date, and the day number (DN). Angle φ and the hours of daylight are shown in the upper right of each figure.

Fig. 5.

(Color online) Top view of the Earth looking down on the North Pole (N), for eight days throughout the year. The circles, in decreasing diameters, represent the Equator, 46.3° N latitude, and the Arctic Circle (66.5° N latitude). Asterisks denote the LDC, which is the boundary between sunlight and darkness. Sxy, defined in Eq. (8), points to the location of noon on the latitude circle. The upper left of each figure lists the angle κ, the date, and the day number (DN). Angle φ and the hours of daylight are shown in the upper right of each figure.

Close modal

Many seasonal features can be understood using Fig. 5, such as: (1) the equator always has 12 h of daylight and 12 h of darkness; (2) for the Arctic Circle ($L=90°–23.5°=66.5°$), there is total darkness on December 21 and total sunlight on June 21; (3) there are 12 h of daylight everywhere on Earth on the equinoxes (March 21 and September 21); (4) in the northern hemisphere, the hours of daylight increase from December 21 (winter solstice) to June 21 (summer solstice), and then decrease until December 21; (5) for any date, as you travel to lower latitudes, the hours of daylight increase. Notice (6) also that from December 21 to January 19, the hours of daylight increase only 0.6 h, while from February 20 to March 21 they increase by 1.6 h. Thus, to an observer on Earth, the Sun appears to travel slowly at the solstices and faster at the equinoxes.

Figure 5 of Ref. 2 displays the number of hours of daylight, calculated as above, vs. the date (or day number), and compares them to information obtained from the National Oceanic and Atmospheric Administration (NOAA) website.10 The precise NOAA values are about 0.2 h longer than those calculated here, a difference of 2.4% on December 21 and 1.3% on June 21. The discrepancy is due to small effects that are not considered here, such as the refraction of sunlight in the Earth's atmosphere.

In what direction does an observer on Earth look to view the Sun? To answer this question, we introduce a new coordinate system Q that rotates with the observer about the planet's N–S axis. As shown in Fig. 6, the + $y Q$ axis points north, the $+ x Q$ axis points east, while the $+ z Q$ axis is in the observer's vertical direction. At noon on December 21, the sunlight vector S lies in the plane of the paper, and a straightforward analysis using plane geometry is possible. For other locations of the observer or other times and days of the year, rotation matrices show what an observer sees while watching the path of the Sun. However, keep Fig. 6 in mind, as it illustrates the essential mathematics for determining the components of S in the observer's reference frame.

Fig. 6.

Orientation of the observer and Q-coordinate system at 46.3° N latitude on December 21 at noon. The vector S shows the direction of sunlight striking the Earth, while the oppositely directed vector SOQ shows the direction of the Sun seen by the observer at local noon. The local north–south direction is tangent to a longitude circle, while the east–west direction is tangent to the latitude circle. The +xQ axis is perpendicular to the paper pointing outward.

Fig. 6.

Orientation of the observer and Q-coordinate system at 46.3° N latitude on December 21 at noon. The vector S shows the direction of sunlight striking the Earth, while the oppositely directed vector SOQ shows the direction of the Sun seen by the observer at local noon. The local north–south direction is tangent to a longitude circle, while the east–west direction is tangent to the latitude circle. The +xQ axis is perpendicular to the paper pointing outward.

Close modal
On December 21, sunlight strikes the Earth at an angle of 23.5° with respect to the +xE axis, so
(12)
which agrees with Eq. (7) for κ = 0. The observer views the Sun in the direction of SOQ, which is antiparallel to S, i.e., SOQ = –S. The dashed line in Fig. 6 is the (local) horizontal axis –yQ, which points due south. In the figure, we see that the sum 46.3° + 23.5° + E = 90°, so the elevation angle E of the Sun is equal to 20.2° and SOQ is given by
(13)
Equation (13) tells us what to expect (for December 21 at noon) when rotation matrices are used.

In Fig. 6, the observer's Q-coordinate system (called the “upper” Q system) has its origin at the location of the observer. However, when using a rotation matrix to transform from one coordinate system to another, the origins of the two systems must coincide. Therefore, we define another Q-coordinate system (called the “lower” Q system) whose axes are parallel to the upper Q system but whose origin is located at the Earth's center. The components of S are identical in the upper and lower Q systems. Similarly, the xQ and yQ values of the LDC are identical in the two systems, but the zQ values differ by the radius of the Earth $R E$. While we visualize events in the upper Q-coordinate system, we employ the rotation matrices using the lower Q system.

Two rotations are required to transform from the E-system to the Q-system. Recall that all coordinate systems have their origin at the center of the Earth. The first is shown in Fig. 7(a) for December 21: The –yE axis is rotated by angle θ about the zE axis so that it passes directly BELOW the location of the observer. This rotation yields an intermediate coordinate system, which we will call the M-coordinate system and label the rotation matrix RME. To obtain the value of θ, we use Fig. 5(a), noting that angle $β=125.9°$ and also that the sunlight is symmetrical about the xE axis. The positions of sunrise and sunset are then $180° ∓ β / 2=(117.05°, 242.95°)$. To rotate the −yE axis to the sunrise point, $θ=117.05−270°=−152.95°,$ where the negative angle denotes a CW rotation. The rotation matrix RME is obtained by substituting θ = −152.95° in Eq. (A8) of the  Appendix.

Fig. 7.

(a) Top view of the 46.3° N latitude circle. A rotation (RME) of –152.95° about the zE axis, which is perpendicular to the paper pointing outward, places the –yM axis below the observer at sunrise. The rotated axes are called the M-coordinate system. The + symbols indicate equally spaced locations between sunrise and sunset. (b) Side view of the Earth showing the second rotation (RQM) in which the +zM axis is rotated by +43.7° (the complement to the 46.3° N latitude) about the +xM axis to pass through the sunrise location. The rotated axes are called the Q-coordinate system. The +xM and +xQ axes are perpendicular to the plane of the paper, pointing outward. The “upper” Q-coordinate system is drawn at the surface of the Earth and has the arrangement of axes shown in Fig. 6.

Fig. 7.

(a) Top view of the 46.3° N latitude circle. A rotation (RME) of –152.95° about the zE axis, which is perpendicular to the paper pointing outward, places the –yM axis below the observer at sunrise. The rotated axes are called the M-coordinate system. The + symbols indicate equally spaced locations between sunrise and sunset. (b) Side view of the Earth showing the second rotation (RQM) in which the +zM axis is rotated by +43.7° (the complement to the 46.3° N latitude) about the +xM axis to pass through the sunrise location. The rotated axes are called the Q-coordinate system. The +xM and +xQ axes are perpendicular to the plane of the paper, pointing outward. The “upper” Q-coordinate system is drawn at the surface of the Earth and has the arrangement of axes shown in Fig. 6.

Close modal
The second rotation is shown in Fig. 7(b), where the plane contains the zM axis, the rotated axis yM, and the sunrise point. After rotating by $+43.7°$ about the xM axis, the zM axis passes through the sunrise point. The rotation matrix RQM is obtained by substituting +43.7° for angle θ in Eq. (A6) of the  Appendix. This yields the Q-coordinate system shown in Fig. 7(b), where the +xQ (and +xM) axes are perpendicular to the plane of the paper, pointing outward. The “upper” Q-system can then be constructed at the location of the observer, with the desired arrangement of axes shown in Fig. 6. For a different latitude L, the rotation angle used in RQM is 90° – L. Of course, the locations of sunrise and sunset would need to be determined to find angle θ. The transformation from E- to Q-coordinates may, thus, be written in matrix notation as $S Q= R QM R ME S E$, and the generalized equation for the direction to the Sun, according to the observer, is
(14)
where L is the observer's latitude, θ is the angular rotation of the –yE axis about the zE axis, 90 – L is the rotation about the xM axis, and κ depends on the day of the year.

We can now use Eq. (14) to find SOQ when the observer views the Sun at sunrise on December 21 (κ = 0). Evaluating Eq. (14) for L = 46.3° and θ = –152.95° yields $SO Q= 0.8167 , − 0.5770 , 0.$ The positive sign of SOQx and the negative sign of SOQy show that the observer is facing southeast (components to the south and east) to view the Sun. Specifically, the angle ϕ between the +xQ axis and the vector SOQ is given by $tan − 1 − 0.5770 / 0.8167 = − 35.2 °$. As expected, the zQ component is zero at sunrise.

At noon on December 21, the observer is midway between sunrise and sunset, or on the –xE axis, as shown in Fig. 7(a). Therefore, the –yE axis must be rotated CW by $θ=−90°$ about the zE axis. The values for L, κ, and S remain the same as those for sunrise. Evaluating Eq. (14) yields $SO Q= 0 , − 0.9385 , 0.3453,$ indicating that the vector points due south. The elevation angle E is the angle between the Sun and the horizon and is given by $E= tan − 1 0.3453 / 0.9385=20.2°.$ This is the value that was found in Fig. 6 using geometry.

#### 1. Path of the Sun on December 21

The apparent path of the Sun on December 21 is obtained by evaluating Eq. (14) throughout the day, as the observer moves from sunrise to sunset. Since $θ=−152.95°$ at sunrise, and $θ=−27.05°$ at sunset (Fig. 7(a)), the value of θ at equally spaced intervals during the day is easily found. A three-dimensional plot of (SOQx, SOQy, SOQz) for 17 times (shown as + signs) between sunrise and sunset is displayed in Fig. 8(a).

Fig. 8.

The observer, located at the origin, sees these three-dimensional paths of the Sun, which are obtained by plotting the vector SOQ for (a) December 21, (b) March 21, and (c) June 21, where the noon elevation angles E are 20.2°, 43.7°, and 67.2°, respectively. The bearing angle ϕ is defined in (b). At sunrise, ϕ is (a) – 35.2° (southeast), (b) 0° (due east), and (c) +35.2° (northeast). The time between each asterisk is approximately 30 min, so the largest number of asterisks appears for June 21 (summer solstice).

Fig. 8.

The observer, located at the origin, sees these three-dimensional paths of the Sun, which are obtained by plotting the vector SOQ for (a) December 21, (b) March 21, and (c) June 21, where the noon elevation angles E are 20.2°, 43.7°, and 67.2°, respectively. The bearing angle ϕ is defined in (b). At sunrise, ϕ is (a) – 35.2° (southeast), (b) 0° (due east), and (c) +35.2° (northeast). The time between each asterisk is approximately 30 min, so the largest number of asterisks appears for June 21 (summer solstice).

Close modal

#### 2. Path of the Sun on March 21 and June 21

The calculations of SOQ for March 21 and June 21 (Figs. 8(b) and 8(c)) are similar to those for December 21, but since the sunrise and sunset locations are different for each date, the value of θ in Eq. (14) also changes. In Figs. 8(b) and 8(c), the path of the Sun is shown for these dates. Table II summarizes the calculations for the vector SOQ on the three dates shown in Fig. 8. Elevations obtained from the NOAA website10 for the three dates are as follows: 20.31°, 44.17°, and 67.15°. The values obtained here are, thus, in excellent agreement with those computed by NOAA.

Table II.

Results of calculations for the vector SOQ. The second column is calculated from Eq. (14). The third column lists the sunrise angle ϕ or the elevation angle at noon. The term “south-east” indicates the vector SOQ has components in the south and east directions, and so on.

Date and time SOQ Angle ϕ or E
Dec 21 sunrise  (0.8167, −0.5770, 0)  ϕ = −35.2°, south-east
Dec 21 noon  (0, −0.9385, 0.3453)  E = 20.2°, due south
March 21 sunrise  (1.0000, 0, 0)  ϕ = 0, due east
March 21 noon  (0.−0.7230, 0.6909)  E = 43.7°, due south
June 21 sunrise  (0.8167, 0.5770, 0)  ϕ = 35.2°, north-east
June 21 noon  (0, −0.3875, 0.9210)  E = 67.2°, due south
Date and time SOQ Angle ϕ or E
Dec 21 sunrise  (0.8167, −0.5770, 0)  ϕ = −35.2°, south-east
Dec 21 noon  (0, −0.9385, 0.3453)  E = 20.2°, due south
March 21 sunrise  (1.0000, 0, 0)  ϕ = 0, due east
March 21 noon  (0.−0.7230, 0.6909)  E = 43.7°, due south
June 21 sunrise  (0.8167, 0.5770, 0)  ϕ = 35.2°, north-east
June 21 noon  (0, −0.3875, 0.9210)  E = 67.2°, due south

In this section, we shall compare sunrise photographs taken by the author, several of which are shown in Fig. 1, to the theory developed above. A simplified bird's eye view of the area near the author's home is shown in Fig. 9. Distances along the horizon were obtained by comparing features on Google Maps to the same features visible in the photographs, as discussed in Ref. 1. Figure 9 depicts seven sunrises chosen about one month apart. The asterisks show where the Sun appears at dawn along the bluff of the river. The dotted lines are drawn using the calculated values of ϕ for these dates. They converge at point C, which is ≈1320 m west of line AB and is taken to be the origin of the Q-coordinate system. The distance of each sunrise from that on the equinox (March 22) at point D is given by $D photo=DC tan ϕ$.

Fig. 9.

(Color online) Top view of the area surrounding the Columbia River, drawn to scale. The distance AB = 1860 m, AD = DB = 930 m, DE= 110 m, EF = 855 m, and DC = 1320 m. Point G, indicated by the x-symbol, is the location of the author's home, where the photographs were taken. It is 45 m from the bank of the Columbia River and 74 m north of point C. The number in parentheses, below each date, indicates the DN and the calculated angle ϕ at sunrise.

Fig. 9.

(Color online) Top view of the area surrounding the Columbia River, drawn to scale. The distance AB = 1860 m, AD = DB = 930 m, DE= 110 m, EF = 855 m, and DC = 1320 m. Point G, indicated by the x-symbol, is the location of the author's home, where the photographs were taken. It is 45 m from the bank of the Columbia River and 74 m north of point C. The number in parentheses, below each date, indicates the DN and the calculated angle ϕ at sunrise.

Close modal

To explore this further, we analyzed 55 sunrise photos taken between June 22 and December 21, measuring the distance of sunrise in each photo from that on the equinox (point D, March 22).

Using the calculated value of ϕ for each date, a linear least-squares fit11 to the data (Fig. 10) yields
(15)
Fig. 10.

Distance between the sunrise location on a given date and that at the equinox on March 21 is plotted on the vertical axis. For each date, angle ϕ is calculated and tan (ϕ) plotted on the horizontal axis. The solid line and the dashed line are the least squares fits to Eqs. (15) and (16). Negative values on the vertical axis were obtained from December 21 to March 21, while positive values were obtained from March 21 to June 21.

Fig. 10.

Distance between the sunrise location on a given date and that at the equinox on March 21 is plotted on the vertical axis. For each date, angle ϕ is calculated and tan (ϕ) plotted on the horizontal axis. The solid line and the dashed line are the least squares fits to Eqs. (15) and (16). Negative values on the vertical axis were obtained from December 21 to March 21, while positive values were obtained from March 21 to June 21.

Close modal
When the intercept is forced to be zero, the best fit straight line is given as follows:
(16)
The slope of either line is very close to DC, demonstrating that the photographic data are self-consistent and in good agreement with the mathematical analysis. Ideally, point G in Fig. 9 (the camera location) would be coincident with point C. The fact that it is not suggests that the geometry of Fig. 9 is oversimplified.

When I began this project, my first calculations predicted 8.4 h of daylight for December 21, while the local newspaper reported a time of 8.63 h. This suggested that the simplified model of the Earth-Sun system presented here, and the introduction of the light-dark circle (LDC), could provide a straightforward explanation of the seasonal variations of sunlight. Using rotation matrices rather than vector analysis, diagrams such as Fig. 5 were generated to explain the duration of daylight in the northern hemisphere throughout the year. I have not seen diagrams similar to Fig. 5 in any textbook, so perhaps they introduce a new way to visualize seasonal phenomena. Further analysis using rotation matrices allows one to explain the observed path of the Sun from any location on Earth, and the results are in good agreement with sunrise photographs taken throughout the year.

The author has no conflicts of interest to disclose.

Figure 11 shows two coordinate systems: System B, represented by axes with solid lines, in which vector r has known components (xB, yB, zB); and system C, with dashed lines, in which r has unknown components (xC, yC, zC). In Fig. 11(a), system C is obtained by rotating the XB and ZB axes counterclockwise by a (positive) angle θ. Because the YB and YC axes are parallel, $y B= y C$. The components of zB along the XC and ZC axes, shown with thick dashed lines, are xC = −sin θ zB and zC = cos θ zB. The components of xB along the XC and ZC axes, shown with prominent dotted lines, are xC = cos θ xB and zC = sin θ xB. Summing the components on the XC and ZC axes, we find the following:
(A1)
(A2)
(A3)
Expressing these equations in the matrix form yields the following:
(A4)
or
(A5)
where rC and rB are column vectors and RCB is the 3 × 3 rotation matrix.
Fig. 11.

(a) Rotation of the XB and ZB axes by a positive angle θ about the +YB axis, which is perpendicular to the paper pointing inward. (b) Rotation of the YB and ZB axes by a positive (counterclockwise) angle θ about the XB axis, which points outward from the paper. (c) Rotation of the XB and YB axes by a positive angle θ about the ZB axis, which points outward from the paper. If a rotation is opposite to that in a specified figure then θ has a negative value.

Fig. 11.

(a) Rotation of the XB and ZB axes by a positive angle θ about the +YB axis, which is perpendicular to the paper pointing inward. (b) Rotation of the YB and ZB axes by a positive (counterclockwise) angle θ about the XB axis, which points outward from the paper. (c) Rotation of the XB and YB axes by a positive angle θ about the ZB axis, which points outward from the paper. If a rotation is opposite to that in a specified figure then θ has a negative value.

Close modal
In Fig. 11(b), coordinate system C is obtained by rotating system B by angle θ about the XB axis. The results are given in matrix form by
(A6)
or
(A7)
Similarly, the results for rotation by angle θ about the ZB axis, shown in Fig. 11(c), are given by
(A8)
For additional information, the reader is referred to Refs. 7 and 8.
1.
Margaret Stautberg
Greenwood
, “
95 Sunrises along the Columbia riverbank
,” <https://skyandtelescope.org/observing/stargazers-corner/95-Sunrises> accessed on May 28, 2020.
2.
Margaret Stautberg
Greenwood
, “
Geometry and the cause of the seasons: The changing hours of daylight and elevation angle of the Sun at noon
,”
Phys. Teach.
60
,
694
698
(
2022
).
3.
Jean
Meeus
,
Astronomical Algorithms
, 2nd ed. (
Willmann-Bell
,
Richmond
,
1998
).
4.
Peter
Duffett-Smith
and
Jonathan
Zwart
,
, 4th ed. (
Cambridge U. P
.,
Cambridge
,
2017
).
5.
J. L.
Lawrence
,
Celestial Calculations: A Gentle Introduction to Computational Astronomy
(
MIT Press
,
Cambridge
,
2019
).
6.
Jonathan
Levine
, “
Seasonal changes in the apparent position of the Sun as elementary applications of vector operations
,”
Eur. J. Phys.
35
,
065020
(
2014
).
7.
See <https://en.wikipedia.org/wiki/Rotation_matrix#in_three_dimensions> for a detailed discussion of rotation matrices.
8.
Mary L.
Boas
,
Mathematical Methods in the Physical Sciences
, 2nd ed. (
John Wiley
,
New York
,
1983
), p.
410
.
9.
See supplementary material online for a second method to find the E coordinates of sunrise and sunset algebraically.
10.
See <https://gml.noaa.gov/grad/solcalc/> or search for NOAA Solar Calculator using computer browser.
11.
Search the internet for “Physics 305 Hints: Using LINEST in Excel” for information about error analysis.