mirror of
https://github.com/cosinekitty/astronomy.git
synced 2026-02-05 05:01:55 -05:00
cbcacc4b57
Before making these changes, I had the following discrepancies
between the calculations made by the different programming
language implementations of Astronomy Engine:
C vs C#: 5.55112e-17, worst line number = 6
C vs JS: 2.78533e-12, worst line number = 196936
C vs PY: 1.52767e-12, worst line number = 159834
Now the results are:
Diffing calculations: C vs C#
ctest(Diff): Maximum numeric difference = 5.55112e-17, worst line number = 5
Diffing calculations: C vs JS
ctest(Diff): Maximum numeric difference = 1.02318e-12, worst line number = 133677
Diffing calculations: C vs PY
ctest(Diff): Maximum numeric difference = 5.68434e-14, worst line number = 49066
Diffing calculations: JS vs PY
ctest(Diff): Maximum numeric difference = 1.02318e-12, worst line number = 133677
Here is how I did this:
1. Use new constants HOUR2RAD, RAD2HOUR that directly convert between radians and sidereal hours.
This reduces tiny roundoff errors in the conversions.
2. In VSOP longitude calculations, keep clamping the angular sum to
the range [-2pi, +2pi], to prevent it from accumulating thousands
of radians. This reduces the accumulated error in the final result
before it is fed into trig functions.
The remaining discrepancies are largely because of an "azimuth amplification" effect:
When converting equatorial coordinates to horizontal coordinates, an object near
the zenith (or nadir) has an azimuth that is highly sensitive to the input
equatorial coordinates. A tiny change in right ascension (RA) can cause a much
larger change in azimuth.
I tracked down the RA discrepancy, and it is due to a different behavior
of the atan2 function in C and JavaScript. There are cases where the least
significant decimal digit is off by 1, as if due to a difference of opinion
about rounding policy.
My best thought is to go back and have a more nuanced diffcalc that
applies less strict tests for azimuth values than the other calculated values.
It seems like every other computed quantity is less sensitive, because solar
system bodies tend to stay away from "poles" of other angular coordinate
systems: their ecliptic latitudes and equatorial declinations are usually
reasonably close to zero. Therefore, right ascensions and ecliptic longitudes
are usually insensitive to changes in the cartesian coordinates they
are calculated from.
7 lines
460 B
Plaintext
7 lines
460 B
Plaintext
Moon horizontal position: azimuth = 274.486, altitude = 52.101
|
|
Moon check: (1.0, -2.072588670409689e-16, -1.0362943352048446e-16, 2021-03-22T02:45:00.000Z)
|
|
Sun vector: (-0.08853032065573567, -0.31939127849027193, -0.9396631499852881, 2021-03-22T02:45:00.000Z)
|
|
Tilt angle of sunlit side of the Moon = -108.773 degrees counterclockwise from up.
|
|
Moon magnitude = -10.27, phase angle = 84.22 degrees.
|
|
Angle between Moon and Sun as seen from Earth = 95.64 degrees.
|