Pipeline design

The configuration record

Every run starts by fixing, and every output ends by printing, the following. These fields specify the inputs needed for a reproducible comparison. Algorithm and solver choices also matter ^{7 8}.

Sources for each choice are in the stage notes below and in Earth model and time and lunar and solar model. The edge-mode radius is an explicit modelling option from the cited flash-spectrum analysis. It is not a settled physical constant. Keep the 959.63 arcsecond convention available for reproducing the corresponding tables. Detector and limb definitions affect other radius determinations ^{10 11 12}. The frame kernel is named moon_de440_200625.tf in Wright and Young's appendix; NAIF now serves the 2025 revision of the same frame as moon_de440_250416.tf ^{6 13}.

Stage 1: ephemeris, orientation and time

Inputs. de440s.bsp (1849 to 2150, 31 MB) for the Sun, Earth and Moon; the DE440 lunar PCK and frame kernel; IERS finals.all or USNO deltat.data and deltat.preds; a leap-second table. Reject dates outside the chosen kernel's coverage. Catalogue work over wider epochs needs a full kernel covering those dates ^{14 15 16 17}.

Contract. The layer exposes, at a TT instant: apparent geocentric places of the Sun and Moon on the true equator and equinox of date, with light-time and aberration applied to both; apparent topocentric places for a site; Greenwich apparent sidereal time; the Moon-to-observer vector in the Moon mean-Earth frame for libration. The Moon is read as body 301 relative to 399 and the Sun as body 10, never the barycentre, since the Sun wanders up to 1.6 solar radii from it ¹⁸. TT may stand in for TDB, the difference being under 2 m on the shadow ¹⁹.

The trap. Mixing a geometric Sun with an apparent Moon displaces the shadow axis by 20″, which is 38 km on the ground. Using TT where UT belongs displaces it by 32 km. These two errors dwarf everything else in the budget and are the first things to test ^{20 21}.

ΔT. Three regimes: measured IERS or USNO values for 1962 onward, deltat.preds with its error column for a few years ahead, and the Stephenson, Morrison and Hohenkerk 2016 spline and parabola beyond that, with the tidal-acceleration adjustment matched to the ephemeris ^{22 23}. ΔT is stored as metadata and applied once, when μ is converted to a longitude in stage 4 or to a local hour angle in stage 6, at 15.041″ of rotation per second ^{3 24}.

Stage 2: enumeration and provisional type

Method. Step through mean new moons. Reject a lunation when $|\sin F| > 0.36$. Compute gammagammaThe distance of the shadow axis from the Earth's centre in Earth equatorial radii at greatest eclipse, positive when the axis passes north of the centre. Values of $||$ below about 0.9972 give a central eclipse. and $u$ from Meeus's short series. Classify: central if $|\gamma| < 0.9972$; no eclipse if $|\gamma| > 1.5433 + u$; non-central umbral if $0.9972 < |\gamma| < 0.9972 + |u|$; total if $u < 0$; annular if $u > 0.0047$; for $0 \le u \le 0.0047$, hybrid if $u < 0.00464\sqrt{1 - \gamma^2}$ and otherwise annular; partial magnitude $(1.5433 + u - |\gamma|)/(0.5461 + 2u)$ ²⁵.

Then refine. For each candidate, compute the elements of stage 3 and take γ as $\sqrt{x^2+y^2}$ at its minimum, with the sign of $y$ to preserve the catalogue's north/south convention ²⁶. Decide the final type on the umbral radius at the surface point under the axis, at greatest eclipse and at both ends of the central line, so that hybrids and their three classes fall out. Keep a non-central class: 94 eclipses in five millennia have a one-limit track and no central line ^{27 28}.

Saros. Kluepfel's closed form from the lunation number reproduces the Saros number of all 11,898 catalogued eclipses ²⁹.

Validation. The 221 eclipses of 1951 to 2050 in NASA's ASCII catalogue: zero type disagreements, greatest-eclipse times within 1.2 min from the series, then within 0.1 s from the elements ³⁰.

Stage 3: Besselian elements

Inputs. Apparent geocentric places from stage 1, $k_1$, $k_2$, $s_0$, $\pi_0 = 8.794″$, $a_e$, ΔT.

Equations, in order ^{3 31}. Form the dimensionless distance ratio $b$ before choosing separate output units. In the Supplement's parallax form, $b=\sin\pi_0/(R\sin\pi_m)$ with the Sun's distance $R$ in au. For the coordinates below, $r_m=1/\sin\pi_m$ is the lunar distance in Earth equatorial radii. The cone equations retain $R$ in au. Convert all angular inputs, including $s_0$ and $\pi_0$ supplied in arcseconds, to radians before evaluating trigonometric functions. The shadow axis direction $(a,d)$ and scale $g$ come from $g\cos d\cos a = \cos\delta_\odot\cos\alpha_\odot - b\cos\delta_m\cos\alpha_m$ and its two companions. The Moon's coordinates in the fundamental system are

