How to get planet centric orbits into REBOUND using sim.add class to JPL Horizons?

[ijiraq]

Why doesn't this return the inclination of Titan's orbit relative to the Saturn equator?

import rebound 
sim = rebound.Simulation()
sim.units = ('AU', 'days', 'Msun')
sim.add("699", hash='Saturn')
sim.add("606", hash='Titan', primary=sim.particles['Saturn'])
print(sim.particles['Titan'].calculate_orbit())

<rebound.Orbit instance, a=0.008168242278988228 e=0.02862878656499152 inc=0.48349150583342243 Omega=2.950727393359049 omega=3.0584123778020955 f=-3.2800042295711433>

Should it?

[hannorein]

By default REBOUND using the ecliptic coordinate plane when initializing particles with NASA Horizon data (see NASA Horizons for the precise definition). You need to find the appropriate coordinate transformation from the ecliptic to the Saturnian equatorial plane.

Note that since REBOUND treats all particles as point particles (unless you implement anything else, e.g. using REBOUNDx), it doesn't keep track of a particle's spin and therefore can't keep track of a particle's equator either.

[ijiraq]

Thanks for the reply and confirming that the orbit coming back is in ecliptic coordinates rather than Saturn Equatorial, which makes sense in some ways. I'll adjust my inputs to account for that. (as you noted I could specify that Saturn is the primary for the orbit computation but that is also the default in that particular case). Cheers.

[hannorein]

Yes. That's the right approach (I think NASA Horizons will do the transformation for you). REBOUND is agnostic about the coordinate system, so whatever you use for the initial conditions is also what you get back out.