Corrected a mistake in the explanation of the
C function Astronomy_GravSimInit: the `bodyStates`
parameter is NOT barycentric -- it is relative to the
originBody parameter.
Python had improperly formatted documentation for
Time.FromTerrestrialTime parameter `tt`.
The Python markdown generator `pydown` did not
correctly handle links to compound symbols like
`#GravitySimulator.Update`. It also was trying
to link to `StateVector[]` instead of `StateVector`.
Removed unnecessary and unhelpful documentation
for C# internal class constructors. They do not appear
in the generated markdown documentation anyway.
Other minor wording revisions in the documentation.
Finished coding the Python version of the gravity simulator.
No unit tests have been written yet.
Cleaned up documentation in the other languages.
Made some functions static that did not need to be members.
Added the following functions:
GravSimTime
GravSimNumBodies
GravSimOrigin
GravSimSwap
The GravSimSwap function allows an instantaneous "undo"
of a simulation step, which required refactoring the internal
data structures of the simulator. I did this because I realized
I needed a way to undo exploration of time steps near a fixed
current time. This will make it easier to implement a
light travel time solver.
I have decided the ability to select different
collections of gravitating bodies causes far
more complexity in the code than it is worth.
So now the gravity simulator always calculates
the Sun and all planets except Pluto.
This greatly simplifies the core code, gets
a good balance between efficiency and accuracy,
and makes the test matrix much simpler.
The GravSimOriginState function was doing too
much work for bodies that weren't already calculated
and cached. Instead of always calling Astronomy_BaryState,
it is possible to do a lot less work for most bodies
by recycling the SSB calculation that has already
been done.
Other minor code cleanup -- mostly stuff that makes
the code easier to read and understand.
Starting implementation of a generalized gravity simulator.
Already calculating the movement of Ceres, but with less
accuracy than I had hoped. I don't know if the lack of modeling
pull of the other asteroids has a larger effect than I expected,
or there is just something wrong with the implementation.
It makes more sense to report Jupiter's moons with
individually named structure fields rather than an array.
It reduces the overall code and documentation size,
and outside of unit testing, there are few cases
where iterating over an array of moons is more
lucid than using the names of the moons.
This is a breaking change, but hopefully very few
developers are using this function yet.
Fixing the breakage is very simple.
The existing lunar libration functions in the
other languages (C, C#, Python, JavaScript) were
calculating the Moon's ecliptic latitude and longitude
in radians, not degrees as intended. They have been fixed.
Implemented the libration function for Kotlin.
For years before 1582 or years after 3668, the Seasons functions
were unable to find many equinoxes and/or solstices.
The problem was that over time, the Earth's axis precesses
enough that the calendar dates of these events drifts outside
the fixed search ranges I had provided for them.
I expanded the search ranges so all season changes can be found
for a much wider range of years, as verified by unit tests:
C/C++: -2000..9999
C#: 1..9999
JavaScript: -2000..9999
Python: 1..9999
Kotlin: 1..9999
Note: C#, Python, and Kotlin currently do not allow
years values below +1. In fact, I discovered we were not
noticing when an invalid year was passed into the Kotlin code.
I updated that code to throw an exception when the year does
not match what was expected. It is disturbing that the
GregorianCalendar class silently ignores invalid years!
Constricted the search tolerance from 1 second to 0.01
seconds for the seasons search, to ensure more consistent
behavior.
Fixed a bug in the Kotlin search() function's
quadratic interpolation that was causing the convergence
to be slower than it should have been.
Implemented Kotlin functions and code generator
for calculating the state vectors of Jupiter's
largest 4 moons.
Added cautionary comments about needing to correct
Jupiter's moons for light travel time.
This is the first pass to get everything needed
for the AstroCheck tests. I tried comparing
C output to Kotlin output, and there are some
serious problems to debug:
$ ./ctest diff 2.8e-16 temp/{c,k}_check.txt
First file: temp/c_check.txt
Second file: temp/k_check.txt
Tolerance = 2.800e-16
lnum a_value b_value factor diff name
FAIL 137746 4.2937184148112564e+01 4.2944101081740065e+01 0.03364 2.327e-04 helio_x
FAIL 373510 1.4197190315274938e+01 1.4193716564905307e+01 0.03364 1.168e-04 helio_y
FAIL 137746 -6.5897675150466091e+00 -6.5929481589493522e+00 0.03364 1.070e-04 helio_z
FAIL 59150 1.8035183339348251e+01 1.8035909197904104e+01 0.01730 1.255e-05 sky_j2000_ra
FAIL 137747 -8.1222057639092533e+00 -8.1250990689970894e+00 0.00556 1.607e-05 sky_j2000_dec
FAIL 137747 4.8436159305823310e+01 4.8441487614058218e+01 0.03481 1.855e-04 sky_j2000_dist
FAIL 322846 8.7596368704201495e+01 2.6760770774700188e+02 0.00278 4.995e-01 sky_hor_az
FAIL 405828 -6.5075824596574279e+01 5.6922941329250996e+01 0.00556 6.778e-01 sky_hor_alt
OK 92717 4.1268347083494783e-03 4.1268347083494774e-03 223.21429 1.936e-16 jm_x
OK 45091 -8.0149190392649894e-03 -8.0149190392649929e-03 79.42812 2.756e-16 jm_y
OK 135377 1.5470777280065808e-03 1.5470777280065804e-03 223.21429 9.680e-17 jm_z
OK 216836 4.5725777238332412e-03 4.5725777238332394e-03 126.58228 2.196e-16 jm_vx
OK 351647 5.1351566793199944e-03 5.1351566793199962e-03 126.58228 2.196e-16 jm_vy
OK 351647 2.5217607180929289e-03 2.5217607180929298e-03 126.58228 1.098e-16 jm_vz
Score = 6.778e-01
ctest(Diff): EXCEEDED ERROR TOLERANCE.
So I'm checking this in as work-in-progress.
Ported the following types to the Kotlin code:
GlobalSolarEclipseInfo
EclipseEvent
LocalSolarEclipseInfo
TransitInfo
ShadowInfo
IllumInfo
AxisInfo
NodeEventKind
NodeEventInfo
Made some wording fixes in the documentation for the
other languages.
Converting between radians and degrees.
Clamping angles to a desired range of degrees.
Converting between vector, spherical, horizontal.
Refraction and inverse refraction.
Implemented most of the RotationMatrix functions.
Added unit tests for combining rotation matrices, using a
rotation matrix to rotate a vector, and pivoting a rotation
matrix around its axes.
Replaced AstroVector operator '*' with infix function 'dot',
because it removes ambiguity between vector dot products
and vector cross products.
Later I will add a 'cross' infix function too.
Corrected minor typo in documentation for Python, C, C#, JavaScript.
"trasnform" -> "transform"
There was already an internal function for calculating
Greenwich Apparent Sidereal Time (GAST). By request,
I have exposed this function for outside users.
Added a minimal unit test to verify the function is
callable and returns the correct result for one case.
This function is already exhaustively tested by unit
tests that verify other functions that already called
this function when it was internal, so minimal testing
is sufficient in this case.
Added the following new functions to all 4 languages:
MassProduct: find the GM product for all Solar System bodies.
LagrangePoint: calculate L1..L5 state vectors for a pair of bodies.
LagrangePointFast: calculate L1..L5 state vectors given
state vectors and GM products of a pair of bodies.
In languages that support it, using hypot(x,y) is a little
easier to read than sqrt(x*x + y*y). Some documentation
(e.g. the man page for the C function) leads me to believe
hypot might also be better behaved than sqrt in some cases.
The JavaScript Math.hypot() is especially nice because it works
for any number of dimensions, so I can use it in 2D and 3D cases.
C only allows 2D usage, as does Python 3.7. Python 3.8 added
support for any number of dimensions, but I don't want to break
compatibility with Python 3.7 just yet. Therefore, in C and Python,
I am only using hypot for 2D cases.
C# does not appear to have any kind of hypot function,
so no changes were made to the C# code.
Thanks to https://github.com/ebraminio for this suggestion.
My custom Markdown documentation generator for C had
a bug when emitting the listing of a #define.
It is not valid to try to hyperlink to other symbols,
because the Markdown syntax gets listed literally inside
the context of a C code block.
In most cases, people calculating Lagrange points will just
want to pass in the bodies and not have to worry about calculating
their state vectors and masses.
Renamed Astronomy_LagrangePoint to Astronomy_LagrangePointFast.
Added new function Astronomy_LagrangePoint that accepts body enum
values instead of state vectors and masses. It knows to optimize
the precision of the calculation by calling GeoMoonState for the
Earth/Moon case.
It is conceptually simpler to take cross products to
generate 3 coordinate axes (essentially a rotation matrix)
that represent radial, tangential, and normal directions
with respect to the major and minor bodies.
Now correctly calculating L4 and L5 positions, but
there is a large error in their velocity vectors.
Refactored ctest.c LagrangeTest() to be a lot easier
to understand and modify. A new function VerifyStateLagrange()
allows passing test parameters in a more function-oriented way.
Confirmed that L4 and L5 always lie in the same plane with
the position vector and velocity vector.
Use the formulas I already had to calculate first
approximations for L1, L2, L3 distances.
Then use Newton's Method to home in on the positions
where centrifugal acceleration balances with net
gravitational acceleration.
I realized there was a small mistake in how I was
calculating the distance scaling factor for L1 and L2.
It was relative to the distance between the minor body
and the barycenter, not the minor body and the major body.
This significantly improves the accuracy for Earth/Moon
Lagrange points, but still has more error compared
to JPL Horizons than I currently understand.
The Microsoft C compiler is oddly picky about declaring a const variable.
Apparently it cannot do math with other const variables in its
initializer expression, unlike other C compilers.
So I had to change MOON_GM from a const to a #define.
The Lagrange point calculation is still not finished,
but L1 and L2 are working. L3 is probably correct, but there
is no test data for it.
I replaced the test data with new JPL Horizons output that
is centered on the primary body instead of the Solar System Barycenter.
This allows Astronomy_LagrangePoint() to be agnostic about
the coordinate systems of the state vectors handed to it.
I still need to get L4 and L5 calculations to match JPL Horizons
data, but it is not yet clear how to do that.
Implemented a pair of C functions for finding a series of
Moon nodes:
Astronomy_SearchMoonNode
Astronomy_NextMoonNode
Finished the C unit test "moon_nodes" that verifies
my calculations against Fred Espenak's test data.
This is a thin wrapper function for the internal
function CalcMoon, which has already been extensively
validated. It will enable outside users to search
for ascending and descending nodes of the Moon,
or to calculate ecliptic spherical coordinates for the Moon
for any other useful purpose.
Changed the documentation for the GeoMoon and GeoMoonState
functions to make it explicit that they calculate coordinates
oriented with respect to the Earth's J2000 equator (EQJ).
This is because I will soon add ecliptic (ECL) counterparts
for the GeoMoon function, to more directly search for ascending
and descending nodes of the Moon.
See this discussion:
https://github.com/cosinekitty/astronomy/issues/150
For the case of calculating a map, where each pixel
on the map represents a different location on the Earth,
it is more efficient to factor out expensive calculation
of sidereal times, assuming the entire map represents
some phenomenon at a single moment in time.
For example, to determine whether the Moon is visible
at different places on the Earth, the following
functions can be calculated across thousands of
different (lat, lon) geographic coordinates around
the world:
ObserverVector
Rotation_EQD_HOR
Before iterating over the map pixels, a program
can call GeoMoon, then convert EQJ coordinates to EQD.
Then by passing the same time value in a loop to
ObserverVector and Rotation_EQD_HOR, the program
can calculate a vector from the observer to the Moon
in EQD coordinates, then convert EQD to HOR.
The z-coordinate of the horizontal coordinates
determines whether the Moon is above or below the
observer's horizon at that point on the Earth.
This calculation pattern performed redundant
sidereal time calculations for each pixel on the map.
I changed the code for all 4 languages to cache
sidereal time so that it only needs to be calculated
once.
In the C version of Astronomy Engine, this resulted
in a speedup factor of about 2.3 in the above use case.
(See the function MapPerformanceTest in generate/ctest.c.)
Reduce the number of redundant Earth nutation calculations
by passing astro_time_t values as pointers in more functions.
Nutation values can then be cached in the time parameter
and passed to other functions that can then avoid calculating
the same nutation again.
Nutation is an expensive calculation, so reducing this overhead
can dramatically speed up certain use cases.
This was only needed in C, because this is the only language
in which times are passed by value. In Python, C#, and JavaScript,
times are objects that are already passed by reference, and
they already benefit from this nutation recyling approach.
The following functions have had their parameters changed.
This is a breaking change, but in every case, the caller
usually just needs to change `time` to `&time`.
Astronomy_Rotation_EQD_EQJ
Astronomy_Rotation_EQD_ECL
Astronomy_Rotation_EQD_HOR
Astronomy_Rotation_EQJ_EQD
Astronomy_Rotation_EQJ_HOR
Astronomy_Rotation_ECL_EQD
Astronomy_Rotation_ECL_HOR
Astronomy_Rotation_HOR_EQD
Astronomy_Rotation_HOR_EQJ
Astronomy_Rotation_HOR_ECL
Astronomy_RotationAxis
Astronomy_VectorObserver
The code generator was creating slightly different numeric
values for the Pluto state tables and the Jupiter rotation matrix.
I decreased the output precision by one decimal digit.
This should allow the code generator to produce identical
source code on both Linux and macOS.
Added radius data for the Sun, Moon, and remaining planets.
Test the raytracer for all other bodies except the Earth and Sun.
There is a problem with Pluto that I still need to figure out.
Fixed an issue in the doxygen-to-markdown translator I wrote
(hydrogen.js): it did not handle when one #define referred
to another #define. Created a more generic markdown expansion
that works in all cases, and creates embedded hyperlinks.
Refactored the Jupiter imager to be a generic planet imager.
Added support for drawing an image of Venus.
Verified that its extreme crescent phase looks correct
for the current date.
I will add radius constants to astronomy.h for each body I support.
