KeplerOrbit¶
- class twobody.KeplerOrbit(elements=None, elements_type='kepler', barycenter=None, **kwargs)[source]¶
Bases:
object
Represents a bound Kepler orbit.
- Parameters
- elements
twobody.OrbitalElements
subclass instance Either pass in an
OrbitalElements
object, e.g., an instance oftwobody.KeplerElements
, or pass in the element names themselves. If the latter, anything passed in as kwargs gets passed to the elements class specified byelements_type
. The element names for the defaultelements_type
are included below for convenience.- elements_typestr (optional)
Ignore if you pass in an instantiated
OrbitalElements
object. This argument controls the class that thekwargs
are passed to. The default is'kepler'
, meaning all keyword arguments get passed to thetwobody.KeplerElements
class.- barycenter
twobody.Barycenter
(optional) Parameters that control specification of the barycenter of the orbit.
- elements
Examples
As described above, you can either create an
Elements
object and then pass this toKeplerOrbit
, e.g.,>>> import astropy.units as u >>> from astropy.time import Time >>> from twobody import KeplerElements >>> t0 = Time(2459812.641, format='jd') # reference epoch >>> elem = KeplerElements(a=1.5*u.au, e=0.5, P=1.*u.year, ... omega=67*u.deg, i=21.*u.deg, ... Omega=33*u.deg, M0=53*u.deg, t0=t0) >>> orb = KeplerOrbit(elem)
Or, you can pass in the element names as arguments to the
KeplerOrbit
class:>>> orb = KeplerOrbit(a=1.5*u.au, e=0.5, P=1.*u.year, ... omega=67*u.deg, i=21.*u.deg, Omega=33*u.deg, ... M0=53*u.deg, t0=t0)
Attributes Summary
Methods Summary
icrs
(time)Return the ICRS (i.e.
orbital_plane
(time)Compute the orbit at specified times in the two-body barycentric frame aligned with the orbital plane (xyz).
plot_rv
(time[, ax, rv_unit, t_kwargs, ...])Plot the line-of-sight or radial velocity at the specified times.
radial_velocity
(time[, anomaly_tol, ...])Compute the radial velocity of the body at the specified times relative to the barycenter or reference point, i.e. in the reference plane system not in a solar system barycentric frame.
reference_plane
(time)Compute the orbit at specified times in the two-body barycentric frame aligned with the reference plane coordinate system (XYZ).
unscaled_radial_velocity
(time[, ...])Compute the unscaled radial velocity of the body at the specified times relative to the barycenter or reference point, i.e. in the reference plane system not in a solar system barycentric frame.
Attributes Documentation
- barycenter¶
Methods Documentation
- icrs(time)[source]¶
Return the ICRS (i.e. reference plane) position and velocity of the orbit at the specified time(s).
- Parameters
- timearray_like,
Time
Time array. Either in BMJD or as an Astropy time.
- timearray_like,
- orbital_plane(time)[source]¶
Compute the orbit at specified times in the two-body barycentric frame aligned with the orbital plane (xyz).
- Parameters
- timearray_like,
astropy.time.Time
Array of times as barycentric MJD values, or an Astropy
Time
object containing the times to evaluate at.
- timearray_like,
- plot_rv(time, ax=None, rv_unit=None, t_kwargs=None, plot_kwargs=None)[source]¶
Plot the line-of-sight or radial velocity at the specified times.
- Parameters
- timearray_like,
Time
Time array. Either in BMJD or as an Astropy time.
- ax
Axes
, optional The axis to draw on (default is to grab the current axes using
gca
).- rv_unit
UnitBase
, optional Units to plot the radial velocities in (default is km/s).
- t_kwargsdict, optional
Keyword arguments passed to
astropy.time.Time
with the input time array. For example,dict(format='mjd', scale='tcb')
for Barycentric MJD.- plot_kwargsdict, optional
Any additional arguments or style settings passed to
matplotlib.pyplot.plot()
.
- timearray_like,
- Returns
- ax
Axes
The matplotlib axes object that the RV curve was drawn on.
- ax
- radial_velocity(time, anomaly_tol=None, anomaly_maxiter=None)[source]¶
Compute the radial velocity of the body at the specified times relative to the barycenter or reference point, i.e. in the reference plane system not in a solar system barycentric frame.
This should always be close (in a machine precision sense) to the
z
velocity oforbit.reference_plane(time).
When the barycenter is assumed to be at rest with respect to tangential motion relative to the observer, this should be equivalent to
orbit.icrs(time).radial_velocity
As mentioned above and in Compute and plot a radial velocity curve given orbital elements, the radial velocity computed this way assumes that the barycenter does not move tangentially between epochs and thus ignores spherical projection effects. For sources with large proper motions, the true observable line-of-sight velocity will change slightly over time. We can visualize the expected differences given a full specification of the position and motion of the barycenter:
>>> import astropy.units as u >>> import astropy.coordinates as coord >>> from astropy.time import Time >>> from twobody import Barycenter, KeplerOrbit >>> origin = coord.ICRS(ra=170.8743*u.deg, dec=-71.34*u.deg, ... distance=57.134*u.pc, ... pm_ra_cosdec=-206.718*u.mas/u.yr, ... pm_dec=301.82*u.mas/u.yr, ... radial_velocity=41.84*u.km/u.s) >>> baryc = Barycenter(origin=origin, t0=Time('J2000')) >>> orb = KeplerOrbit(P=1.5*u.year, e=0.67, a=1.77*u.au, ... omega=17.14*u.deg, i=65*u.deg, Omega=0*u.deg, ... M0=35.824*u.deg, t0=Time('J2015.0'), ... barycenter=baryc) >>> t = Time('J2000') + np.linspace(0, 15, 10000) * u.year >>> true_rv = orb.icrs(t).radial_velocity >>> approx_rv = orb.radial_velocity(t) >>> (true_rv - approx_rv).to(u.m/u.s).max() <Quantity 3.315196290913036 m / s>
In this case, the maximum difference is only ~3 m/s.
- Parameters
- timearray_like,
astropy.time.Time
Array of times as barycentric MJD values, or an Astropy
Time
object containing the times to evaluate at.- anomaly_tolnumeric (optional)
Tolerance passed to
eccentric_anomaly_from_mean_anomaly
for solving for the eccentric anomaly. See default value in that function.- anomaly_maxiternumeric (optional)
Maximum number of iterations to use in
eccentric_anomaly_from_mean_anomaly
for solving for the eccentric anomaly. See default value in that function.
- timearray_like,
- reference_plane(time)[source]¶
Compute the orbit at specified times in the two-body barycentric frame aligned with the reference plane coordinate system (XYZ).
- Parameters
- timearray_like,
astropy.time.Time
Array of times as barycentric MJD values, or an Astropy
Time
object containing the times to evaluate at.
- timearray_like,
- unscaled_radial_velocity(time, anomaly_tol=None, anomaly_maxiter=None)[source]¶
Compute the unscaled radial velocity of the body at the specified times relative to the barycenter or reference point, i.e. in the reference plane system not in a solar system barycentric frame.
See the docstring of
radial_velocity
for more information and caveats.- Parameters
- timearray_like,
astropy.time.Time
Array of times as barycentric MJD values, or an Astropy
Time
object containing the times to evaluate at.- anomaly_tolnumeric (optional)
Tolerance passed to
eccentric_anomaly_from_mean_anomaly
for solving for the eccentric anomaly. See default value in that function.- anomaly_maxiternumeric (optional)
Maximum number of iterations to use in
eccentric_anomaly_from_mean_anomaly
for solving for the eccentric anomaly. See default value in that function.
- timearray_like,
- Returns
- rvnumeric [m/s]
Relative radial velocity - does not include systemtic velocity!