## Introduction

The rotation of the Earth and its orbit around the Sun are the result of many complex factors, so the mathematical equations required to accurately model them are not trivial. Firstly, the shape of the Earth’s orbital path around the Sun varies from close to circular to being more elliptical every 100,000 years or so. Secondly, the tilt angle of its rotational axis relative to the orbital plane varies periodically from around 22° to 25° and back over 41,000 years. Plus, the direction of the tilt axis itself slowly rotates relative to the Sun every 23,000 years or so. Even tidal actions within the oceans and seas act to slightly increase the speed of the Earth’s axial rotation by about 1 second every 463,000 years (ref: 1). All of these have to be accounted for in order to predict the Sun’s current orbit.

Moreover, what we refer to as one day (24 hours) is the time it takes for the Sun to pass over the same line of longitude on each sequential rotation of the Earth. However, the Earth has to rotate around its axis about 361° before the Sun is back over the same line of longitude, given that the Earth has also completed close to 1° of its orbit around the Sun during that same period (it takes 365.25 days to travel the full 360°, so in 1 day it orbits 0.985626283367556°) (read more). Thus, the equations also have to compensate for the fact that our perception of time is based on a unit that is around 4 minutes more than the actual time required for a full 360° sidereal rotation.

Also, the most commonly required information from any solar position calculation is the relative azimuth and altitude of the Sun in the sky when viewed from a particular geographic location. Thus, in addition to all of the celestial calculations required to locate the Earth’s position and rotation in space relative to the Sun, we also have to reorient the perspective so that the final trigonometry is taken relative to a horizontal plane that is tangential to that location on the Earth’s surface, as shown in Figure 1.

## Earth-Sun Relationship

To help visually understand the orbital relationship of the Earth around the Sun, you are referred to our Earth and Sun online web app. Use it to interactively investigate and experiment with different dates, times and locations. You can also animate the Earth’s orbit and rotation, as well as highlighting the significant latitudes and times of day to see how and why they are relevant.

## Algorithms and Calculations

To reiterate, even calculating just a single solar position involves modelling the Earth’s orbit around the Sun and accounting for factors such as elliptical eccentricity, axial obliquity, right ascension and declination, as well as apparent and mean solar time. All of these involve complex trigonometry, which makes each calculation quite computationally intensive.

If you want another solar position the following day, all of these calculations must be done again. However, if you want one at a different time on the same day, then actually only some of them really need to be done again as you can either resuse or interpolate between some of the previously calculated values with the same level of accuracy. Similar situations occur in many other solar radiation and sky distribution calculations.

Thus, the design and structure of the solar calculations in the Performative Design libraries are almost entirely governed by the need to recycle and reuse any previously done calculations in the smartest and most optimised way possible. Because each recalculation is so computationally intensive, and most solar analysis require many thousands if not millions of them, any opportunity to cache, interpolate or simplify can have enormous speed benefits.

For example, a typical annual cumulative sky calculation will require 43,800 solar positions (assuming 6 minute increments over an average of 12 hours per day for 365 days) as well as 3,622,260 sky radiance calculations (assuming 827 sky segments in 1 hour increments over 12 hours per day for 365 days), with each single one involving more than a hundred arithmetic operations and between 12 and 26 trigonometric function calls. In fact, the most realistic and best method for interpolating sky conditions also samples the Sun position and resulting sky distribution in 6 minute intervals, which actually requires closer to 36,222,600 individual sky radiance calculations over a typical year.

Thankfully, the calculations themselves offer many opportunities for optimisation and the author has been working on them for nearly 30 years in all sorts of applications, infrastructures and programming environments. To the extent that the JavaScript version of the solar calculations used by some of our interactive web applications allows them to generate real-time 3D annual Sun-path diagrams together with equally dynamic full sky radiance distributions, even on a mid-range mobile phone or tablet.

Unfortunately, the algorithms implemented in the Performative Design libraries and online web apps are proprietry and not generally available. Whilst they are all based on the freely available astronomical algorithms published by Jean Meeus (ref: 2), they include extensive contextual optimisations that only really make sense as part of the overall library. All of the knowledge, experience and thousands of hours of iterative work and validation that has gone into these optimisations has to be paid for, and that is by the professional and commercial edge that it gives to the author’s own applications and online tools.

## Available Equations and Code

Thankfully there are lots of different implementations of solar position algorithms available on the web. Some are highly simplified and relatively quick to calculate, but which may not be sufficiently accurate for more demanding applications. Others are much more accurate and, as a consequence, more complex and computationally intensive. The following are some links to publically available solar position calculations and equations, separated into simplified and detailed categories.

#### Simplified Equations

