Now the C function Astronomy_Ecliptic returns ecliptic
coordinates in true equinox of date instead of the
J2000 mean equinox. I'm doing this because it is a
better fit for physical phenomena that ecliptic
coordinates are often used for. For example, lunar nodes,
eclipses, phase angles, and oppositions make more sense
with true latitude and longitude of date.
I will port these changes to the other languages also.
More work standardizing the nomenclature of the
orientation systems across all language documents.
Added C functions to calculate rotation matrices
for EQJ/ECT and ECT/EQJ.
I'm taking gcc 12.2 for a test drive today.
It reports a few warnings that slipped through earlier versions.
None of the warnings concern me for actual code safety,
but I went ahead and resolved them to keep the build clean.
Also provide a hook for a CPP environment variable to override
the C++ compiler to use, instead of forcing g++.
Define ECT = True Ecliptic of Date in the documentation.
I will soon convert the Ecliptic() functions to return ECT instead of
ECL, but I will retain ECL support via rotation matrix functions.
Somehow gcc didn't warn me that a `void` function
was trying to return a `void` value. It makes sense
from a functional language "unit type" perspective,
but it wasn't intentional, and it is weird. Fixed it.
Added EclipticGeoMoon as output to the temp/*_check.txt files as 'm' lines.
This ensures that all the languages calculate nearly identical values.
Optimized EclipticGeoMoon a little more by eliminating a redundant
call to mean_obliq.
Updated the C function Astronomy_EclipticGeoMoon to
correct ecliptic coordinates for nutation.
This means that the returned value is expressed in
true equinox of date instead of mean equinox of date.
This results in the moon_ecm test decreasing the max
longitude error from 24 arcseconds to 6 arcseconds.
EclipticGeoMoon is now about 40% slower, but it still
runs in about 0.4 microseconds per call.
I bootstrapped based on the pretty good optimizations that
codegen did for the Python version of the (now truncated)
IAU2000B nutation formula. I will do the same for the other
nutation formulas.
I realize now that the URLs I was using to download stuff
from GitHub Actions were redirects. So I need to use 'curl -L'
to follow the redirects. But I also removed the redirect by
using the ultimate URLs.
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.