Last active
September 30, 2021 22:10
-
-
Save bamford/a55afccbb68d0ed12b1dc5ebc5af37cd to your computer and use it in GitHub Desktop.
Eclipsing Binaries
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
#!/bin/python | |
from astropy import time, coordinates as coord, units as u, table as tab | |
from astropy.table import Table | |
from astropy.time import Time, TimeDelta | |
import numpy as np | |
path = 'eclipses' | |
t = Table.read('eclipsing_binaries.fits') | |
t.write('{}/binaries.html'.format(path), overwrite=True, | |
htmldict={'cssfiles': ['./eclipseTable.css'], | |
'table_class': 'eclipseTable'}) | |
# Now calculate times of next minima, and select objects suitable for observation... | |
t['coord'] = coord.SkyCoord(t['ra'], t['dec'], unit=('hour', 'deg')) | |
# compute next minima | |
Mnow = Time.now().jd | |
t['epoch'] = (np.ceil((Mnow - t['M0']) / t['P'])).astype(np.int) | |
t['Mnext'] = t['M0'] + t['epoch'] * t['P'] | |
# compute local utc of next minima (correct from heliocentric time) | |
nottingham = coord.EarthLocation(lon=-1.19201, lat=+52.94166, height=35.0) | |
times_hjd = time.Time(t['Mnext'], format='jd', scale='utc', location=nottingham) | |
ltt_helio = times_hjd.light_travel_time(t['coord'], 'heliocentric') | |
times_jd = times_hjd - ltt_helio | |
t['Mnext_UTC'] = times_jd.iso | |
# This is a worse-case-scenario error, reality may not be so bad: | |
t['Mnext_err_min'] = ((t['M0_err'] + t['epoch'] * t['P_err']) * 24 * 60).round(1) | |
# compute sky position | |
altaz = t['coord'].transform_to(coord.AltAz(obstime=times_jd, location=nottingham)) | |
t['alt'], t['az'] = altaz.alt.round(1), altaz.az.round(1) | |
# Get all EA/B with primary eclipse between 5pm and 10pm UTC tonight with duration less than 3 hours, altitude above 30 deg and brighter than V=12... | |
Tearly = Time(Time.now().value.replace(hour=17, minute=0, second = 0)) | |
Tlate = Time(Time.now().value.replace(hour=22, minute=0, second = 0)) | |
sel = t['minimum_type'] != 'SEC' | |
sel &= t['Mnext'] > Tearly.jd | |
sel &= t['Mnext'] < Tlate.jd | |
sel &= [('EA' in x)|('EB' in x) for x in t['var_type']] | |
sel &= (t['D_day'] > 0.01) & (t['D_day'] < 3/24.) | |
sel &= t['alt'] > 30 | |
sel &= t['V'] < 12 | |
useful_columns = ['constellation', 'id', 'var_type', 'V', 'Mnext_UTC', | |
'Mnext_err_min', 'P', 'D_min', 'alt', 'az', 'ra', 'dec', 'url'] | |
tout = t[sel][useful_columns] | |
# Use the url to see the O-C diagram. A positive O-C means the minimum is later than predicted. Use the left axis to find the O-C in days. 1 hour = 0.04 days. | |
fn = '{}/tonight'.format(path) | |
tout.write(fn+'.fits', overwrite=True) | |
tout.write(fn+'.html', overwrite=True, | |
htmldict={'cssfiles': ['./eclipseTable.css'], | |
'table_class': 'eclipseTable'}) |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
eclipsing_binaries.py
generates the nightly table of eclipses here.