Introduction

Today there was a solar eclipse that was not visible on my side of the
Atlantic, but was seen on the European side, either as a partial eclipse,
towards the south, or a total one, towards the north [1]. Eclipses being rare
and solar power being a new thing, this event caused unprecedented reduction of
solar power [2].

A good spot for viewing the total eclipse was Longyearbyen, on Svalbard, and
readers seeking good images might want to include that name in web searches. (I
could not find open-source images at the time of writing, but of course that
was not really my goal in this blog posting…)

I thought it might be interesting to see whether the sun and moon functions in
the oce package could reproduce the eclipse timing, so I constructed a function
to measure the mismatch between sun and moon position in the sky, and set up an
optimization problem to find the time of least mismatch.

The oce functions sunAngle() and moonAngle() are at the heart of the
work. Each returns a list that contains, among other things, altitude and
azimuth. I set up a mismatch function to calculate a combination of these,
with a scale factor (line 11 of the code) to account for a sort of
converging-meridian effect of azimuth. This function is unity at the horizon,
where a degree change in azimuth is the same angular distance as a degree
change of altitude, and it falls to 0 at the zenith.

Below is the code, and the graph it makes. The black line is the estimated time
when sun and moon centres are nearest each other in the sky, and the time is
written in black near the top of the graph. The red line is the estimate of
eclipse maximum from e.g. [3] The red line is the estimate of eclipse maximum
from [3].

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27

library(oce)
angle <- function(t, lon=15+40/60, lat=78+12/60)
{
    ## Default location is Longyearbyen, Svalbard
    sa <- sunAngle(t, lon, lat)
    ma <- moonAngle(t, lon, lat)
    saz <- sa$azimuth
    sal <- sa$altitude
    maz <- ma$azimuth
    mal <- ma$altitude
    C <- cos(0.5*(ma$altitude+sa$altitude) * pi / 180)
    sqrt((C*(saz-maz))^2 + (sal-mal)^2)
}
# time is from reference 3
nasa <- as.POSIXct("2015-03-20 9:45:39", tz="UTC")
times <- nasa + seq(-1800, 1800, 30)
misfit <- NULL
for (i in seq_along(times)) {
    misfit <- c(misfit, angle(times[i]))
}
oce.plot.ts(times, misfit, ylab="Angular misfit [deg]")
o <- optimize(function(t) angle(nasa+t), lower=-3600, upper=3600)
eclipse <- nasa + o$minimum
abline(v=eclipse, col='black')
abline(v=nasa, col='red')
mtext(sprintf("Here: %s ", format(eclipse)), line=-1, adj=1, col="black")
mtext(sprintf("NASA: %s ", format(nasa)), line=-2, adj=1, col="red")

Results and discussion

The present method suggests the maximum eclipse occured 3 minutes after the NASA estimate (see [3]).

I’m not too sure why this is, so I’ve put some thoughts in the exercises.

Exercises

1. Verify or correct the “C” factor in line 11 of the code.

2. Determine why the times disagree by a few minutes.

3. Determine whether oce code is now out of synch with UTC. I notice on the
   NASA site [3] they show TD (terrestrial dynamical time) being about a minute
   after UTC [4]. I wonder whether TD may be a replacement for the time I used
   in the oce code, which was based on algorithms from the 1970s?

References and resources

1. Overview of eclipse (wikipedia).

2. Effect of eclipse on power grids (reuters).

3. Image from
   NASA
   showing eclipse details, including timing.

4. Time types at NASA

5. R source code used here: 2015-03-20-eclipse.R.

6. Jekyll source code for this blog entry: 2015-03-20-eclipse.Rmd