Solar Angle

The Earth

Your browser does not support the canvas element.

Ever been stuck in an office where you couldn't see outside? Or wanted to know if the sun was up somewhere else on the globe? You might find this little page useful.

History

While working at AFGWC on a radiation package in the '80s for a climate / weather forecast model, I wanted something that could give me the angle of the sun from a grid point on the planet. This algorithm is the result.

It's not terribly fast, but the original was extremely accurate for the 50 years' worth of astronomical data I could find. The requirement wasn't for accuracy, but speed. Oh, well. Over the decades I lost the code and any notes I'd made, and this is the recreation from memory.

Method

The dot-product of two normalized vectors in cartesian space is the cosine of the angle of separation (see the box below). The apparent latitude of the sun doesn't change much over the course of a day so I use an approximation. I simply did the dot-product and solved for the unknown: the longitude of the sun. This is all observer-centric (no offset from the center of the planet - though I think the original used a translation of coordinates).

The two constants I calculated for the two waves for the Equation Of Time are from memory. The only part I vaguely remembered was the math.

Derivation of the formula for θs:
s⃗ * o⃗ = cos( δ )
s⃗ = ( xs,ys,zs )     o⃗ = ( xo,yo,zo )
x = r * cos( φ ) * cos( θ )
y = r * cos( φ ) * sin( θ )
z = r * sin( φ )
normalize: r ≡ 1

cos( δ ) = xs * xo +ys * yo +zs * zo
cos( δ ) = a * cos( θs ) + b * sin( θs ) + c
a = cos( φo ) * cos( θo ) * cos( φs )
b = cos( φo ) * sin( θo ) * cos( φs ) which drops if θo ≡ 0
c = sin( φo ) * sin( φs )

d = a * cos( θs ) + b * sin( θs )
using the double angle identity:
d = e * sin( α + θs )
d = e * cos( α ) * cos( θs ) - e * sin( α ) * sin( θs )
a = e * cos( α )
b = -e * sin( α )
tan( α ) = sin( α ) / cos( α )
tan( α ) = ( -b/e ) / ( a/e ) = -b/a
α = tan-1( -b/a )

a2 = e2 * cos2( α )
b2 = e2 * sin2( α )
a2 + b2 = e2 * ( cos2( α ) +sin2( α ) )
a2 + b2 = e2
e = √( a2 + b2 )
d = e * cos( α + θs )
d/e = cos( α + θs )
α = cos-1( d/e )

θs = cos-1( d/e )

To see what's actually implemented... best to look at the code.
Equation of Time approximation (eccentricity +obliquity):
eot( date ) = effectEccentricity +effectObliquity
effectEccentricity = EOT_DEVIATION_PERIGEE *sin( 2 *π *( date -APSIS_NEAR_PERIHELION) /YEAR_LENGTH_ANOMALISTIC )
effectObliquity = EOT_DEVIATION_SOLSTICE *sin( 4 *π *( date -SOLSTICE_EVENT) /YEAR_LENGTH_TROPICAL )
EOT_DEVIATION_PERIGEE ≊ 480 (from memory)
EOT_DEVIATION_SOLSTICE ≊ 600 (from memory)

Description

Derivation of the double angle identity:
e± i * α = cos( α ) ± i * sin( α ) (by Euler)
e± i * β = cos( β ) ± i * sin( β )
ei * α * ei * β = ei * ( α + β )
ei * ( α + β ) = cos( α + β ) + i * sin( α + β )

cos( α + β ) + i * sin( α + β ) =
( cos( α ) + i * sin( α )) * ( cos( β ) + i * sin( β ))
or
cos( α + β ) = cos( α ) * cos( β ) - sin( α ) * sin( β )
sin( α + β ) = cos( α ) * sin( β ) + sin( α ) * cos( β )

The yellow dot is the position where the sun is directly overhead (or within a half of a degree). The blue and yellow circles are where the sun is at 30° and 60° from the vertical (take note if your skin is sensitive to harsh sunlight).

The three stripes show where it's civil, nautical, or astronomical twilight.

The ECMAScript code starts and creates however many worker threads it needs to have one worker per core on the machine. It then divides the problem up into vertical pixel lines and calculates the angle at each pixel.

The page should be updated every minute, unless the checkbox is unselected.

φs approximation:
φs ≊ 23.44 * cos(( date -SOLSTICE_EVENT ) *2π / YEAR_LENGTH_TROPICAL )
(yes, this could be greatly improved)

A Point On The Earth


latitude
longitude
Noon

Event GMT
Begins Ends
Daylight
Twilight - Civil
Twilight - Nautical
Twilight - Astronomical

To get the times for sunrise, sunset, and the twilights, enter the latitude and longitude into the input fields. The code will display these values and the time of noon (when the sun is over your longitude).

All times are given in GMT. This makes it easy (or less confusing) when inspecting other latitudes and longitudes.

I haven't checked the accuracy of any of this at present, and my constants for my version of the Equation of Time are only approximations. And since I only work on Linux, and use Firefox, I suspect this might not work on other environments.