Constrain the search for partial eclipse semiduration to
within what we already found for the penumbral eclipse.
Same for total/partial. This is a very small improvement because
narrowing the search window does not improve quadratic interpolation
very much. But it is an extremely cheap and safe optimization.
It turns out that searching plus or minus 0.03 days around the
full moon is ample for finding minimum shadow distance.
This reduces CalcMoon() call count from 127155 to 112827.
Performance ratio with original algorithm = 5.13.
When the full moon's ecliptic latitude is larger than 1.8 degrees,
even a penumbral eclipse is not possible. Thus there is no need
to search for the minimum shadow distance in that case.
This decreased unit test CalcMoon() count to 127155.
Improvement ratio over original algorithm = 4.55.
Greatly reduced the number of CalcMoon() calls needed to find
the time of the minimum shadow distance, when searching for a lunar eclipse.
Use Astronomy_Search() instead of dumb search.
Added undocumented global variable for counting how manyh times CalcMoon()
is called.
The call count went from 578569 down to 207186 (ratio = 2.79).
Execution time likewise decreased from 2.9 seconds to 1.1.
We don't need to generate Delta T tables any more because
none of the Astronomy Engine supported languages use it any more.
Now they all use Espenak/Meeus polynomials instead.
I'm in the process of replacing how Astronomy Engine calculates
Delta T. Instead of a series of line segments based on canned data,
I'm switching over to use the Espenak/Meeus piecewise polynomials.
Also allowing the user to change the Delta T function to match
an external reference. I will use this in the unit tests that
reference JPL Horizons data, so that I can greatly tighten the
test tolerances.
I downloaded the airless_Moon_dt.txt file from JPL Horizons so that
I could reverse engineer their DeltaT function. But keeping it in
the horizons directory caused unit test breakage. Moved to delta_t
directory so dtplot can still use it to make graphs.
I had to increase certain error tolerances in the unit tests.
Reworked the unit tests to make more sense by waiting until
each language step is done to check against each other.
That way I can run a single language step independently.
This bug has been there a long time, but I didn't notice it
until I started experimenting with larger Delta-T values.
Fortunately, the bug was in the unit test, not Astronomy Engine
release code. Calling NOVAS place(), I need to pass in Delta-T
expressed in seconds, not days.
Write eclipse/c_le_stats.csv that contains all the error
values for each lunar eclipse calculation. This helped
me figure out the main source of lunar eclipse error was
discrepancy with Delta-T values in the future.
Using some trial and error, I found that using 85 km instead of 65.4 km
for the thickness of the Earth's atmosphere results in better overall
fit with the test data.
Adding support for lunar eclipse calculations in the C code
caused me to tweak the values of the Sun and Moon radii.
This in turn caused unit test failures.
Made slight changes to the unit tests to get things passing again.
I found that lunar eclipse data is available for many centuries.
I downloaded the data for the years 1701..2200.
Wrote norm.py to extract and convert the parts I care about
into a format that will be much easier to parse in the unit
tests for all four languages.
Regenerate the normalized data from the 'run' script.
This way, I have documentation for where the data came from.
Using geocentric Moon instead of heliocentric Moon
gives more floating point precision for determining
the distance between the Moon and the Earth's shadow ray.
I figured out a formula that determines how far away
the Moon is from the center of the Earth's shadow.
This confirms the formula makes sense for a known
total lunar eclipse on May 26, 2021.
Increase type safety by making the enumerated type Body
derive from Enum rather than IntEnum, as recommended by
https://www.python.org/dev/peps/pep-0435/
Fixed places where I was treating Body values as integers.
Now when a Time object is evaluated and represented in
the Python interpreter, it results in a string of the form:
astronomy.Time(ut)
where ut is the numeric representation of the ut field.
This mimics the exact way such a Time value could be constructed.
That is, eval(repr(t)) results in a time value equal to t.
It turns out I was off by nearly 18 hours in the B1875 epoch.
This has a tiny effect on the orientation of the Earth's axis.
Instead of: ut = 1875-01-01T12:00:00.000Z
the correct epoch is: ut = 1874-12-31T18:12.21.950Z
See the comments in the Constellation functions in
each of the source files for more info.
