My GitHub actions automated tests kept failing today
because commit_hook.bat kept failing to download
the Doxygen zipped binaries for Windows.
So I have downloaded it myself and mirrored it in
my `ephemeris` repo where I keep other large files that
are only needed by Astronomy Engine tests, not end users.
Downloading from GitHub to a GitHub Actions worker
should be much more reliable (probably faster too).
Another advantage is the version will be stable,
so I don't have to keep fixing things every time they
make a change to Doxygen. It's good enough the way it is!
Instead of stopping at the first failure, diffcalc
scripts (Linux and Windows) will now continue to run
all tests before reporting pass/fail. This allows
fixing multiple issues all in one GitHub Actions pass.
Also fixed the immediate issue where I needed to
increase tolerances slightly, as a follow-up fix
for nutation formula changes.
While trying to convert ecliptic coordinates from mean
equinox of date to true equinox of date, I ran into excessive
overhead from the IAU2000B nutation model. The fact that it
uses 77 trigonometric terms made the calculations a lot slower.
https://apps.dtic.mil/sti/pdfs/AD1112517.pdf
Page 4 in the above document mentions a shorter series
“NOD version 2” that has 13 terms instead of 77 as used in IAU2000B.
I had not noticed NOD2 before, because it appears only in
the FORTRAN version of NOVAS 3.x, not the C version.
After reading the FORTRAN code, I realized NOD2 is the same
as IAU2000B, only it keeps the first 13 of 77 terms.
The terms are already arranged in descending order of
significance, so it is easy to truncate the series.
Based on this discovery, I realized I could achieve all of
the required accuracy needed for Astronomy Engine by
keeping only the first 5 terms of the nutation series.
This tremendously speeds up nutation calculations while
sacrificing only a couple of arcseconds of accuracy.
It also makes the minified JavaScript code smaller:
Before: 119500 bytes.
After: 116653 bytes.
So that's what I did here. Most of the work was updating
unit tests for accepting slightly different calculation
results.
The nutation formula change did trigger detection of a
lurking bug in the inverse_terra functions, which convert
a geocentric vector into latitude, longitude, and elevation
(i.e. an Observer object). The Newton's Method loop in
this function was not always converging, resulting in
an infinite loop. I fixed that by increasing the
convergence threshold and throwing an exception
if the loop iterates more than 10 times.
I also fixed a couple of bugs in the `demotest` scripts.
Modified the code generator and source templates to allow
using fewer than 77 terms from the IAU2000B nutation model.
The number of terms causes calculations to be far slower than
I would like, plus most of these terms provide far less than
one arcsecond of difference in the output.
I want to experiment with truncated versions of the series.
https://apps.dtic.mil/sti/pdfs/AD1112517.pdf
Page 4 in the above document references a shorter series
"NOD version 2" that has only 13 terms instead of 77 as used in IAU2000B.
I found that code and confirmed the 13 terms are the first 13
consecutive terms from the original 77.
This commit does not change the number of terms calculated;
it only enables doing so by changing a single #define in
generate/codegen.c.
Once I change the nutation model, I'm sure there will be multiple
tweaks to unit tests to get everything working again.
Astronomy_SphereFromVector was not checking its vector
argument for having a bad status. Now it does.
Added comments that clarify exactly what nutation and precession
functions do.
Added a performance test to measure how fast
a call to EclipticGeoMoon is. Results:
Over 3 trials, mean time = 24.064 seconds.
Number of calls = 14,609,700.
Time per call = 1.647 microseconds.
That is faster than I thought!
Also replaced conditional compilation with a runtime flag
to indicate whether a test should be skipped by the "all"
command. This allows me to run performance tests without
hacking a #define in the code, plus it ensures that even
though I don't run the performance tests every time I change
Astronomy Engine, they at least compile without errors.
This was another test where I was trying to figure out
an apparent anomaly that was actually based on my
flawed assumption about why obliquity changes over time.
I don't need the MoonLatitudes test any more.
Its purpose was to investigate what I thought was an
ecliptic latitude discrepancy. Now I understand that
changes in obliquity are caused by changes in the ecliptic
plane, not the equatorial plane.
I'm still trying to understand a discrepancy
between ecliptic latitudes calculated in two
different equinoxes: mean equinox of date and
mean J2000 equinox. I keep thinking the two
latitude values should be the same, because they
are measured against an essentially unchanging
ecliptic plane. I understand that ecliptic longitudes
should change as the Earth's equator precesses
by 47 arcseconds per century.
Added new test MoonEcliptic(), that compares
JPL Horizons calculations of the Moon's ecliptic
coordinates relative to equator of date versus
Astronomy_EclipticGeoMoon():
C MoonEcliptic: PASS: count=73050, max lat=1.825 arcsec, max lon=23.462 arcsec.
It is so weird that the maximum latitude discrepancy is less than 2 arcseconds.
Updated MoonVector test to measure maximum errors in addition to RMS errors.
Display distance errors in kilometers instead of dimensionless units.
Defined constant JD_2000 so all my tests that
convert Julian Dates from JPL Horizons can share it.
There is a discrepancy between ecliptic latitudes measured
using two different mean equinoxes: 1900 and 2000.
The discrepancy is very close to the change of obliquity
expected for a 100-year period: 47 arcseconds.
But ecliptic latitudes should remain constant as
the Earth's axis precesses!
MoonLatitudes now generates a CSV file for me to study the
ecliptic latitude discrepancy as a function of time.
I can see that it increases away from the year 2000,
and that it is periodic with the Moon's orbit.
This indicates an error in correcting for precession.
Generated new JPL Horizons test data for MoonVector that
covers the years 1900..2100 instead of 1980..2020.
Confirmed that its maximum error is still 1.16 arcsec.
That is error comparing EQJ vectors and ECL vectors.
The C unit test program ctest.c now includes a MoonLatitudes
function invoked by the name `moon_lat`.
This function demonstrates that there is a 47 arcsecond
discrepancy between ecliptic latitudes using the mean
equator of date and ecliptic latitudes using the mean
equator of J2000. This should not happen! I'm investigating.
Also, when running ctest without command-line arguments,
print a list of all available test names, so I don't have
to go open up the file ctest.c and read the code.
The documentation for SearchRiseSet and SearchAltitude needed
clarification about refraction and the part of the body solved
for (center versus limb). The JavaScript version was especially
lacking compared to documentation for the other languages.
Also documented SearchAltitude's limitations; it does not
work at or near maximum/minimum altitude.
Mention that user-defined stars are allowed for
SearchRiseSet, SearchAltitude, and SearchHourAngle.
Fixed a couple places where the Kotlin documentation had
broken links to other functions.
I had a copy-n-paste typo in the `dec` parameters
for all of the DefineStar functions. Fixed it.
The TypeScript version of HelioState did not handle
user-defined stars. Added support there.
Because we instantly know the heliocentric
distance of a user-defined star, there is no
need to convert it into a vector and then take
the length of the vector.
All of the HelioDistance functions now return
the distance directly, as an optimization.
Also, I decided it didn't make sense to have a
default definition for user-defined stars.
If the caller doesn't define a star, it should
be treated as an invalid body.
Added Python support for user-defined stars.
Defined new StateVector methods: Position and Velocity.
Defined division operator: Vector / float.
Bumped version number to 2.1.12.
SearchRiseSet and SearchHourAngle now work with user-defined stars.
Fixed a bug in InternalSearchAltitude where I failed to return an error
code when its call to MaxAltitudeSlope failed.
DefineStar now requires passing in the heliocentric
distance of the star expressed in light-years.
That way, I can directly support returning vectors
to a star from HelioVector, GeoVector, etc.
SearchRiseSet and SearchHourAngle now work with user-defined stars.
I'm starting to implement the ability to define
up to 8 distinct points in the sky as "stars"
that will be allowed as a `body` parameter to
some Astronomy Engine functions, to be determined.
Made sure all the altitude search functions
verify that the geographic latitude and target altitude
are valid numbers in the range [-90, +90].
Reworked the C version of the code to be clearer:
eliminated goofy ALTDIFF macro, split out max
altitude derivative into its own function MaxAltitudeSlope,
just like the other language implementations do.
Minor rewording of comments in MaxAltitudeSlope functions.
Python InvalidBodyError now includes the invalid body
in the diagnostic message.
Changes to the altitude search caused a tiny
change in rise/set times in all the demo programs.
I remembered to update the expected output for
Kotlin, but forgot to do it for Java.
The C# version of the altitude searches now work correctly
near the Earth's poles. This is a port of the C version.
Cleaned up a little code in the C version.
Do not recurse deeper when any solution would
have a rise+set or a set+rise consecutively
within one second. This would risk confusing the
search, which has a 0.1 second tolerance for convergence.
For practical reasons, such events are dubious anyway,
being swamped by unpredictable atmospheric refraction.
Plus it's just faster to not have to recurse so deeply.
Add a safety valve that fails the search if recursion
gets deeper than the above tolerance would ever allow
to happen in the first place.
This change improves performance significantly:
C RiseSet: passed 5909 lines: time errors in minutes: rms=0.2904, max=1.1541, recur=15, altcount=156881
Added helper function AscentError for cleaner code in error cases.
