The test build failed because diffcalc reported a small
discrepancy between the C and C# output.
So I made the threshold more lenient for now.
I want to come back later and figure out if I can get back
to exact agreement between C and C# code.
Told wget not to output rediculous progress bar stuff
that eats thousands of lines of log output.
It turns out different Node.js versions do math differently,
which caused a Travis CI build failure.
Scale topocentric distance the same way I scale heliocentric distance.
Adjusted diffcalc bash script and diffcalc.bat Windows batch file accordingly.
The differ now prints the final "score" so I'm less likely to make
a mistake spotting the correct maximum difference.
Removed unused variable in ctest.c DiffLine(): maxdiff.
When comparing calculations of body vectors, scale
the size of the difference by the minimum orbital
radius (or typical radius in the case of the Solar
System Barycenter).
This concludes my investigations of discrepancies between
the various language calculations. I have done as much
as I can without implementing my own trig functions,
which is not worth the effort (or the loss of efficiency
in JavaScript).
Scaling the errors relative the measurement units reveals
that the discrepancies are reasonable for the 16-digit
precision one expects from 64-bit floating point numbers.
The worst case is C vs JavaScript, with a scaled error
of about 7.2e-15. I can live with that.
A given amount of error in an angle measured in
sidereal hours is 15 times more important than the
same numeric error in an angle measured in degrees.
Scale angular errors by the range of values they
could take on. Longitude-like angles in degrees
have a range of 360, while latitude-like angles
range over 180 degrees (-90 to +90).
Split out separate Windows batch file diffcalc.bat,
just like I already split out bash script diffcalc.
I was able to get the diffcalc test to confirm that
the C and C# algorithms are producing absolutely identical
output values. I just needed to make the C code print its
output in scientific notation with 18 significant figures.
For angular results that are longitude-like, the severity
of the error decreases as the location approaches the
poles of that coordinate system. So scale the longitude
difference by the cosine of the latitude.
For example, in horizontal coordinates, the closer an object
is to the zenith, the less significant an error in its azimuth
becomes. At the zenith, azimuth loses all meaning.
Also, now we pass in the tolerance on the command line, because
there are known issues for each language. I expect C# and C
to be in very close agreement, and I want to know if they start
drifting apart. However, JS has known issues with atan2(), so I'm more
lenient with C vs JS.
Now track and report max difference independently for each
type of calculated value, instead of merging them all
together in one big soup.
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.
