Did you by any chance, happen to look at my code for reference? Or are they somewhat similar just out of coincidence and similar naming conventions? If you did, in fact, look at my project for reference, I'm glad to have helped!
Just for comparison, here's the orbit determination from initial conditions code from a newer project of mine. It doesn't compute the classical orbital directly, because that tends to add numerical error (ie. first you compute inclination using acos() and a moment later you compute the orientation matrix using cos(inclination)).
I don't remember looking at this project. I think I did read other implementations of this algorithm, and they all end up very similar, so it's probably due to algorithm itself.
It might be interesting to link to all codebases that might be relevant. Maybe I will add a list to https://airsick.guide/tools.html.
---
I only searched briefly for more stable methods. Thanks for the reference, I'll have a look!
Yeah, there aren't too many ways of detemining orbital elements from initial conditions so it is very likely that the code implementing that ends up looking the same. The similarity in variable names and naming_conventions made the two functions look rather similar.
What comes to stability: given just a single pair of (position, velocity) vectors, it's hard to make something that is really stable. If you have several pairs, you can use the method of least squares to try to fit an orbit into observations.
Not using too many arcus trig functions (acos, asin), will give a bit better results near zero.
Additionally, not having to solve the Kepler's equation would be an improvement for near-parabolic trajectories. In other words: it might be preferable to store (true anomaly at epoch, timestamp of epoch) rather than time of periapsis passage or mean anomaly at epoch to avoid this. The issue is that E - e*sin(E) loses numerical precision when e is near 1 (parabolic) and E is near zero. Same applies for the hyperbolic equivalent of Kepler's equation.
Anyway, nice to find people with similar interests! Good luck with your projects.
Regarding the naming convention, I am trying to stick to PEP8 as strictly as possible (sometimes, even in C) and to use overly explicit variable names. I had not noticed that you did that too for the structure member, with `longitude_of_ascending_node` and such. My objective is to have the code be as readable as possible to a newcomer (hence Python, too).
You should consider using my libtwobody project as a "backend" for your Python library. That is actually one of the design goals for the library, to provide low level primitive operations for higher level tasks and be usable from "scripting" languages.
Advantages over a simple implementation such as yours: proper handling of parabolic and near parabolic orbits, state of the art solvers for time of flight (Kepler's equation) and universal variable formulation. Additionally I've done some work with higher level algorithms for interplanetary trajectory planning.
The library may seem a bit cryptic if you're not too familiar with the literature on this subject, but it should not be too difficult to get up to speed with it. It might take a little bit of trivial software engineering to make it work to get it working as a dll but shouldn't be too much work. I'm happy to help, drop me an email (address in commit log) or contact me on irc (freenode, same nick as here on hn).
Looks somewhat similar to what I did a few years ago in a now-abandoned project of mine: https://github.com/rikusalminen/libkepler/blob/master/src/ke...
Did you by any chance, happen to look at my code for reference? Or are they somewhat similar just out of coincidence and similar naming conventions? If you did, in fact, look at my project for reference, I'm glad to have helped!
Just for comparison, here's the orbit determination from initial conditions code from a newer project of mine. It doesn't compute the classical orbital directly, because that tends to add numerical error (ie. first you compute inclination using acos() and a moment later you compute the orientation matrix using cos(inclination)).
https://github.com/rikusalminen/twobody/blob/master/include/...
You can, of course, get the classical elements from the orbit, they are useful for presenting human-readable numbers.