hannorein / rebound

💫 An open-source multi-purpose N-body code.
https://rebound.readthedocs.io/
GNU General Public License v3.0
862 stars 221 forks source link

JPL's HORIZONS query #530

Closed AstroEloy closed 3 years ago

AstroEloy commented 3 years ago

I have noticed that if I do a query with rebound for a specific date and the whole solar system, the position and velocity of the earth are different from those returned directly by Horizons.

What is the reason for this? Is it because of the reference system? Is there any way to fix it and have Rebound get exactly the same answer from Horizons?

from astroquery.jplhorizons import Horizons
import rebound

time = '2015-03-30T21:33:52.000'
JD = Time(time).jd

#JPL'S HORIZONS
obj = Horizons(id="Geocenter", location="@sun", epochs=JD, id_type='id').vectors()
print(obj['x']*1.496e+8, obj['y']*1.496e+8, obj['z']*1.496e+8) #km

#REBOUND
sim = rebound.Simulation()
date = "JD"+str(JD)
sim.units = ('s', 'km', 'kg')

sim.add("Sun", date=date)
sim.add("Mercury", date=date)
sim.add("Venus", date=date)
sim.add('Geocenter', date=date)
sim.add("Luna", date=date)
sim.add("Mars", date=date)
sim.add('Jupiter', date=date)
sim.add('Saturn', date=date)
sim.add('Uranus', date=date)
sim.add('Neptune', date=date)

sim.N_active = sim.N
print(sim.status())
hannorein commented 3 years ago

Hm. Both astroquery and REBOUND use J2000 as the default reference frame. But REBOUND uses @0 instead of @sun. Any chance that's the difference?

AstroEloy commented 3 years ago

Thanks! that works!

Can I use @sun with Rebound instead of @0 ?

hannorein commented 3 years ago

Regarding @sun. You'd either need to change this manually in the horizons.py file (just search for @0), or once you've imported the particle data, just move to the heliocentric frame with sim.move_to_hel().

AstroEloy commented 3 years ago

It works! Thank you very much!