The hour angle is $\mu = \theta - a$ with θ the apparent sidereal time evaluated at the TT instant as if it were UT, the ephemeris sidereal time of the 1961 Supplement. NASA's tables, its JavaScript calculator, Eclipse-Engine and the 2027 visualiser all tabulate μ this way, so a table's μ is free of ΔT and the shift is applied downstream ^{3 5 24}. The cone half-angles are

the vertex heights $c_1 = z + k_1/\sin f_1$ and $c_2 = z - k_2/\sin f_2$, and the shadow radii on the plane $l_1 = c_1\tan f_1$ and $l_2 = c_2\tan f_2$. $l_2 < 0$ is total, $l_2 > 0$ annular. Publish $l_2$ with its sign.

Output. A function elements_at(t_tt) returning $(x, y, d, \mu, l_1, l_2, \tan f_1, \tan f_2, z)$, and for export the NASA polynomial form: cubic in $x$ and $y$, quadratic in $d$, $l_1$, $l_2$, linear in μ, fitted to five samples over six hours about $t_0$ ³². Metadata: ephemeris, $t_0$, ΔT, $k_1$, $k_2$, $s_0$, ellipsoid.

Validation. The 1961 Supplement's worked example for 1961 February 15 08h ET, and the NASA 2024 April 8 polynomial at $t = 0$, expecting agreement at $10^{-5}$ Earth radii once ephemeris, $k$ and ΔT match. A difference above $10^{-4}$ Earth radii in $x$ or $y$, about 640 m, means a frame or aberration inconsistency ^{3 33}.

Stage 4: global circumstances, smooth Moon

Frame. Bessel's substitution turns the ellipsoid into a unit sphere with $\rho_1 = \sqrt{1 - e^2\cos^2 d}$, $\rho_2 = \sqrt{1 - e^2\sin^2 d}$ and a rotated declination $d_1$. In that frame the central line at time $t$ is $\xi = x$, $\eta_1 = y/\rho_1$, $\zeta_1 = \sqrt{1 - x^2 - \eta_1^2}$ ³⁴. For west-positive longitude, $\lambda_W = \mu_{UT} - H$, where $H$ is the local hour angle recovered from the observer geometry. With angles in degrees and ΔT in seconds, $\mu_{UT} = \mu - \frac{1.002738\times15}{3600}\,\Delta T$. This correction is for ephemeris-meridian elements. Omit it if μ was already computed from actual UT1 sidereal time ²⁴.

