Barycentric velocity correction¶
The baryvel()
and baryCorr()
functions allow to
calculate Earth’s helio- and barycentric motion and project it
onto a given direction toward a star. The helcorr()
includes a barycentric correction including the effect caused
by the rotating Earth.
PyA’s baryvel function is a port of its IDL Astrolib’s counterpart. The helcorr function has been ported from the REDUCE package (see Piskunov & Valenti 2002, A&A 385, 1095).
- PyAstronomy.pyasl.baryvel(dje, deq)¶
Calculate helio- and barycentric velocity.
Note
The “JPL” option present in IDL is not provided here.
- Parameters
- djefloat
Julian ephemeris date
- deqfloat
Epoch of mean equinox of helio- and barycentric velocity output. If deq is zero, deq is assumed to be equal to dje.
- Returns
- dvelharray
Heliocentric velocity vector [km/s].
- dvelbarray
Barycentric velocity vector [km/s].
Notes
Note
This function was ported from the IDL Astronomy User’s Library.
- IDL - Documentation
pro baryvel, dje, deq, dvelh, dvelb, JPL = JPL
- NAME:
BARYVEL
- PURPOSE:
Calculates heliocentric and barycentric velocity components of Earth.
- EXPLANATION:
BARYVEL takes into account the Earth-Moon motion, and is useful for radial velocity work to an accuracy of ~1 m/s.
- CALLING SEQUENCE:
BARYVEL, dje, deq, dvelh, dvelb, [ JPL = ]
- INPUTS:
DJE - (scalar) Julian ephemeris date. DEQ - (scalar) epoch of mean equinox of dvelh and dvelb. If deq=0
then deq is assumed to be equal to dje.
- OUTPUTS:
DVELH: (vector(3)) heliocentric velocity component. in km/s DVELB: (vector(3)) barycentric velocity component. in km/s
The 3-vectors DVELH and DVELB are given in a right-handed coordinate system with the +X axis toward the Vernal Equinox, and +Z axis toward the celestial pole.
- OPTIONAL KEYWORD SET:
- JPL - if /JPL set, then BARYVEL will call the procedure JPLEPHINTERP
to compute the Earth velocity using the full JPL ephemeris. The JPL ephemeris FITS file JPLEPH.405 must exist in either the current directory, or in the directory specified by the environment variable ASTRO_DATA. Alternatively, the JPL keyword can be set to the full path and name of the ephemeris file. A copy of the JPL ephemeris FITS file is available in
- PROCEDURES CALLED:
Function PREMAT() – computes precession matrix JPLEPHREAD, JPLEPHINTERP, TDB2TDT - if /JPL keyword is set
- NOTES:
Algorithm taken from FORTRAN program of Stumpff (1980, A&A Suppl, 41,1) Stumpf claimed an accuracy of 42 cm/s for the velocity. A comparison with the JPL FORTRAN planetary ephemeris program PLEPH found agreement to within about 65 cm/s between 1986 and 1994
If /JPL is set (using JPLEPH.405 ephemeris file) then velocities are given in the ICRS system; otherwise in the FK4 system.
- EXAMPLE:
- Compute the radial velocity of the Earth toward Altair on 15-Feb-1994
using both the original Stumpf algorithm and the JPL ephemeris
IDL> jdcnv, 1994, 2, 15, 0, jd ;==> JD = 2449398.5 IDL> baryvel, jd, 2000, vh, vb ;Original algorithm
==> vh = [-17.07243, -22.81121, -9.889315] ;Heliocentric km/s ==> vb = [-17.08083, -22.80471, -9.886582] ;Barycentric km/s
- IDL> baryvel, jd, 2000, vh, vb, /jpl ;JPL ephemeris
==> vh = [-17.07236, -22.81126, -9.889419] ;Heliocentric km/s ==> vb = [-17.08083, -22.80484, -9.886409] ;Barycentric km/s
IDL> ra = ten(19,50,46.77)*15/!RADEG ;RA in radians IDL> dec = ten(08,52,3.5)/!RADEG ;Dec in radians IDL> v = vb[0]*cos(dec)*cos(ra) + $ ;Project velocity toward star
vb[1]*cos(dec)*sin(ra) + vb[2]*sin(dec)
- REVISION HISTORY:
Jeff Valenti, U.C. Berkeley Translated BARVEL.FOR to IDL. W. Landsman, Cleaned up program sent by Chris McCarthy (SfSU) June 1994 Converted to IDL V5.0 W. Landsman September 1997 Added /JPL keyword W. Landsman July 2001 Documentation update W. Landsman Dec 2005
- PyAstronomy.pyasl.baryCorr(jd, ra, dec, deq=0.0)¶
Calculate barycentric correction.
This function uses the
baryvel()
function to calculate the helio- and barycentric motion of the Earth and projects it onto the direction to the star.Note
Positive return values indicate that the Earth moves toward the star.
- Parameters
- jdfloat
The time at which to calculate the correction.
- rafloat
Right ascension in degrees.
- decfloat
Declination in degrees.
- deqfloat, optional
The mean equinox of barycentric velocity calculation (see
bryvel()
). If zero, it is assumed to be the same as jd.
- Returns
- Projected heliocentric velocityfloat
Heliocentric velocity toward star [km/s]
- Projected barycentric velocityfloat
Barycentric velocity toward star [km/s]
- PyAstronomy.pyasl.helcorr(obs_long, obs_lat, obs_alt, ra2000, dec2000, jd, debug=False)¶
Calculate barycentric velocity correction.
This function calculates the motion of an observer in the direction of a star. In contrast to
baryvel()
andbaryCorr()
, the rotation of the Earth is taken into account.The coordinates (ra2000, dec2000) are precessed to the epoch defined by jd. These coordinates are used in the calculation.
Note
This function was ported from the REDUCE IDL package. See Piskunov & Valenti 2002, A&A 385, 1095 for a detailed description of the package and/or visit http://www.astro.uu.se/~piskunov/RESEARCH/REDUCE/
Warning
Contrary to the original implementation the longitude increases toward the East and the right ascension is given in degrees instead of hours. The JD is given as is, in particular, nothing needs to be subtracted.
- Parameters
- obs_longfloat
Longitude of observatory (degrees, eastern direction is positive)
- obs_latfloat
Latitude of observatory [deg]
- obs_altfloat
Altitude of observatory [m]
- ra2000float
Right ascension of object for epoch 2000.0 [deg]
- dec2000float
Declination of object for epoch 2000.0 [deg]
- jdfloat
Julian date for the middle of exposure.
- Returns
- Barycentric correctionfloat
The barycentric correction accounting for the rotation of the Earth, the rotation of the Earth’s center around the Earth-Moon barycenter, and the motion of the Earth-Moon barycenter around the center of the Sun [km/s].
- HJDfloat
Heliocentric Julian date for middle of exposure.
Notes
- IDL REDUCE - Documentation
Calculates heliocentric Julian date, barycentric and heliocentric radial velocity corrections from:
INPUT: <OBSLON> Longitude of observatory (degrees, western direction is positive) <OBSLAT> Latitude of observatory (degrees) <OBSALT> Altitude of observatory (meters) <RA2000> Right ascension of object for epoch 2000.0 (hours) <DE2000> Declination of object for epoch 2000.0 (degrees) <JD> Julian date for the middle of exposure [DEBUG=] set keyword to get additional results for debugging
OUTPUT: <CORRECTION> barycentric correction - correction for rotation of earth,
rotation of earth center about the earth-moon barycenter, earth-moon barycenter about the center of the Sun.
<HJD> Heliocentric Julian date for middle of exposure
Algorithms used are taken from the IRAF task noao.astutils.rvcorrect and some procedures of the IDL Astrolib are used as well. Accuracy is about 0.5 seconds in time and about 1 m/s in velocity.
History: written by Peter Mittermayer, Nov 8,2003 2005-January-13 Kudryavtsev Made more accurate calculation of the sidereal time.
Conformity with MIDAS compute/barycorr is checked.
2005-June-20 Kochukhov Included precession of RA2000 and DEC2000 to current epoch
Example: Applying baryvel and baryCorr¶
from __future__ import print_function, division
from PyAstronomy import pyasl
jd = 2.476468576e6
heli, bary = pyasl.baryvel(jd, deq=2000.0)
print("Earth's velocity at JD: ", jd)
print("Heliocentric velocity [km/s]: ", heli)
print("Barycentric velocity [km/s] : ", bary)
# Coordinates of Sirius
ra = 101.28715535
dec = -16.71611587
vh, vb = pyasl.baryCorr(jd, ra, dec, deq=2000.0)
print("Barycentric velocity of Earth toward Sirius: ", vb)
Example: Obtaining a barycentric correction¶
from __future__ import print_function, division
from PyAstronomy import pyasl
# Coordinates of European Southern Observatory
# (Coordinates of UT1)
longitude = 289.5967661
latitude = -24.62586583
altitude = 2635.43
# Coordinates of HD 12345 (J2000)
ra2000 = 030.20313477
dec2000 = -12.87498346
# (Mid-)Time of observation
jd = 2450528.2335
# Calculate barycentric correction (debug=True show
# various intermediate results)
corr, hjd = pyasl.helcorr(longitude, latitude, altitude, \
ra2000, dec2000, jd, debug=True)
print("Barycentric correction [km/s]: ", corr)
print("Heliocentric Julian day: ", hjd)