Astro-Algo: Astronomical Algorithms for Clojure

At The Climate Corporation, we have “sprintbaticals”, two-week projects where we can work on something a bit different. This post is about work done by Matt Eckerle during his recent sprintbatical.

Sometimes I want to know where the sun will be at a certain time. I don’t want to know approximately where the sun will be, in a close-enough-for-a-Facebook-app kind of way. I want to know precisely where the sun will be, in a which-picnic-table-at-the-beer-garden-will-be-in-the-sun-at-exactly-5pm-on-my-34th-birthday kind of way. Or, more to the point of my employment at The Climate Corporation, what zenith angle should I use to correct radiation measurements taken at location x and time t? (It turns out that sunlight is important for crop growth. Who knew?) For this sort of endeavor, I keep a copy of Jean Meeus’ Astronomical Algorithms at my desk at all times.

IMG_20140331_150813Using a book on astronomy to compute solar position may sound like overkill, but using various solar calculators over the years has given me a lot of perspective on what is considered “accurate” by non-astronomers. Even the NOAA solar calculator,  which uses methods described in AA, doesn’t bother to consider the difference between dynamical time to universal time. This discrepancy of about 60 seconds seems to be obscured by the fact that they round output times to the nearest minute. And NOAA’s solar calculator is much, much, much better than most solar calculators out there.

Introducing Astro-Algo, Astronomical Algorithms for Clojure

Naturally, I wanted a better solar calculator to use with Clojure. So I wrote one in pure clojure. And it seemed silly to write just a solar calculator library when I had a whole astronomy book at my disposal, so I wrote a library for astronomy and implemented the Sun as the first celestial body of the library. I hope to get more bodies (starting with the Moon) into the library eventually, but here it is. To use Astro-Algo, add the following dependency to your Leiningen project: [astro-algo "0.1.2"]

The fine print

Astro-Algo is a Clojure library designed to implement computational methods described in Astronomical Algorithms by Jean Meeus. In this release, the following is complete:

  • date utils with reduction of time scales table from 1890 to 2016
  • celestial Body protocol to get equatorial coordinates and standard altitude
  • implementation of Body protocol for the Sun
  • function to get local coordinates for a body
  • function to get passages (rising, transit, setting) for a body

The solar position computations implemented in this library are considered low accuracy. They account for some effects of nutation with simplifications described in AA Ch. 25. These calculations should be accurate to within 0.01 degrees of geometric position, however larger differences may be observed due to atmospheric refraction, which is not accounted for other than the standard values for the effect of atmospheric refraction at the horizon used in the calculation of passage times. Passages are calculated according to AA Ch. 15 by interpolating the time that the celestial body crosses the local meridian (for transit) or the standard altitude (for rising and setting). Standard altitudes for rising and setting are taken from the celestial body unless the passages for twilight are desired, in which case the standard altitude for the specified twilight definition is used.

The core API

com.climate.astro-algo.bodies

protocol Body

equatorial-coordinates

Calculate apparent right ascension and declination of a body.
Input:
  this (Body instance)
  dt (org.joda.time.DateTime UTC)
Returns a map with keys:
  :right-ascension (radians west of Greenwich)
  :declination (radians north of the celestial equator)

English translation: Interface Body has method equatorial-coordinates, which returns the orientation of the Body relative to the earth at any time.

standard-altitude

Returns the "standard" altitude, i.e. the geometric altitude of the
center of the body at the time of apparent rising or setting
Input:
  this (Body instance)
Returns standard altitude (degrees)

English translation: Interface Body has method standard-altitude, which returns the angle above the horizon at which the body is considered to have risen or set. (This can be different for different kinds of celestial bodies.)

local-coordinates

Calculate local coordinates of a celestial body
Input:
  body (Body instance)
  dt (org.joda.time.DateTime UTC)
  lon - longitude (degrees East of Greenwich)
  lat - latitude (degrees North of the Equator)
Returns a map with keys:
  :azimuth (radians westward from the South)
  :altitude (radians above the horizon)

English translation: Function local-coordinates returns the orientation of a Body relative to an observer on earth at a specified longitude and latitude.

passages

