hannorein / rebound

đŸ’« An open-source multi-purpose N-body code.
https://rebound.readthedocs.io/
GNU General Public License v3.0
824 stars 218 forks source link

omega and pomega values and others that are dependent on inclination angles #603

Closed zhexingli closed 2 years ago

zhexingli commented 2 years ago

Hello,

I'm trying to plot out some orbital elements of secondary bodies using output from REBOUND but some of the parameter values returned are not what I'd expect. For example, my pomega and omega values appears to be mirrored by the Omega value when I set inclination to be edge-on (i=90, or pi/2) and Omega = 0.

Looking at the REBOUND code, it appears that there's a sign flip for angles like this (other angles as well such as for mean anomaly and mean longitude etc): when math.cos(i) > 0, we have pomega = omega + Omega, otherwise, we have pomega = Omega - omega. However, with i = pi/2, my laptop gives math.cos(i) a near zero but positive value, so I'd expect REBOUND to take the first equation when calculating pomega (and also other angles that have conditions like this being checked) and gives pomega = omega, but from my plot it appears it's using the second one. This happens to any of my case that uses i=pi/2. May I ask why the edge-on inclination case is not included as a condition for the first expression? And how should I treat edge-on orbit cases like this?

Thanks

hannorein commented 2 years ago

These edge cases are tricky because some quantities might not be well defined. The implementation in REBOUND should be stable and self-consistent (i.e. you can go back and forth from cartesian to orbital elements and get the same result). But in many cases it might make sense to change your coordinate system. For example, you could just rotate your coordinate system so that all "edge on" cases don't have i=90 but i=0. Alternatively, you could use the ix, iy, h, k variables for this.

dtamayo commented 2 years ago

Agreed! Pomega only really makes sense when nearly coplanar, and you want it to approach the angle from the x axis to pericenter, so for prograde near coplanar it should be Omega + omega but for retrograde near i=180 it needs to be Omega - omega. This confused me a lot with retrograde irregular satellites. I agree it's strange that it's discontinuous at i=90, but I think that's always been the traditional definition of pomega

On Thu, Jul 7, 2022, 1:16 PM Hanno Rein @.***> wrote:

These edge cases are tricky because some quantities might not be well defined. The implementation in REBOUND should be stable and self-consistent (i.e. you can go back and forth from cartesian to orbital elements and get the same result). But in many cases it might make sense to change your coordinate system. For example, you could just rotate your coordinate system so that all "edge on" cases don't have i=90 but i=0. Alternatively, you could use the ix, iy, h, k variables for this.

— Reply to this email directly, view it on GitHub https://github.com/hannorein/rebound/issues/603#issuecomment-1177516572, or unsubscribe https://github.com/notifications/unsubscribe-auth/AA2ABFYAEDKQLE6SM46U3LLVS3DBXANCNFSM523VPFFA . You are receiving this because you are subscribed to this thread.Message ID: @.***>

zhexingli commented 2 years ago

Thank you Hanno and Dan for the clarification! One more question, how's i defined in REBOUND? It's not quite clear from the documentation page, is i=90 edge-on for face-on? It probably doesn't matter for a lot of the cases but does it when dealing with inclination sensitive works? Is the observer at the reference x direction or at z direction if we use the orbit plot in the documentation page? Or is it really the matter or picking the right axis to work with?

hannorein commented 2 years ago

If i=0, then the orbit is in the xy plane. Whether you call this edge-on or face-on just depends on where you place your observer.

zhexingli commented 2 years ago

Ok thanks!