One of the best sets of simplified equations known to the author are the General Solar Position Calculations provided by the NOAA Global Monitoring Division. The linked PDF file lays out each algorithmic step quite clearly and in a way that can be very easily translated into code in pretty well any programing language.

Another interesting set of simplified but reasonably accurate equations was developed by Plataforma Solar de Almeria (PSA) for use by low-power microprocessors in the controllers of individual solar concentrating reflectors (ref: 3). These have to be sufficiently accurate for each reflector in a solar farm to focus their beams effectively onto the central furnace, whilst not being overly taxing on limited computational resources. These equations are available without license within their published paper as well as downloadable source code in C++.

#### Detailed Equations

The Position of the Sun topic on Wikipedia is a pretty good place to start when trying to understand in more detail the fundamental steps involved. However, translating the equations given there into code is not particularly straightforward.

A much more detailed and extensive explanation is given in the Solar Position Algorithm for Solar Radiation Applications PDF file provided by NREL. This walks through each step in the process and gives all of the definitive equations required. Whilst these equations are slightly more suitable for translation into code, and can be downloaded as source files in C, they are covered by terms that require a license for commercial use and prohibit redistribution.

The NOAA provide the equations used in their well known Solar Position Calculator as downloadable Excel and OpenOffice spreadsheets. Mike Bostock (the creator of D3.js) includes an excellent JavaScript implementation of these equations in some of his solar demonstration projects, such as the Solar Terminator example. To view the solar position code, scroll down to near the end of the example web page. He also makes this same code available via unpkg.com. The code on his example page is available under the GNU General Public License, version 3, but the licensing of the code on unpkg.com is unclear.

It is important to note that the NOAA code and Mike Bostock’s implementation are detailed and accurate, but are not in any way optimised. All of the component values are recalculated on each invocation and several of them are converted from degrees to radians and then back again in the final results. However, for most situations these are likely to be more than sufficient.

## A Brief History

The progenitor for the solar position calculations used in the current Performative Design building analysis library were a set of equations published by the CSIRO (ref: 4), implemented by the author in early 1990. Soon after, a local government agency published more detailed equations with accompanying computer code intended for use throughout Western Australia (W.A. Dept. of Agriculture) (ref: 5) where the author was studying at the time. The author’s supervisor had just completed a research project that validated these equations against historical almanac data. A comprehensive regression analysis had suggested slight modifications to some of the coefficients, so this published code together with the suggested modifications served as the basis for the next iteration in 1993.

The early 2000’s saw significant industry activity in this area, triggered by the release of CIBSE TM33 and ANSI/ASHRAE 140 - both providing tests for validating building analysis software and their solar position calculations. These were part of preparations for implementing the EPBD in Europe and the new ASHRAE 90.1 in the United States. In 2004, NREL published their own Solar Position Algorithms (ref: 6), which then became something of a defacto standard against which any other calculations could be compared.

Part of the author’s ongoing development of Ecotect involved extensive optimisation and continuous validation of these calculations, which included collaboration with NREL as part of a project to integrate them within Revit after the acquisition by Autodesk.

After leaving Autodesk in late 2011, the need to develop an entirely new and unencumbered codebase led to a new implementation based on equations published by Jean Meeus (ref: 2) and influenced by the Earth System Research Laboratory at the US National Oceanic & Atmospheric Administration (NOAA). These were first implemented concurrently in Java and C# to support some Processing applications being developed, and then reimplemented in JavaScript with significant additional optimisations to support the author’s online web applications that followed.

## References

1. Muneer, T. (2004), Solar Radiation and Daylight Models. Elsevier, ISBN 0 7506 5974 2. (View online)

2. Meeus, J.H. (1991), Astronomical Algorithms. Willmann-Bell, Incorporated. (View online)

3. Blanco-Muriel, M., Alarcón-Padilla, D.C., López-Moratalla, T., and Lara-Coira, M.Í. (2001), Computing the solar vector, Solar Energy, vol. 70, no. 5, pp. 431 - 441.

4. Spencer, J.W. (1965), Calculation of Solar Position for Building Purposes, CSIRO Division of Building Research Technical Paper No. 14.

5. Roderick, M.L. (1992), Methods for calculating solar position and day length including computer programs and subroutines. Department of Agriculture and Food, Western Australia, Perth. Report 137. (View online)

6. Reda, I., Andreas, A. (2003), Solar Position Algorithm for Solar Radiation Applications. 55 pp.; NREL Report No. TP-560-34302, Revised January 2008. (View online)

7. The United States Naval Observatory (USNO), The Seasons and the Earth’s Orbit, Published by the Astronomical Applications Department. (View online)

8. Mee, N. (2009), How long is a day?, Plus Magazine. (View online)