Quick snippet to use astropy to convert from heliocentric times to barycentric and vice-versa

helio_bary.py

from astropy.coordinates import SkyCoord, EarthLocation
from astropy import units as u
from astropy.time import Time


def helio_to_bary(coords, hjd, obs_name):
    helio = Time(hjd, scale='utc', format='jd')
    obs = EarthLocation.of_site(obs_name)
    star = SkyCoord(coords, unit=(u.hour, u.deg)) 
    ltt = helio.light_travel_time(star, 'heliocentric', location=obs)
    guess = helio - ltt
    # if we assume guess is correct - how far is heliocentric time away from true value?
    delta = (guess + guess.light_travel_time(star, 'heliocentric', obs)).jd  - helio.jd
    # apply this correction
    guess -= delta * u.d

    ltt = guess.light_travel_time(star, 'barycentric', obs)
    return guess.tdb + ltt


def bary_to_helio(coords, bjd, obs_name):
    bary = Time(bjd, scale='tdb', format='jd')
    obs = EarthLocation.of_site(obs_name)
    star = SkyCoord(coords, unit=(u.hour, u.deg))
    ltt = bary.light_travel_time(star, 'barycentric', location=obs) 
    guess = bary - ltt
    delta = (guess + guess.light_travel_time(star, 'barycentric', obs)).jd  - bary.jd
    guess -= delta * u.d
    
    ltt = guess.light_travel_time(star, 'heliocentric', obs)
    return guess.utc + ltt

Hi Stuart,

Good to meet you yesterday and I hope you don't mind me asking a question based on your helio-bary-helio code.

In the Scheduler we've assumed that t0 (the time anchor point provided by observers in phase 2 for creating the observation windows) will come to as in TDB format. We then convert this to UTC as our other datetimes are all in this format (rise/set times, dawn/dusk, etc).

Astropy has a method for converting UTC to TDB, but not vice versa as far as I can see. So I've done the conversion as if the formats are the other way round (t0 comes as UTC and needs to be converted to TDB), and then used the answer to actually convert TDB to UTC. This isn't axact, but it works quite well when I test the answer against the following website: http://astroutils.astronomy.ohio-state.edu/time/UTC2bjd.html.

def convert_tdb_to_utc(dt, targ):
    """
    Convert TDB to UTC
    'targ' is a SkyCoord object of the target star
    """

    t = Time(dt, scale='utc', location=(str(INT_Lon) + 'd', str(INT_Lat) + 'd'))

    td = t.light_travel_time(targ)

    # deduct the offset from the original UTC time to get the new UTC time. This has the effect
    #   of treating the original times as barycentric time and reverse engineering the UTC time.
    dt = dt - td.datetime

    return dt

So, my question is, does this approach sound reasonable, or is there a much easier, more accurate (but quick) way of converting TDB to UTC?

Best regards

Ian

Hi Ian,

What you've done is slightly incorrect. There are two issues; one is that the observer-provided t0 is on the TDB timescale, and you should initialise as such:

t = Time(dt, scale='tdb',  location=(str(INT_Lon) + 'd', str(INT_Lat) + 'd'))

The second is that the call to light_travel_time will not quite get the correct time, as it is meant to be used the other way around. This is what the second iteration in my functions is meant to correct for; guess the light_travel_time, use the guess to correct to UTC, and see if the corrected UTC converts to the given TDB. A version of your code using these fixes would look like:

def convert_tdb_to_utc(dt, targ):
    """
    Convert TDB to UTC
    'targ' is a SkyCoord object of the target star
    """
    t = Time(dt, scale='tdb',  location=(str(INT_Lon) + 'd', str(INT_Lat) + 'd'))
    td = t.light_travel_time(targ)
    guess = t - td
    # if guess is correct, what's the difference between the inferred BJD and the given BJD?
    delta = (guess + guess.light_travel_time(targ)).jd  - t.jd
    # apply this correction
    guess -= delta * u.d
    # return on UTC scale as date time
    return guess.utc.datetime

However, the correction involves a second call to light_travel_time and only provides a minor change in precision - around 1/100th of a second. So a cruder but faster version (x2) would just be

def convert_tdb_to_utc(dt, targ):
    """
    Convert TDB to UTC
    'targ' is a SkyCoord object of the target star
    """
    t = Time(dt, scale='tdb',  location=(str(INT_Lon) + 'd', str(INT_Lat) + 'd'))
    td = t.light_travel_time(targ)
    correct = t - td
    return correct.utc.datetime

Thanks Stuart, I like the idea of checking the guess and correcting it slightly. We don’t have to run this process too many times in the code, only in setup in fact, so I might put this in for completeness. Thanks for your help, it’s much appreciated. Ian