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).

propertymeaningKepOrbitVisual{KepOrbit}ThieleInnesOrbitRadialVelocityOrbitCartesianOrbitVisual{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).

FunctionKepOrbitVisual{KepOrbit}ThieleInnesOrbitRadialVelocityOrbitCartesianOrbitVisual{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")).

FunctionKepOrbitVisual{KepOrbit}ThieleInnesOrbitRadialVelocityOrbitCartesianOrbitVisual{CartesianOrbit}
meananom
trueanom
eccanom
posx
posy
posz
velx
vely
velz
raoff
decoff
radvel
posangle
pmra
pmdec
accra
accdec

Documentation

PlanetOrbits.AbstractOrbitType
AbstractOrbit

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

source
PlanetOrbits.AbstractOrbitSolutionType
AbstractOrbitSolution

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.

source
PlanetOrbits.AutoType
PlanetOrbits.Auto()

Automatic choice of Kepler solver algorithm. Currently defaults to PlanetOrbits.Markley()

source
PlanetOrbits.GoatType
PlanetOrbits.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.

source
PlanetOrbits.KepOrbitType
KepOrbit(
    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.

source
PlanetOrbits.MarkleyType
PlanetOrbits.Markley()

Kepler solver implementation from AstroLib, based on Markley (1995) Celestial Mechanics and Dynamical Astronomy, 63, 101 (DOI:10.1007/BF00691917).

source
PlanetOrbits.RadialVelocityOrbitType
RadialVelocityOrbit(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.

source
PlanetOrbits.RootsMethodType
PlanetOrbits.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)
source
PlanetOrbits.ThieleInnesOrbitType
ThieleInnesOrbit(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.

source
PlanetOrbits.VisualType
Visual{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.

source
PlanetOrbits.VisualOrbitType
Visual{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.

source
PlanetOrbits.accdecFunction
accdec(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.

source
PlanetOrbits.accraFunction
accra(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.

source
PlanetOrbits.decoffFunction
decoff(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.

source
PlanetOrbits.eccanomMethod
eccanom(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.

source
PlanetOrbits.meananomMethod
meananom(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.

source
PlanetOrbits.mjdMethod
mjd("2020-01-01")

Get the modfied julian day of a date, or in general a UTC timestamp.

source
PlanetOrbits.mjd2dateMethod
mjd2date(modified_julian)

Get a Date value from a modfied julian day, rounded to closest day

Examples

julia> mjd2date(59160.8)
2020-11-08
source
PlanetOrbits.orbitMethod
orbit(...)

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 τ
source
PlanetOrbits.orbitsolveFunction
orbitsolve(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.

source
PlanetOrbits.orbitsolve_νFunction
orbitsolve_ν(elem, ν)

Solve a keplerian orbit from a given true anomaly [rad]. See orbitsolve for the same function accepting a given time.

source
PlanetOrbits.orbitsolve_νFunction
orbitsolve_ν(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.

source
PlanetOrbits.periastronFunction

periastron(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

source
PlanetOrbits.pmdecFunction
pmdec(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.

source
PlanetOrbits.pmraFunction
pmra(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.

source
PlanetOrbits.posangleMethod
posangle(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.

source
PlanetOrbits.posxFunction
posx(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.

source
PlanetOrbits.posyFunction
posy(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.

source
PlanetOrbits.poszFunction
posz(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.

source
PlanetOrbits.projectedseparationMethod
projectedseparation(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.

source
PlanetOrbits.radvelFunction
radvel(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.

source
PlanetOrbits.raoffFunction
raoff(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.

source
PlanetOrbits.trueanomMethod
trueanom(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.

source