Finding lat/lon during an eclipse

[Bernmeister]

Using this example as a starting point, I wanted to find the latitude/longitude to match those in eclipse tables for a given eclipse. All good for lunar eclipses, not so for solar and I'm not sure where I'm going wrong...

#!/usr/bin/env python3
# -*- coding: utf-8 -*-


from skyfield.api import load, wgs84
from skyfield.positionlib import position_of_radec
import datetime


ephemeris = load( "de421.bsp" )
timeScale = load.timescale( builtin = True )

# https://eclipse.gsfc.nasa.gov/5MCLE/5MKLEcatalog.txt
lunarEclipses = [
    "09701   2022 Nov 08  11:00:22     73    282  136   T+  p-   0.2570  2.4143  1.3589  353.9  219.8   85.0   17N  169W",
    "09702   2023 May 05  17:24:05     73    288  141   N   h-  -1.0349  0.9636 -0.0457  257.5    -      -     17S   98E",
    "09703   2023 Oct 28  20:15:18     74    294  146   P   a-   0.9471  1.1181  0.1220  264.6   77.4    -     14N   52E",
    "09704   2024 Mar 25  07:13:59     74    299  113   N   -t   1.0609  0.9557 -0.1325  279.1    -      -      1S  106W",
    "09705   2024 Sep 18  02:45:25     74    305  118   P   -a  -0.9792  1.0372  0.0848  246.3   62.8    -      3S   42W" ]

# https://eclipse.gsfc.nasa.gov/5MCSE/5MKSEcatalog.txt
solarEclipses = [
    "9558  478   2022 Oct 25  11:01:20     73    282  124   P   -t   1.0701  0.8619  61.6N  77.4E   0  244",             
    "9559  478   2023 Apr 20  04:17:56     73    288  129   H   -n  -0.3952  1.0132   9.6S 125.8E  67  334   49  01m16s",
    "9560  478   2023 Oct 14  18:00:41     74    294  134   A   -p   0.3753  0.9520  11.4N  83.1W  68  208  187  05m17s",
    "9561  479   2024 Apr 08  18:18:29     74    300  139   T   n-   0.3431  1.0566  25.3N 104.1W  70  149  198  04m28s",
    "9562  479   2024 Oct 02  18:46:13     74    306  144   A   p-  -0.3509  0.9326  22.0S 114.5W  69   31  266  07m25s" ]


def printLatitudeLongitude( body, dateTimeStringFromCatalog ):
    dateTimeUTC = datetime.datetime.strptime( dateTimeStringFromCatalog, "%Y %b %d %H:%M:%S" ).replace( tzinfo = datetime.timezone.utc )
    now = timeScale.from_datetime( dateTimeUTC )
    ra, dec, distance = ephemeris[ "earth" ].at( now ).observe( body ).radec()
    subpoint = wgs84.subpoint( position_of_radec( ra.hours, dec.degrees, t = now, center = 399 ) ) # NAIF code for the Earth center of mass
    print( subpoint.latitude, subpoint.longitude )


print( "Lunar eclipses: catalog versus moon as body:" )
for eclipse in lunarEclipses:
    dateTime = ' '.join( eclipse.split()[ 1 : 4 + 1 ] )
    latLon = ' '.join( eclipse.split()[ 16 : 17 + 1 ] )
    print( latLon )
    printLatitudeLongitude( ephemeris[ "moon" ], dateTime )
    print()

print()

print( "Solar eclipses: catalog versus moon as body versus sun as body:" )
for eclipse in solarEclipses:
    dateTime = ' '.join( eclipse.split()[ 2 : 5 + 1 ] )
    latLon = ' '.join( eclipse.split()[ 13 : 14 + 1 ] )
    print( latLon )
    printLatitudeLongitude( ephemeris[ "moon" ], dateTime )
    printLatitudeLongitude( ephemeris[ "sun" ], dateTime )
    print()

Results are:

Lunar eclipses: catalog versus moon as body:
17N 169W
16deg 51' 14.6" -169deg 15' 06.6"

17S 98E
-17deg 14' 40.7" 97deg 45' 30.8"

14N 52E
14deg 05' 12.2" 51deg 42' 46.4"

1S 106W
-01deg 12' 13.5" -106deg 33' 51.3"

3S 42W
-02deg 35' 13.9" -42deg 20' 37.3"


