Eclipsing Binaries

eclipsing_binaries.py

#!/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'})

EclipsingBinaries.ipynb

EclipsingBinariesExample.ipynb

{
  "cells": [
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "This notebook uses a [combined version](https://www.dropbox.com/s/jrmhnk6vgezq63f/eclipsing_binaries.fits?dl=0) of the tables prepared by J.M. Kreiner from the [Mt Suhora Astronomical Observatory](http://www.as.up.krakow.pl/ephem/). It provides an example of calculating times of next minima and selecting suitable targets for observation.  Note that these are linear ephemerides: O-C offsets need to be considered."
    },
    {
      "metadata": {
        "collapsed": true,
        "trusted": true
      },
      "cell_type": "code",
      "source": "from astropy import time, coordinates as coord, units as u, table as tab\nfrom astropy.table import Table\nfrom astropy.time import Time, TimeDelta",
      "execution_count": 8,
      "outputs": []
    },
    {
      "metadata": {
        "trusted": true
      },
      "cell_type": "code",
      "source": "t = Table.read('eclipsing_binaries.fits')",
      "execution_count": 9,
      "outputs": []
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "Calculate times of next minima, and select objects suitable for observation..."
    },
    {
      "metadata": {
        "collapsed": true,
        "trusted": true
      },
      "cell_type": "code",
      "source": "Mnow = Time.now().jd\nt['epoch'] = (np.ceil((Mnow - t['M0']) / t['P'])).astype(np.int)\nt['Mnext'] = t['M0'] + t['epoch'] * t['P']",
      "execution_count": 10,
      "outputs": []
    },
    {
      "metadata": {
        "collapsed": true,
        "hide_input": false,
        "trusted": true
      },
      "cell_type": "code",
      "source": "t['coord'] = coord.SkyCoord(t['ra'], t['dec'], unit=('hour', 'deg'))",
      "execution_count": 11,
      "outputs": []
    },
    {
      "metadata": {
        "trusted": true
      },
      "cell_type": "code",
      "source": "# compute local utc of next minima (correct from heliocentric time)\nnottingham = coord.EarthLocation(lon=-1.19201, lat=+52.94166, height=35.0)\ntimes_hjd = time.Time(t['Mnext'], format='jd', scale='utc', location=nottingham)\nltt_helio = times_hjd.light_travel_time(t['coord'], 'heliocentric')\ntimes_jd = times_hjd - ltt_helio\nt['Mnext_UTC'] = times_jd.iso\n# This is a worse-case-scenario error, reality may not be so bad:\nt['Mnext_err_min'] = ((t['M0_err'] + t['epoch'] * t['P_err']) * 24 * 60).round(1)",
      "execution_count": 21,
      "outputs": [
        {
          "name": "stderr",
          "output_type": "stream",
          "text": "WARNING: Tried to get polar motions for times after IERS data is valid. Defaulting to polar motion from the 50-yr mean for those. This may affect precision at the 10s of arcsec level [astropy.coordinates.builtin_frames.utils]\n"
        }
      ]
    },
    {
      "metadata": {
        "trusted": true
      },
      "cell_type": "code",
      "source": "# compute sky position\naltaz = t['coord'].transform_to(coord.AltAz(obstime=times_jd, location=nottingham))\nt['alt'], t['az'] = altaz.alt.round(1), altaz.az.round(1)",
      "execution_count": 22,
      "outputs": [
        {
          "name": "stderr",
          "output_type": "stream",
          "text": "WARNING: Tried to get polar motions for times after IERS data is valid. Defaulting to polar motion from the 50-yr mean for those. This may affect precision at the 10s of arcsec level [astropy.coordinates.builtin_frames.utils]\n"
        }
      ]
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "Get all EA/B with primary eclipse between 6pm and 10pm UTC (7pm-11pm BST) on 2017-10-27 with duration less than 3 hours, altitude above 30 deg and brighter than V=12..."
    },
    {
      "metadata": {
        "trusted": true
      },
      "cell_type": "code",
      "source": "Tearly = Time('2017-10-27T18:00:00')\nTlate = Time('2017-10-27T22:00:00')\nsel = t['minimum_type'] != 'SEC'\nsel &= t['Mnext'] > Tearly.jd\nsel &= t['Mnext'] < Tlate.jd\nsel &= [('EA' in x)|('EB' in x) for x in t['var_type']]\nsel &= (t['D_day'] > 0.01) & (t['D_day'] < 3/24.)\nsel &= t['alt'] > 30 * u.deg\nsel &= t['V'] < 12\nuseful_columns = ['constellation', 'id', 'var_type', 'V', 'Mnext_UTC',\n                  'Mnext_err_min', 'P', 'D_min', 'alt', 'az', 'ra', 'dec', 'url']\nt[sel][useful_columns]",
      "execution_count": 23,
      "outputs": [
        {
          "data": {
            "text/html": "&lt;Table length=5&gt;\n<table id=\"table4605045504\" class=\"table-striped table-bordered table-condensed\">\n<thead><tr><th>constellation</th><th>id</th><th>var_type</th><th>V</th><th>Mnext_UTC</th><th>Mnext_err_min</th><th>P</th><th>D_min</th><th>alt</th><th>az</th><th>ra</th><th>dec</th><th>url</th></tr></thead>\n<thead><tr><th></th><th></th><th></th><th></th><th></th><th></th><th></th><th></th><th>deg</th><th>deg</th><th></th><th></th><th></th></tr></thead>\n<thead><tr><th>str3</th><th>str5</th><th>str10</th><th>float32</th><th>str23</th><th>float64</th><th>float64</th><th>float64</th><th>float64</th><th>float64</th><th>str8</th><th>str9</th><th>str48</th></tr></thead>\n<tr><td>CAS</td><td>CC</td><td>EB/DM</td><td>7.06</td><td>2017-10-27 20:50:00.598</td><td>45.4</td><td>3.366308</td><td>43.2</td><td>56.7</td><td>53.7</td><td>03:14:05</td><td>59:33:49</td><td>http://www.as.up.krakow.pl/minicalc/CASCC.HTM</td></tr>\n<tr><td>CEP</td><td>WW</td><td>EA/SD:</td><td>11.1</td><td>2017-10-27 21:19:08.369</td><td>4.9</td><td>4.600848</td><td>72.0</td><td>70.6</td><td>339.0</td><td>22:18:27</td><td>69:51:40</td><td>http://www.as.up.krakow.pl/minicalc/CEPWW.HTM</td></tr>\n<tr><td>CEP</td><td>XZ</td><td>EB/DM:</td><td>8.0</td><td>2017-10-27 20:29:27.118</td><td>4.4</td><td>5.09725</td><td>72.0</td><td>75.5</td><td>353.3</td><td>22:32:25</td><td>67:09:03</td><td>http://www.as.up.krakow.pl/minicalc/CEPXZ.HTM</td></tr>\n<tr><td>CEP</td><td>AI</td><td>EB/DM</td><td>9.18</td><td>2017-10-27 18:41:04.954</td><td>17.1</td><td>4.225334</td><td>57.6</td><td>82.3</td><td>53.5</td><td>21:46:22</td><td>56:55:02</td><td>http://www.as.up.krakow.pl/minicalc/CEPAI.HTM</td></tr>\n<tr><td>CYG</td><td>KV</td><td>EB/SD</td><td>11.5</td><td>2017-10-27 18:39:01.694</td><td>4.0</td><td>2.8390041</td><td>43.2</td><td>72.2</td><td>209.4</td><td>20:15:38</td><td>36:47:37</td><td>http://www.as.up.krakow.pl/minicalc/CYGKV.HTM</td></tr>\n</table>",
            "text/plain": "<Table length=5>\nconstellation  id  ...   dec                         url                     \n                   ...                                                       \n     str3     str5 ...   str9                       str48                    \n------------- ---- ... -------- ---------------------------------------------\n          CAS   CC ... 59:33:49 http://www.as.up.krakow.pl/minicalc/CASCC.HTM\n          CEP   WW ... 69:51:40 http://www.as.up.krakow.pl/minicalc/CEPWW.HTM\n          CEP   XZ ... 67:09:03 http://www.as.up.krakow.pl/minicalc/CEPXZ.HTM\n          CEP   AI ... 56:55:02 http://www.as.up.krakow.pl/minicalc/CEPAI.HTM\n          CYG   KV ... 36:47:37 http://www.as.up.krakow.pl/minicalc/CYGKV.HTM"
          },
          "execution_count": 23,
          "metadata": {},
          "output_type": "execute_result"
        }
      ]
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "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."
    }
  ],
  "metadata": {
    "_draft": {
      "nbviewer_url": "https://gist.github.com/a55afccbb68d0ed12b1dc5ebc5af37cd"
    },
    "gist": {
      "id": "a55afccbb68d0ed12b1dc5ebc5af37cd",
      "data": {
        "description": "Eclipsing Binaries ",
        "public": true
      }
    },
    "hide_input": false,
    "kernelspec": {
      "name": "python3",
      "display_name": "Python [default]",
      "language": "python"
    },
    "language_info": {
      "name": "python",
      "version": "3.6.1",
      "mimetype": "text/x-python",
      "codemirror_mode": {
        "name": "ipython",
        "version": 3
      },
      "pygments_lexer": "ipython3",
      "nbconvert_exporter": "python",
      "file_extension": ".py"
    },
    "latex_envs": {
      "eqNumInitial": 1,
      "eqLabelWithNumbers": true,
      "current_citInitial": 1,
      "cite_by": "apalike",
      "bibliofile": "biblio.bib",
      "LaTeX_envs_menu_present": true,
      "labels_anchors": false,
      "latex_user_defs": false,
      "user_envs_cfg": false,
      "report_style_numbering": false,
      "autocomplete": true,
      "hotkeys": {
        "equation": "Ctrl-E",
        "itemize": "Ctrl-I"
      }
    }
  },
  "nbformat": 4,
  "nbformat_minor": 1
}

eclipsing_binaries.py generates the nightly table of eclipses here.