Skip to content
This repository was archived by the owner on Aug 28, 2024. It is now read-only.

Commit

Permalink
Improved the horizons query method to accept additional keyword argum…
Browse files Browse the repository at this point in the history
…ents to pass to the astroquery.Horizons method. This allows one to disambiguate names by passing 'id_type' as an optional keyword
  • Loading branch information
MintoDA1 authored and MintoDA1 committed Aug 2, 2023
1 parent 517a2f4 commit f9baf87
Show file tree
Hide file tree
Showing 3 changed files with 27 additions and 16 deletions.
2 changes: 1 addition & 1 deletion examples/Basic_Simulation/basic_simulation.py
Original file line number Diff line number Diff line change
Expand Up @@ -45,7 +45,7 @@
sim.add_solar_system_body(["Sun","Mercury","Venus","Earth","Mars","Jupiter","Saturn","Uranus","Neptune"])

# Add in some main belt asteroids
sim.add_solar_system_body(name=["Ceres","Vesta","Pallas","Hygiea"])
sim.add_solar_system_body(name=["Ceres","Vesta","Pallas","Hygiea"],id_type="smallbody")

# Add in some big KBOs
sim.add_solar_system_body(name=["Pluto","Eris","Sedna","Haumea","Makemake","Quaoar","Orcus","Gonggong","Salacia"])
Expand Down
34 changes: 21 additions & 13 deletions python/swiftest/swiftest/init_cond.py
Original file line number Diff line number Diff line change
Expand Up @@ -25,15 +25,18 @@
Any
)