Limits. A point is on a limit when the eclipse begins and ends at the same instant, which gives $\tan Q = (b' - \zeta d' - a'\sec Q)/(c' - \zeta\mu'\cos d)$ for the position angle $Q$ of the shadow edge, iterated on ζ. Use the 1992 procedure for whole curves: scan $Q$ by degrees, inverse-interpolate, iterate ζ to $10^{-5}$ Earth radii, assign north or south by the sign of $L\cos Q$ ^{34 35}.

Duration and width. Central duration $D = 2|L_2|/n$ with $L_2 = l_2 - \zeta\tan f_2$ and positive relative speed $n$ in Earth radii per hour, so $D$ is in hours. Accurate elapsed duration is C3 minus C2 for the fixed site. Path width uses the positive magnitude of Mikhailov's expression as transcribed in Stellarium ^{35 36}. Rise and set curves from $\cos(\gamma - M) = (m^2 + 1 - l_1^2)/(2m)$ at $\zeta = 0$; curves of maximum eclipse from $\tan Q = -(y' - \eta')/(x' - \xi')$ sweeping ζ ³⁴.

Output. Central line, limits, outlines at stated instants, rise and set curves, maximum-eclipse curves, greatest eclipse, greatest duration, ground speed, as GeoJSON with the time system and ΔT in the properties.

Validation. The 1961 examples 9.6 and 9.7 ($\phi = +44^{\circ} 18'.3$, $\lambda = -29^{\circ} 20'.9$, duration 158.6 s at 08h ET). Then NASA's 2024 path table at three times, expecting 1 km and 0.1 s ^{34 37}.

Stage 5: partial-eclipse maps

Equal-magnitude and equal-obscuration curves are not solved directly. The Supplements interpolate on the maximum-eclipse curves; SVS and the open Eclipse-Engine rasterise instead. For each grid cell, optimise the quantity being contoured over its visible interval. Maximum magnitude and maximum obscuration need not occur at the same instant. Use closest approach as a search seed, and keep apparent radii time-dependent. Evaluate obscuration with the exact lens formula and contour with marching squares. Pad the grid edge and refine vertices near the selected product's maximum ^{34 24 38}.

Stage 6: local circumstances, smooth Moon

Observer. Define $\nu = (\cos^2\phi + (1-f)^2\sin^2\phi)^{-1/2}$ and $S = (1-f)^2\nu$. Then $\rho\sin\phi' = (S+h/a)\sin\phi$ and $\rho\cos\phi' = (\nu+h/a)\cos\phi$, where $h$ and $a$ use the same length unit. For west-positive longitude, the local hour angle is $\theta = \mu(t) - \lambda_W - \frac{1.002738\times15}{3600}\,\Delta T$ with angles in degrees and ΔT in seconds. Apply this correction only to the ephemeris-meridian μ, not to an already-UT1 μ ^{20 39}.

Reduction. Form $(\xi, \eta, \zeta)$ and their rates, $u = x - \xi$, $v = y - \eta$, $m = \sqrt{u^2 + v^2}$, and the observer's-plane radii $L_1' = l_1 - \zeta\tan f_1$, $L_2' = l_2 - \zeta\tan f_2$, starting from the tabulated fundamental-plane radii $l_1$ and $l_2$. Closest approach is a search seed obtained by iterating $t \leftarrow t-(uu'+vv')/n^2$. For maximum magnitude or obscuration, optimise that actual product with time-dependent apparent radii. Contacts are the roots of $m = L_1'$ (C1, C4) and $m = |L_2'|$ (C2, C3), found by the 1992 edition's direct inverse interpolation on $u^2 + v^2 - L^2$ rather than the auxiliary angle ^{39 40}.

Magnitude and obscuration. Magnitude $(L_1' - m)/(L_1' + L_2')$ in the penumbra and $(L_1' - L_2')/(L_1' + L_2')$ inside the umbra or antumbra, switching at $m = |L_2'|$. Obscuration is the lens area $S' = (s^2 A + B - s\sin C)/\pi$ with $s = (L_1' - L_2')/(L_1' + L_2')$ the ratio of the lunar to the solar radius. In this source formula, use $A$ at the lunar centre, $B$ at the solar centre, and $C$ at either intersection of the disc boundaries, with $A = \pi-B-C$. These triangle labels are local to the lens formula and are unrelated to Herald's angular offset below. Magnitude and obscuration depend on the radius ratio as well as separation. At totality the reported magnitude switches to the diameter-ratio convention, so it need not be continuous at the inner contact. Position angle $P$ comes from $\mathrm{atan2}(\pm u,\pm v)$ with the sign reversal for interior contacts of a total eclipse. The vertex angle is $V = P-q$, where $q$ is the parallactic angle and $\tan q = \xi/\eta$ ^{39 41}.

Second path. The direct topocentric method: root-find $\delta(t) = r_s \pm r_m$ on apparent topocentric positions from stage 1, as USNO and the Swiss Ephemeris do. Agreement with the Besselian result within 0.1 s on the same inputs is the regression test for the element generator ^{2 42}.

Validation. NASA's program.js output for any city, which is GPL and self-contained, to the second; the Lusaka 2001 case (C2 13:09:19.3 UT, C3 13:12:32.8 UT before limb correction) to 0.5 s ^{5 43}.

Stage 7: terrain

Why. A limit shifts perpendicular to the path by about $h\cot A\sin D$: 577 m per 1000 m of elevation at a 60° Sun, 1.7 km at 30°, 5.7 km at 10°. SVS measured up to 3 km over the 2017 western states. The umbra widens only by $2h\tan f_2$, 82 m at Everest ^{9 44 45}.

Method. Look up orthometric height $H$ in SRTM GL1 (EGM96) or Copernicus GLO-30 (EGM2008) and add the geoid undulation $N$ to get ellipsoidal $h$. The geoid sits between 106 m below and 85 m above WGS84 and is worth up to 170 m of path at a 30° Sun. Intersect the shadow with the terrain surface rather than shifting a sea-level limit by the elevation factor; keep the factor as a check ^{46 47 6}.

Stage 8: lunar limb

Profile construction ⁶. For an observer at distance $d$ from the Moon's centre: convert each DEM pixel to rectangular coordinates on the 1737.4 km datum; rotate by the topocentric libration from stage 1; take polar $(r, \theta)$ of the transformed $(y, z)$; replace $r$ by the angular radius $\alpha = \tan^{-1}(r/(d - x))$; bin θ into 18,000 elements at 0.02°; keep the maximum α per bin. Rebuild when a libration angle moves by 0.01°. Only DEM longitudes between 75° and 105° east and west and a 450 km polar swath need processing. SLDEM2015 at 128 pixels per degree (240 m) pairs with this resolution. LDEM_128 is the one-file development dataset; the 33 MB LDEM_16 is fit only for a smoke test, since a 16 pixel-per-degree cell is about one arcsecond at the limb, 2 to 3 s of contact time ^{48 49}.

Limb test. Symbols: $n$ the number of elements in $L$, $d_0$ the observer-to-Moon distance at which $L$ was built, $d$ the current observer's distance, $c$ the position angle of the Moon's axis, $\delta$ the Sun-to-Moon apparent distance in solar radii, $\phi$ the position angle of the Sun with respect to the Moon, $r$ the Sun's apparent radius in radians. With $s = (d_0/d)/r$, $a = sL_i$, $\theta = 2\pi i/n + c$, $\rho = a^2 + \delta^2 - 2a\delta\cos(\theta - \phi)$: the eclipse is not total if any $\rho < 1$. Trivial exclusions: outside the umbra if $s\max(L) - \delta < 1$, inside if $s\min(L) - \delta \ge 1$. Reverse the sense for annularity. Antialias with $\epsilon = A(p)/10$. The paper also prints the broken-annular expression $\max(L)-\min(L)>\delta$. This is a quoted criterion, not a directly validated comparison with the units defined here. Angular profile radii and solar-radius-normalised separation must first share a common normalisation. Use the full limb-containment tests for classification, and validate any broken-annular criterion against that geometry ⁶.

Contacts and beads. C2 and C3 at a site are the first and last time steps at which the test passes. Baily's beads are the same test per limb element: each event is the instant the solar limb, drawn with Herald's curve $h = 960\,\text{arcseconds}\,(M-1)(1-\cos C)$, clears a valley floor. Here $C$ is the offset from the nominal contact point, which Herald calls $P$. It is not the contact position angle $P$ in stage 6 or the lens-formula triangle angle. Convert the inward-positive totality curve to the profile's outward-positive radial convention before forming a displacement. A valley delays C2 and advances C3 relative to the same reference limb. Report each event's position angle and depth. Keep Herald's tangency construction as a cross-check because the bulletins' examples are stated in its terms ^{50 51}.

Path products. The umbra is emitted as the polygon the raster produces, not as an ellipse. A limit line becomes an interior and an exterior line with a graze zonegraze zoneThe narrow band, typically 5 to 10 km wide, along each umbral limit where the irregular lunar limb makes the eclipse neither wholly total nor wholly partial and Baily's beads persist. NASA bulletins tabulate its interior and exterior boundaries. between them, 5 to 10 km wide ⁵¹.

What the profile replaces. The centre-of-figure correction of the almanacs is retired: LOLA products are referenced to the centre of mass in the DE421 mean-Earth frame, which is the frame the ephemeris orients. The 0.69 km gap between the 1737.4 km datum and the 1738.09 km sphere of $k = 0.2725076$ is a datum conversion: subtract 0.6908 km from heights above the smaller sphere when quoting them against the almanac mean limb ^{48 13}.

Validation. The Lusaka corrections of +4.0 s at C2 and −1.2 s at C3, within 1 s allowing for Watts against LOLA. The SVS 2024 umbra_hi polygon at one second, reproduced from SLDEM2015. The IOTA bead tables for 2017 at Thermopolis and 2023 at Cape Range, each event within 1 s at the radius correction each paper derived ^{43 52 53 12}.

Stage 9: products and uncertainty

Every product carries the configuration record and an estimate that distinguishes known corrections, unresolved model choices and quantified input uncertainty. The error budget tabulates diagnostic scales, not a universal list of standard deviations. Combine compatible standard uncertainties in quadrature only when independence is established. Otherwise propagate covariance or report model scenarios. Test the resulting band against reference observations before assigning a coverage probability ⁵⁴. Limits are evaluated at $s_0$ and its stated alternatives. The cited 100 to 200 m edge scales are case-specific rather than an accuracy guarantee. Products for dates before 1600 or after 2300 carry a longitude gore for the ΔT standard error. This is stricter than the Five Millennium Canon, which draws gores only where σ exceeds 265 s, before the year 1 and after 2300 ^{7 55 56}.

What this does not cover

This design reports geometric contacts. A common one-to-one angular refraction mapping applied consistently to both complete limbs preserves their contact. A centre-only correction with circular limbs does not provide a validated residual. Atmospheric ray bending and near-horizon visibility need a separate specified model. The cited sources do not establish a general contact-time adjustment ^{9 10}. Lunar eclipses, transits and occultations reuse stages 1 and 8 but are out of scope.