Solar eclipses: catalog versus moon as body versus sun as body:
61.6N 77.4E
-11deg 14' 40.1" 11deg 09' 40.4"
-12deg 10' 25.1" 10deg 41' 30.5"

9.6S 125.8E
11deg 04' 41.3" 115deg 27' 41.1"
11deg 25' 02.2" 115deg 16' 45.5"

11.4N 83.1W
-07deg 56' 43.1" -93deg 29' 46.7"
-08deg 14' 45.5" -93deg 39' 59.6"

25.3N 104.1W
07deg 54' 23.2" -94deg 21' 45.4"
07deg 35' 38.2" -94deg 12' 30.1"

22.0S 114.5W
-04deg 15' 59.7" -104deg 25' 27.7"
-03deg 59' 13.1" -104deg 16' 46.7"

The first lot of results are for lunar eclipse and the latitude/longitude matches nicely. Not so for the second lot of results for solar eclipses...any ideas please?

[brandon-rhodes]

It looks like you are computing the Earth locations for which the Sun is overhead—which I think would only be the location on Earth of the Moon's shadow, and thus of the eclipse totality, if the centers of the three bodies Earth-Moon-Sun were in an exact straight line?

But that's never the case for a solar eclipse. To a smaller or greater degree, the Moon will happen to track above or below the Earth-Sun line, and the sweep of its shadow will pass north or south of the Earth location from which the Sun is overhead.

What Skyfield will need, when someone has time to implement it, is a calculation of where the Moon's shadow falls. I've glanced over the guide to the calculation in the Explanatory Supplement but haven't had time to try turning it into code. A general routine for finding the intersection of a line and an ellipsoid may be needed? I thought I'd written one up to solve another problem, but all I see at the moment is intersect_line_and_sphere() in geometry.py so maybe I never got the more general case written up; there's a SPICE routine that could be helpfully adapted, if I recall.

[Bernmeister]

Thanks @brandon-rhodes. I suspect that after reading:

It looks like you are computing the Earth locations for which the Sun is overhead—which I think would only be the location on Earth of the Moon's shadow, and thus of the eclipse totality, if the centers of the three bodies Earth-Moon-Sun were in an exact straight line?

But that's never the case for a solar eclipse. To a smaller or greater degree, the Moon will happen to track above or below the Earth-Sun line, and the sweep of its shadow will pass north or south of the Earth location from which the Sun is overhead.

that I had/have an incorrect picture in my head of the solar eclipse; perhaps lunar too!

For lunar, I assumed that given the date/time of the eclipse (from tables), that if I compute the RA/Dec for the moon and pass to wgs84.subpoint( position_of_radec() ), that would give me a vector (and subsequent lat/long). That matches the value in the tables and so I did the same for solar (again computing the moon's RA/Dec as it is the moon we're "looking" at).

When the lat/long didn't match that in the tables (for solar), I switched to using the sun's Ra/Dec and again did not match the tables, but the Ra/Dec matched closely with that of the moon.

So if I get the correct lat/long for a lunar eclipse, but incorrect for solar, did I just get lucky for the lunar case and completely misunderstood solar? I did't assume a straight line from Earth-Moon-Sun for any of the cases; simply I assumed (naively or otherwise) that at a date/time of any choice, if I compute the Ra/Dec of a body (moon or sun) then what is the equivalent lat/long.

I would not have thought that we'd need to project the shadow back on to Earth...very confused!

[brandon-rhodes]

I would not have thought that we'd need to project the shadow back on to Earth...very confused!

I have always found the illustrations of Guy Ottewell useful in understanding the scales and geometry of astronomy, so this might be a useful place to start if you want to learn more:

https://www.amazon.com/Under-Standing-Eclipses-Guy-Ottewell/dp/0934546746/

And the cover illustration itself might be helpful:

One point of confusion might be the tables you have found, and are comparing your results against. They are not, alas, telling you where the lunar eclipse is. The lunar eclipse, strictly speaking, is something happening over on the Moon, and its location is determined the same way that the location of a solar eclipse is determined: by projecting the shadow of the eclipsing body to learn where it intersects the spheriod of the being-eclipsed body.

Your table, instead of telling you where the lunar eclipse is happening—"here's the lunar lat/lon where Earth's shadow is deepest"—is merely telling you which Earth observers are positioned to look straight up and see the eclipsed or partially-eclipsed Moon directly overhead. That's a simple subpoint calculation, as the opposite question would be: "where (lat/lon) can I stand on the Moon and see the Earth directly overhead as it's being eclipsed?"