def horizons_get_physical_properties(altid):
def horizons_get_physical_properties(altid,**kwargs):
"""
Parses the raw output from JPL Horizons in order to extract physical properties of a body if they exist
Parameters
----------
altid : string
Raw response from the JPL Horizons query
altid : list of str
List of ids to use for Horizons query
**kwargs: Any
Additional keyword arguments to pass to the query method (see https://astroquery.readthedocs.io/en/latest/jplhorizons/jplhorizons.html)
Returns
-------
Expand Down Expand Up @@ -104,7 +107,7 @@ def get_rotrate(raw_response):
else:
if "Synchronous" in raw_response: # Satellites have this:
rotrate = [s for s in raw_response.split('\n') if 'Orbital period' in s][0]
rotrate = rotrate.split('Orbital period')[1].replace('~',' ').replace('d',' ').strip()
rotrate = rotrate.split('Orbital period')[1].replace('~',' ').replace('d',' ').replace('=',' ').strip()
rotrate = 2*np.pi / (float(rotrate) * swiftest.JD2S)
else:
rotrate = None
Expand All @@ -125,9 +128,12 @@ def get_rotpole(jpl):

rotpole = SkyCoord(ra=RA * u.degree, dec=DEC * u.degree).cartesian
return np.array([rotpole.x.value, rotpole.y.value, rotpole.z.value])

if type(altid) != list:
altid = [altid]

for id in altid:
jpl,idlist,namelist = horizons_query(id=id,ephemerides_start_date='2023-07-26',verbose=False)
jpl,idlist,namelist = horizons_query(id=id,ephemerides_start_date='2023-07-26',verbose=False,**kwargs)
Rpl = get_radius(jpl.raw_response)
if Rpl is not None:
Rpl *= 1e3
Expand All @@ -152,7 +158,7 @@ def get_rotpole(jpl):
return Gmass,Rpl,rot


def horizons_query(id, ephemerides_start_date, exclude_spacecraft=True, verbose=False):
def horizons_query(id, ephemerides_start_date, exclude_spacecraft=True, verbose=False,**kwargs):
"""
Queries JPL/Horizons for a body matching the id. If one is found, a HorizonsClass object is returned for the first object that
matches the passed id string. If more than one match is found, a list of alternate ids is also returned. If no object is found
Expand Down Expand Up @@ -219,7 +225,7 @@ def get_altid(errstr,exclude_spacecraft=True):
try:
jpl = Horizons(id=id, location='@sun',
epochs={'start': ephemerides_start_date, 'stop': ephemerides_end_date,
'step': ephemerides_step})
'step': ephemerides_step},**kwargs)
eph=jpl.ephemerides()
altid = [id]
altname =[jpl.table['targetname'][0]]
Expand Down Expand Up @@ -248,7 +254,8 @@ def get_altid(errstr,exclude_spacecraft=True):
def solar_system_horizons(name: str,
param: Dict,
ephemerides_start_date: str,
ephemeris_id: str | None = None):
ephemeris_id: str | None = None,
**kwargs: Any):
"""
Initializes a Swiftest dataset containing the major planets of the Solar System at a particular data from JPL/Horizons
Expand All @@ -263,7 +270,8 @@ def solar_system_horizons(name: str,
Date to use when obtaining the ephemerides in the format YYYY-MM-DD.
ephemeris_id : string (optional)
If passed, this is passed to Horizons instead of `name`. This can be used to find a more precise body than given by `name`.
**kwargs: Any
Additional keyword arguments to pass to the query method (see https://astroquery.readthedocs.io/en/latest/jplhorizons/jplhorizons.html)
Returns
-------
Expand Down Expand Up @@ -339,7 +347,7 @@ def solar_system_horizons(name: str,

print(f"Fetching ephemerides data for {ephemeris_id} from JPL/Horizons")

jpl,altid,altname = horizons_query(ephemeris_id,ephemerides_start_date)
jpl,altid,altname = horizons_query(ephemeris_id,ephemerides_start_date,**kwargs)
if jpl is not None:
print(f"Found ephemerides data for {altname[0]} ({altid[0]}) from JPL/Horizons")
if name == None:
Expand All @@ -365,16 +373,16 @@ def solar_system_horizons(name: str,
omega = jpl.elements()['w'][0]
capm = jpl.elements()['M'][0]

Gmass,Rpl,rot = horizons_get_physical_properties(altid)
Gmass,Rpl,rot = horizons_get_physical_properties(altid,**kwargs)
# If the user inputs "Earth" or Pluto, then the Earth-Moon or Pluto-Charon barycenter and combined mass is used.
# To use the Earth or Pluto alone, simply pass "399" or "999", respectively to name
if ephemeris_id == "Earth":
print("Combining mass of Earth and the Moon")
Gmass_moon,tmp,tmp = horizons_get_physical_properties(["301"])
Gmass_moon,tmp,tmp = horizons_get_physical_properties(["301"],**kwargs)
Gmass += Gmass_moon
elif ephemeris_id == "Pluto":
print("Combining mass of Pluto and Charon")
Gmass_charon,tmp,tmp = horizons_get_physical_properties(["901"])
Gmass_charon,tmp,tmp = horizons_get_physical_properties(["901"],**kwargs)
Gmass += Gmass_charon

if Gmass is not None:
Expand Down
7 changes: 5 additions & 2 deletions python/swiftest/swiftest/simulation_class.py
Original file line number Diff line number Diff line change
Expand Up @@ -2165,7 +2165,8 @@ def add_solar_system_body(self,
name: str | List[str] | None = None,
ephemeris_id: int | List[int] | None = None,
date: str | None = None,
source: str = "HORIZONS"):
source: str = "HORIZONS",
**kwargs: Any):
"""
Adds a solar system body to an existing simulation Dataset from the JPL Horizons ephemeris service. The JPL Horizons service
will be searched for a body matching the string passed by `name`, or alternatively `ephemeris_id` if passed. Bodies will be
Expand Down Expand Up @@ -2193,6 +2194,8 @@ def add_solar_system_body(self,
source : str, default "Horizons"
The source of the ephemerides.
>*Note.* Currently only the JPL Horizons ephemeris is implemented, so this is ignored.
**kwargs: Any
Additional keyword arguments to pass to the query method (i.e. astroquery.Horizons)
Returns
-------
None
Expand Down Expand Up @@ -2232,7 +2235,7 @@ def add_solar_system_body(self,

body_list = []
for i,n in enumerate(name):
body = init_cond.solar_system_horizons(n, self.param, date, ephemeris_id=ephemeris_id[i])
body = init_cond.solar_system_horizons(n, self.param, date, ephemeris_id=ephemeris_id[i],**kwargs)
if body is not None:
body_list.append(body)

Expand Down

0 comments on commit f9baf87

Please sign in to comment.