Calculate the UTC times of rising, transit and setting of a celestial body
Input:
  body (Body instance)
  date (org.joda.time.LocalDate)
  lon - longitude (degrees East of Greenwich)
  lat - latitude (degrees North of the Equator)
  opts
    :include-twilight (one of "civil" "astronomical" "nautical" default nil)
    :precision - iterate until corrections are below this (degrees, default 0.1)
Returns a map with keys:
  :rising (org.joda.time.DateTime UTC)
  :transit (org.joda.time.DateTime UTC)
  :setting (org.joda.time.DateTime UTC)

English translation: Function passages returns the times that a Body will rise, transit and set given a date, longitude and latitude. You can optionally specify a twilight definition to include, which changes the angle below the horizon used for rising and setting. Optional precision allows you to tune the amount of interpolation before convergence.

Usage examples

(:require
  [com.climate.astro-algo.bodies :refer [local-coordinates passages]]
  [clj-time.core :as ct])
(:import
  [com.climate.astro_algo.bodies Sun])

; Suppose we want to get the cosine of the solar azimuth
; on 10/18/1989 00:04:00 (UTC) in San Francisco 
; (at the moment of the 1989 Loma Prieta earthquake).

; We begin by start by instantiating a Body (the Sun).
(let [sun (Sun.)]
  ; Now we declare the time/space coordinates we need
  ; for local-coordinates...
  (let [dt (ct/date-time 1989 10 18 0 4 0)
        lon -122.42 ; [degrees east of Greenwich]
        lat 37.77   ; [degrees north of the Equator]
        ; ...and get azimuth and altitude in radians.
        {:keys [azimuth altitude]} (local-coordinates sun dt lon lat)]
    ; To get the cosine of the zenith,
    ; we can take the sin of the altitude.
    (Math/sin altitude)))

; => 0.2585521627566685
; This represents the ratio of the number of photons hitting a horizontal
; surface to the number of photons hitting a surface facing directly at
; the sun. A negative value would indicate that the sun is geometrically
; beneath the horizon.

(let [sun (Sun.)]
  ; get civil daylength on 12/21/2013 at Stonehenge
  (let [date (ct/local-date 2012 12 21)
        lon -1.8262 ; [degrees east of Greenwich]
        lat 51.1788 ; [degrees north of the Equator]
        {:keys [rising transit setting]}
        (passages sun date lon lat
                  ; can be none, civil, astronomical, nautical
                  ; none is default
                  :include-twilight "civil"
                  ; set lower precision to sacrifice some accuracy
                  ; and speed up the computation
                  ; 0.1 is default (degrees)
                  :precision 1.0)]
    ; calculate daylength in hours
    (-> (ct/interval rising setting)
      ct/in-millis ; [milliseconds]
      (/ 3600000.))))  ; [hours]

; => 9.214338611111112
; If we're going to be there all day,
; we're going to need more than one bag of chips.

How does it stack up?

I did a spot check to compare the accuracy of Astro-Algo against NOAA’s calculator and java library for sunrise and sunset (called luckycat in the chart below).

chart_3

As you can see, we are making progress!

How can I learn more?

The complete API is documented at the astro-algo GitHub project.

Can I get involved?

Yes! We love pull requests. If you want to get involved, getting a copy of AA is highly recommended.

Thanks!

Special thanks to Bhaskar (Buro) Mookerji, who reviewed my code and gave great feedback. Big thanks also to Buro, Leon Barrett, and Steve Kim for reviewing this post.

On deck

This library could use a refactor for more efficient batch jobs (like calculating solar local coordinates for a bunch of locations at the same time). I’d also like to add the calculations for lunar position and phase soon. Of course, I might have to work on the beer garden picnic table problem first… 🙂

Tagged with: , ,
Posted in Engineering
2 comments on “Astro-Algo: Astronomical Algorithms for Clojure
  1. Super cool! thanks – I’ve been working on implementing the AA’s algorithms myself for my project – http://vijaykiran.com/2013/10/trailer-for-my-clojure-project-yearlight/

  2. Teresa Waters says:

    December 14, 2015

    Dear Mr. Eckerle,
    Thank you for helping us with hour of code. I also liked when you talked about what you do at your job. Now I don’t know If I want to write code, or be a mechanical engineer. Being able to code for a company and it being your job must be really cool. Thanks again.

    Ian, writing for Mrs. Waters’ class

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s

%d bloggers like this: