API Documentation
The following tables show what functions are supported for what kind of orbit. If you're not sure yet what kind of orbit to use, just use the orbit
function! ✅ indiciates that a function is available, and ❌ indicates it is not due to the orbit not storing sufficient information. ⚠️ indicates that it could be supoprted, but is not yet implemented.
Required Parameters
The following table specifies what properties are required to construct each orbit type. Based on this information, different orbit types have different capabilities (described in following tables).
property | meaning | KepOrbit | Visual{KepOrbit} | ThieleInnesOrbit | RadialVelocityOrbit | CartesianOrbit | Visual{CartesianOrbit} |
---|---|---|---|---|---|---|---|
M | ✔️ | ✔️ | ✔️ | ✔️ | ✔️ | ✔️ | |
τ | ✔️ | ✔️ | ✔️ | ✔️ | ✔️ | ✔️ | |
tref | ✔️ | ✔️ | ✔️ | ✔️ | ✔️ | ✔️ | |
e | ✔️ | ✔️ | ✔️ | ✔️ | |||
i | ✔️ | ✔️ | |||||
ω | ✔️ | ✔️ | ✔️ | ||||
Ω | ✔️ | ✔️ | |||||
A | ✔️ | ||||||
B | ✔️ | ||||||
F | ✔️ | ||||||
G | ✔️ | ||||||
x | ✔️ | ✔️ | |||||
y | ✔️ | ✔️ | |||||
z | ✔️ | ✔️ | |||||
vx | ✔️ | ✔️ | |||||
vy | ✔️ | ✔️ | |||||
vz | ✔️ | ✔️ |
Properties of Orbits
You can use these functions like totalmass(orbit)
.
Function | KepOrbit | Visual{KepOrbit} | ThieleInnesOrbit | RadialVelocityOrbit | CartesianOrbit | Visual{CartesianOrbit} |
---|---|---|---|---|---|---|
totalmass | ✅ | ✅ | ✅ | ✅ | ✅ | ✅ |
period | ✅ | ✅ | ✅ | ✅ | ✅ | ✅ |
distance | ✅ | ❌ | ✅ | ❌ | ❌ | ✅ |
meanmotion | ✅ | ✅ | ✅ | ✅ | ✅ | ✅ |
eccentricity | ✅ | ✅ | ✅ | ✅ | ✅ | ✅ |
inclination | ✅ | ✅ | ✅ | ❌ | ✅ | ✅ |
semimajoraxis | ✅ | ✅ | ✅ | ✅ | ✅ | ✅ |
periastron | ✅ | ✅ | ✅ | ✅ | ✅ | ✅ |
semiamplitude | ✅ | ✅ | ⚠️ | ✅ | ✅ | ✅ |
Properties of Orbit Solutions
You can use these functions like sol = orbitsolve(orbit,mjd("2020-01")); posx(sol)
, or like posx(sol, mjd("2020-01"))
.
Function | KepOrbit | Visual{KepOrbit} | ThieleInnesOrbit | RadialVelocityOrbit | CartesianOrbit | Visual{CartesianOrbit} |
---|---|---|---|---|---|---|
meananom | ✅ | ✅ | ✅ | ✅ | ✅ | ✅ |
trueanom | ✅ | ✅ | ✅ | ✅ | ✅ | ✅ |
eccanom | ✅ | ✅ | ✅ | ✅ | ✅ | ✅ |
posx | ✅ | ✅ | ✅ | ❌ | ✅ | ✅ |
posy | ✅ | ✅ | ✅ | ❌ | ✅ | ✅ |
posz | ✅ | ✅ | ✅ | ❌ | ✅ | ✅ |
velx | ✅ | ✅ | ✅ | ❌ | ✅ | ✅ |
vely | ✅ | ✅ | ✅ | ❌ | ✅ | ✅ |
velz | ✅ | ✅ | ✅ | ✅ | ✅ | ✅ |
raoff | ❌ | ✅ | ✅ | ❌ | ❌ | ✅ |
decoff | ❌ | ✅ | ✅ | ❌ | ❌ | ✅ |
radvel | ✅ | ✅ | ✅ | ✅ | ✅ | ✅ |
posangle | ✅ | ✅ | ✅ | ❌ | ✅ | ✅ |
pmra | ❌ | ✅ | ✅ | ❌ | ❌ | ✅ |
pmdec | ❌ | ✅ | ✅ | ❌ | ❌ | ✅ |
accra | ❌ | ✅ | ❌ | ❌ | ❌ | ❌ |
accdec | ❌ | ✅ | ❌ | ❌ | ❌ | ❌ |
Documentation
PlanetOrbits.AbstractOrbit
— TypeAbstractOrbit
Represents a orbit. Contains all the information to calculate the location of a planet at a given time, true anomaly, eccentric anomaly, or mean anomaly. Different concrete implementations of AbstractOrbit contain varying amounts of information.
Basic information about the orbit can be queried using functions like period(orbit)
.
Orbits can be solved using functions like orbitsolve(orb)
.
See: RadialVelocityOrbit
, KepOrbit
, VisualOrbit
PlanetOrbits.AbstractOrbitSolution
— TypeAbstractOrbitSolution
Represents the solution of an orbit. Contains all the information of an AbstractOrbit, plus information necessary to uniquely locate a planet.
The solution can be queried using a variety of functions such as radvel(solution)
.
The API for creating orbit solutions it not considered public as the fields may change between minor versions. Instead, create solutions only through the public orbitsolve
and orbitsolve_...
functions.
PlanetOrbits.Auto
— TypePlanetOrbits.Auto()
Automatic choice of Kepler solver algorithm. Currently defaults to PlanetOrbits.Markley()
PlanetOrbits.CartesianOrbit
— MethodConvert an existing orbit object to a CartesianOrbit.
PlanetOrbits.Goat
— TypePlanetOrbits.Goat()
Kepler solver implementation from https://arxiv.org/abs/2103.15829 and https://github.com/oliverphilcox/Keplers-Goat-Herd
It is here for comparison purposes only. In general, Markley() is more performant and accurate.
PlanetOrbits.KepOrbit
— TypeKepOrbit(
a, # semi-major axis [AU]
e, # eccentricity
i, # inclination [rad]
ω, # argument of periapsis [rad]
Ω, # longitude of ascending node [rad]
τ, # epoch of periastron passage at MJD=0
M, # mass of primary [M⊙]
)
Represents the Keplerian elements of a secondary body orbiting a primary. Use the traditional Campbell parameterization. Values can be specified by keyword argument or named tuple for convenience.
PlanetOrbits.Markley
— TypePlanetOrbits.Markley()
Kepler solver implementation from AstroLib, based on Markley (1995) Celestial Mechanics and Dynamical Astronomy, 63, 101 (DOI:10.1007/BF00691917).
PlanetOrbits.OrbitSolutionCartesian
— TypeRepresents a CartesianOrbit
evaluated to some position.
PlanetOrbits.OrbitSolutionKep
— TypeRepresents a KepOrbit
evaluated to some position.
PlanetOrbits.OrbitSolutionRadialVelocity
— TypeRepresents a RadialVelocityOrbit
evaluated to some position.
PlanetOrbits.OrbitSolutionThieleInnes
— TypeRepresents a ThieleInnesOrbit
evaluated to some position.
PlanetOrbits.RadialVelocityOrbit
— TypeRadialVelocityOrbit(a, e, ω, τ, M)
Represents an orbit of a planet with only the information retrievable from radial velocity measurements. That is, without inclination, longitude of ascending node, or distance to the system.
PlanetOrbits.RootsMethod
— TypePlanetOrbits.RootsMethod(method::Roots.PlanetOrbits.Roots.AbstractUnivariateZeroMethod, kwargs...)
Wraps a root finding method from Roots.jl. Requires Roots to be loaded first. You can also pass keyword arguments that will be forwarded to Roots to control the tolerance.
Examples:
method = PlanetOrbits.RootsMethod(Roots.Newton())
method = PlanetOrbits.RootsMethod(Roots.Thukral5B())
method = PlanetOrbits.RootsMethod(Roots.Bisection())
method = PlanetOrbits.RootsMethod(Roots.A42())
method = PlanetOrbits.RootsMethod(Roots.Newton(), rtol=1e-3, atol=1e-3)
PlanetOrbits.ThieleInnesOrbit
— TypeThieleInnesOrbit(e, τ, M, plx, A, B, F, G)
Represents a visual orbit of a planet using Thiele-Innes orbital elements. Convertable to and from a VisualOrbit. This parameterization does not have the issue that traditional angular parameters have where the argument of periapsis and longitude of ascending node become undefined for circular and face on orbits respectively.
PlanetOrbits.Visual
— TypeVisual{OrbitType}(..., plx=...)
This wraps another orbit to add the parallax distance field plx
, thus allowing projected quantities to be calculated. It forwards everything else to the parent orbit.
For example, the KepOrbit type supports calculating x and y positions in AU. A Visual{KepOrbit} additionally supports calculating projected right ascension and declination offsets.
The ThieleInnesOrbit
type does not need to be wrapped in Visual
as it the Thiele-Innes constants are already expressed in milliarcseconds and thus it always requires a plx
value.
PlanetOrbits.VisualOrbit
— TypeVisual{OrbitType}(..., plx=...)
This wraps another orbit to add the parallax distance field plx
, thus allowing projected quantities to be calculated. It forwards everything else to the parent orbit.
For example, the KepOrbit type supports calculating x and y positions in AU. A Visual{KepOrbit} additionally supports calculating projected right ascension and declination offsets.
Note: the ThieleInnesOrbit
type does not need to be wrapped in Visual
as it the Thiele-Innes constants are already expressed in milliarcseconds and thus it always requires a plx
value.
PlanetOrbits.accdec
— Functionaccdec(orbit, t)
Get the instantaneous acceleration [mas/year^2] in the right-ascension direction of the secondary at the time t
[days].
accdec(o)
Get the instantaneous acceleration [mas/year^2] in the right-ascension direction of the secondary from an instance of AbstractOrbitSolution
.
accdec(elem, t, M_planet)
Get the instantaneous acceleration [mas/year^2] in the right-ascension direction of the primary in at the time t
[days]. The units of M_planet
and elem.M
must match.
accdec(o)
Get the instantaneous acceleration [mas/year^2] in the right-ascension direction of the primary from an instance of AbstractOrbitSolution
. The units of M_planet
and elem.M
must match.
PlanetOrbits.accra
— Functionaccra(orbit, t)
Get the instantaneous acceleration [mas/year^2] in the right-ascension direction of the secondary at the time t
[days].
accra(o)
Get the instantaneous acceleration [mas/year^2] in the right-ascension direction of the secondary from an instance of AbstractOrbitSolution
.
accra(elem, t, M_planet)
Get the instantaneous acceleration [mas/year^2] in the right-ascension direction of the primary in at the time t
[days]. The units of M_planet
and elem.M
must match.
accra(o)
Get the instantaneous acceleration [mas/year^2] in the right-ascension direction of the primary from an instance of AbstractOrbitSolution
. The units of M_planet
and elem.M
must match.
PlanetOrbits.astuple
— Methodastuple(elements)
Return the parameters of a KepOrbit value as a tuple.
PlanetOrbits.decoff
— Functiondecoff(orbit, t)
Get the offset [mas] from the primary body in Declination at the time t
[days].
decoff(orbit, t)
Get the offset [mas] from the primary body in Declination from an instance of AbstractOrbitSolution
.
PlanetOrbits.distance
— Functiondistance(orbit)
Distance to the system [pc].
PlanetOrbits.eccanom
— Methodeccanom(orbit, t)
Get the eccentric anomaly [radians] of the secondary at the time t
[days].
eccanom(o)
Get the eccentric anomaly [radians] of the secondary from an instance of AbstractOrbitSolution
.
PlanetOrbits.eccentricity
— Functioneccentricity(orbit)
Eccentricity of an orbit, between 0 and 1.
PlanetOrbits.inclination
— Functioninclination(orbit)
Inclination of an orbit, if available [rad].
PlanetOrbits.meananom
— Methodmeananom(orbit, t)
Get the mean anomaly [radians] of the secondary at the time t
[days].
meananom(o)
Get the mean anomaly [radians] of the secondary from an instance of AbstractOrbitSolution
.
PlanetOrbits.meanmotion
— Functionmeanmotion(orbit)
Mean motion [rad/year].
PlanetOrbits.mjd
— Methodmjd("2020-01-01")
Get the modfied julian day of a date, or in general a UTC timestamp.
PlanetOrbits.mjd
— Methodmjd()
Get the current modified julian day of right now.
PlanetOrbits.mjd2date
— Methodmjd2date(modified_julian)
Get a Date value from a modfied julian day, rounded to closest day
Examples
julia> mjd2date(59160.8)
2020-11-08
PlanetOrbits.orbit
— Methodorbit(...)
Construct an orbit from the provided keyword arguments. Will automatically select a subclass of AbstractOrbit based on the information provided. This is a convenience function that is not type stable and should not be used in performance sensitive contexts. Instead, call one of the concrete constructors KepOrbit
, VisualOrbit
, or RadialVelocityOrbit
directly. This function logs the kind of elements created so that it's easy to select the correct constructor.
Required arguments:
- a: semi-major axis [AU]
- M: mass of primary [M⊙]
Optional arguments:
- τ: epoch of periastron passage at MJD=0, default=0
- e: eccentricity, default=0
- ω: argument of periapsis [rad], default=0
- i: inclination [rad]
- Ω: longitude of ascending node [rad]
- plx: parallax [mas]; defines the distance to the primary
- tref=58849 [mjd]: reference epoch for τ
PlanetOrbits.orbitsolve
— Functionorbitsolve(orbit, t, method=Auto())
Given an orbit object and a time t
in days, get the position and velocity of the secondary body (e.g. planet around a star).
This will output a struct that is a subtype of AbstractOrbitSolution
which we can then query with raoff
, decoff
, radvel
, etc.
You can also calculate those quanitities individually (see their docstrings) but if you need more than one, it is most efficient to save the orbit solution once.
Note: these calculations use the small angle approximation, so are only accurate when the star is much further way from the observer than the secondary is from the primary.
See also: orbitsolve_ν
, orbitsolve_meananom
, orbitsolve_eccanom
, projectedseparation
, raoff
, decoff
, radvel
, propmotionanom
.
PlanetOrbits.orbitsolve_eccanom
— Methodorbitsolve_eccanom(elements, EA)
Same as orbitsolve
, but solves orbit for a given eccentric anomaly instead of time.
PlanetOrbits.orbitsolve_meananom
— Methodorbitsolve_meananom(elements, MA)
Same as orbitsolve
, but solves orbit for a given mean anomaly instead of time.
PlanetOrbits.orbitsolve_ν
— Functionorbitsolve_ν(elem, ν)
Solve a keplerian orbit from a given true anomaly [rad]. See orbitsolve for the same function accepting a given time.
PlanetOrbits.orbitsolve_ν
— Functionorbitsolve_ν(elem, ν, EA)
Solve an orbit from a given true anomaly [rad]. See orbitsolve
for the same function accepting a given time. Can optionally pass eccentric anomaly (EA) if already computed.
PlanetOrbits.periastron
— Functionperiastron(elements)
Compute the MJD of periastron passage most recently after the reference epoch tref specified in the orbit. N.B. mjd of 58849 = 2020-01-01
PlanetOrbits.period
— Functionperiod(orbit)
Period of an orbit [days].
PlanetOrbits.pmdec
— Functionpmdec(orbit, t)
Get the instantaneous proper motion anomaly [mas/year] in declination of the secondary at the time t
[days].
pmdec(o)
Get the instantaneous proper motion anomaly [mas/year] in declination of the secondary from an instance of AbstractOrbitSolution
.
pmdec(elem, t, M_planet)
Get the instantaneous proper motion anomaly [mas/year] in declination of the primary in at the time t
[days]. The units of M_planet
and elem.M
must match.
pmdec(o, M_planet)
Same as above, but from an orbit solution.
PlanetOrbits.pmra
— Functionpmra(orbit, t)
Get the instantaneous proper motion anomaly [mas/year] in right-ascension of the secondary at the time t
[days].
pmra(o)
Get the instantaneous proper motion anomaly [mas/year] in right-ascension of the secondary from an instance of AbstractOrbitSolution
.
pmra(elem, t, M_planet)
Get the instantaneous proper motion anomaly [mas/year] in right-ascension of the primary in at the time t
[days]. The units of M_planet
and elem.M
must match.
pmra(o, M_planet)
Same as above, but from an orbit solution.
PlanetOrbits.posangle
— Methodposangle(orbit, t)
Calculate the position angle [rad] of the secondary about its primary from our perspective at the time t
[days].
posangle(o)
Calculate the position angle [rad] of the secondary about its primary from our perspective from an instance of AbstractOrbitSolution
.
posangle(elem, t, M_planet)
Calculate the position angle [rad] of the secondary about its primary from our perspective at the time t
[days]. In this case only, the value of M_planet can be arbitrary.
posangle(o, M_planet)
Calculate the position angle [rad] of the primary from our perspective from an instance of AbstractOrbitSolution
. In this case only, the value of M_planet can be arbitrary.
PlanetOrbits.posx
— Functionposx(orbit, t)
Get the offset [AU] from the primary body at the time t
[days].
posx(orbit, t)
Same as above, but from an instance of AbstractOrbitSolution
.
PlanetOrbits.posx
— MethodGet the position in the x direction in astronomical units.
PlanetOrbits.posy
— Functionposy(orbit, t)
Get the offset [AU] from the primary body at the time t
[days].
posy(o)
Same as above, but from an instance of AbstractOrbitSolution
.
PlanetOrbits.posy
— MethodGet the position in the y direction in astronomical units.
PlanetOrbits.posz
— Functionposz(orbit, t)
Get the offset [AU] from the primary body at the time t
[days].
posz(o)
Same as above, but from an instance of AbstractOrbitSolution
.
PlanetOrbits.posz
— MethodGet the position in the z direction in astronomical units.
PlanetOrbits.projectedseparation
— Methodprojectedseparation(orbit, t)
Calculate the projected separation [mas] of the secondary from its primary at the time t
[days].
projectedseparation(o)
Calculate the projected separation [mas] of the secondary from its primary from an instance of AbstractOrbitSolution
.
PlanetOrbits.radvel
— Functionradvel(orbit, t)
Get the radial velocity [m/s] of the secondary along the line of sight at the time t
[days].
radvel(o)
Get the radial velocity [m/s] of the secondary along the line of sight from an instance of AbstractOrbitSolution
.
radvel(elem, t, M_planet)
Get the radial velocity [m/s] of the primary along the line of sight at the time t
[days]. The units of M_planet
and elem.M
must match.
radvel(o, M_planet)
Get the radial velocity [m/s] of the primary along the line of sight from an AbstractOrbitSolution
. The units of M_planet
and elem.M
must match.
PlanetOrbits.raoff
— Functionraoff(orbit, t)
Get the offset [mas] from the primary body in Right Ascension at the time t
[days].
raoff(o)
Get the offset [mas] from the primary body in Right Ascension from an instance of AbstractOrbitSolution
.
PlanetOrbits.semiamplitude
— Functionsemiamplitude(orbit)
Radial velocity semiamplitude [m/s].
PlanetOrbits.semimajoraxis
— Functionsemimajoraxis(orbit)
Semi-major axis of an orbit, if available [au].
PlanetOrbits.totalmass
— Functiontotalmass(orbit)
Total mass of the system in solar masses
PlanetOrbits.trueanom
— Methodtrueanom(orbit, t)
Get the true anomaly [radians] of the secondary at the time t
[days].
trueanom(o)
Get the true anomaly [radians] of the secondary from an instance of AbstractOrbitSolution
.
PlanetOrbits.velx
— MethodGet the velocity in the x direction in astronomical units / year.
PlanetOrbits.vely
— MethodGet the velocity in the y direction in astronomical units / year.
PlanetOrbits.velz
— MethodGet the velocity in the z direction in astronomical units / year.
PlanetOrbits.years2mjd
— Methodyears2mjd()
Convert from fractional years (e.g. 1995.25) into modified julian date.