I tried more distant objects like Jupiter ... Neptune.
This revealed that at increasing distances, the convergence
threshold in inverse_terra needed to increased also.
So now I use 1 AU as a baseline, and scale up linearly
for more distant objects.
Asking the latitude and longitude directly beneath
the Sun causes inverse_terra not to converge, because the
convergence increment `W` never got below 1.48e-8, but the
convergence limit was 1.0e-8. I increased the limit to 2.0e-8
in all programming language versions.
I'm hoping that is a big enough tolerance for all cases now,
but I will do more testing to see if further fixes are required
for even more distant bodies than the Sun.
It was possible for certain starting times to have a search
failure in this demo, because it could place more than
one zero-crossing of the function in the same time interval.
So now we iterate over an interval of 10 days at a time
until we find the solution.
This demo shows how to search for the next time
the Moon reaches extreme ecliptic latitude or
extreme declination. In other words, it finds
when the Moon reaches the farthest north or south,
expressed in either ecliptic coordinates or equatorial
coordinates.
Both angles are measured using the Earth's equator of date.
Replace the abstract class with a parameter of function type.
This allows the documentation to fully explain how to use
`CorrectLightTravel` without having to look at the code.
The documentation was missing a mention of
the `time` parameter in the following TypeScript
functions:
* `Rotation_ECT_EQD`
* `Rotation_EQD_ECT`
Likewise, the `time` parameter was not documented in
the corresponding Kotlin functions:
* `rotationEctEqd`
* `rotationEqdEct`
These mistakes have been corrected.
The Windows build does not run the demo tests like
the Linux one does. This resulted in not copying
the updated astronomy.ts to the calendar demo, thus
missing the copyright year update for 2023.
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.
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.
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.
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.
This is a whole new algorithm that efficiently finds
all rise/set events, even near the poles.
It uses a recursive bisection search that limits
recursion depth by knowing the maximum possible
|da/dt| = change in altitude with respect to time.
The following Python functions now support searching
in forward or reverse chronological order:
SearchRiseSet
SearchAltitude
SearchHourAngle
Made some minor performance improvements to the
other implementations: return sooner if we
go past time window.
The quadratic interpolator used by `Search` was returning
an unused output: `x`. Search does not need this dimensionless
value; it only cares about the time solution `t` and the slope
of the function at `t`. Removed `x` from the return value
to make slightly smaller/faster code.
Enhanced the JavaScript function Astronomy.SearchMoonPhase
to allow searching forward in time when the `limitDays`
argument is positive, or backward in time when `limitDays`
is negative.
Added unit test "moon_reverse" to verify this new feature.
