KIC08462852-single.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Photometry using astropy photutils and AAVSO comparison stars\n",
    "\n",
    "## 1 Import dependencies and setup matplotlib"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import requests, math, glob\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "from photutils import DAOStarFinder\n",
    "from astropy.stats import mad_std\n",
    "from astropy.io import fits\n",
    "from astropy.coordinates import SkyCoord\n",
    "from astropy.wcs import WCS\n",
    "import astropy.units as u\n",
    "import matplotlib.pyplot as plt\n",
    "from photutils import aperture_photometry, CircularAperture\n",
    "from astroquery.simbad import Simbad\n",
    "from datetime import datetime\n",
    "import warnings\n",
    "warnings.filterwarnings('ignore')\n",
    "\n",
    "%matplotlib inline\n",
    "plt.style.use('seaborn')\n",
    "pd.options.mode.chained_assignment = None "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2  Define input file, star name and comparason star range¶"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "FITS_FILE = '/home/dokeeffe/Pictures/CalibratedLight/2017-10-02/KIC08462852-2017-10-02-20-56-22Light_005.fits'\n",
    "STAR_NAME = 'KIC08462852'\n",
    "BRIGHTEST_COMPARISON_STAR_MAG = 11.0\n",
    "DIMMEST_COMPARISON_STAR_MAG = 13.0"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "## 3 Define a function to download comparison stars from AAVSO. \n",
    "\n",
    "This will return a tuple containing an array and the chart ID."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def get_comp_stars(ra,dec,filter_band='V',field_of_view=18.5):\n",
    "    result = []\n",
    "    vsp_template = 'https://www.aavso.org/apps/vsp/api/chart/?format=json&fov={}&maglimit=18.5&ra={}&dec={}'\n",
    "    print(vsp_template.format(field_of_view, ra, dec))\n",
    "    r = requests.get(vsp_template.format(field_of_view, ra, dec))\n",
    "    chart_id = r.json()['chartid']\n",
    "    print('Downloaded Comparison Star Chart ID {}'.format(chart_id))\n",
    "    for star in r.json()['photometry']:\n",
    "        comparison = {}\n",
    "        comparison['auid'] = star['auid']\n",
    "        comparison['ra'] = star['ra']\n",
    "        comparison['dec'] = star['dec']\n",
    "        for band in star['bands']:\n",
    "            if band['band'] == filter_band:\n",
    "                comparison['vmag'] = band['mag']\n",
    "                comparison['error'] = band['error']\n",
    "        result.append(comparison)\n",
    "    return result, chart_id"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4 Download comparison stars and search simbad for our target.\n",
    "\n",
    "Use astroquery to locate the RA/DEC of our target star.\n",
    "\n",
    "Here we also define `results` which will be a collection of data progressivly enriched as we go\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "https://www.aavso.org/apps/vsp/api/chart/?format=json&fov=18.5&maglimit=18.5&ra=20 06 15.4553&dec=44 27 24.793\n",
      "Downloaded Comparison Star Chart ID X21328CWX\n",
      "5 comp stars found\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "[{'auid': '000-BLS-551',\n",
       "  'dec': '44:22:48.0',\n",
       "  'error': 0.054,\n",
       "  'ra': '20:06:48.09',\n",
       "  'vmag': 11.263},\n",
       " {'auid': '000-BML-045',\n",
       "  'dec': '44:27:05.0',\n",
       "  'error': 0.036,\n",
       "  'ra': '20:06:36.68',\n",
       "  'vmag': 12.415},\n",
       " {'auid': '000-BLS-555',\n",
       "  'dec': '44:30:54.0',\n",
       "  'error': 0.05,\n",
       "  'ra': '20:06:21.25',\n",
       "  'vmag': 12.789},\n",
       " {'auid': '000-BML-046',\n",
       "  'dec': '44:35:41.8',\n",
       "  'error': 0.048,\n",
       "  'ra': '20:06:53.14',\n",
       "  'vmag': 12.983},\n",
       " {'auid': '000-BML-047',\n",
       "  'dec': '44:25:53.7',\n",
       "  'error': 0.038,\n",
       "  'ra': '20:06:00.39',\n",
       "  'vmag': 13.327},\n",
       " {'auid': 'target', 'dec': '44 27 24.793', 'ra': '20 06 15.4553'}]"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "astroquery_results = Simbad.query_object(STAR_NAME)\n",
    "TARGET_RA = str(astroquery_results[0]['RA'])\n",
    "TARGET_DEC = str(astroquery_results[0]['DEC']).replace('+','').replace('-','')\n",
    "results, chart_id = get_comp_stars(TARGET_RA, TARGET_DEC)\n",
    "print('{} comp stars found'.format(len(results)))\n",
    "results.append({'auid': 'target', 'ra': TARGET_RA, 'dec': TARGET_DEC})\n",
    "\n",
    "results"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 5 Extract all sources from image"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<Table length=115>\n",
       "<table id=\"table140277282068856\" class=\"table-striped table-bordered table-condensed\">\n",
       "<thead><tr><th>id</th><th>xcentroid</th><th>ycentroid</th><th>sharpness</th><th>roundness1</th><th>roundness2</th><th>npix</th><th>sky</th><th>peak</th><th>flux</th><th>mag</th></tr></thead>\n",
       "<thead><tr><th>int64</th><th>float64</th><th>float64</th><th>float64</th><th>float64</th><th>float64</th><th>float64</th><th>float64</th><th>float64</th><th>float64</th><th>float64</th></tr></thead>\n",
       "<tr><td>1</td><td>645.55658029</td><td>7.69908327045</td><td>0.45424855519</td><td>0.174645032246</td><td>0.355539480239</td><td>25.0</td><td>0.0</td><td>70365.2937781</td><td>28.269908551</td><td>-3.62831100903</td></tr>\n",
       "<tr><td>2</td><td>1309.16598231</td><td>11.5995843996</td><td>0.44672872522</td><td>0.253701614478</td><td>0.497475903102</td><td>25.0</td><td>0.0</td><td>7287.61409036</td><td>2.85692255025</td><td>-1.13974616759</td></tr>\n",
       "<tr><td>3</td><td>1375.65792899</td><td>11.9995956879</td><td>0.496571246411</td><td>0.185979847486</td><td>0.525493037437</td><td>25.0</td><td>0.0</td><td>3368.49664403</td><td>1.05632128534</td><td>-0.0594900777583</td></tr>\n",
       "<tr><td>4</td><td>1122.94782151</td><td>16.7158994909</td><td>0.460102848306</td><td>-0.0896697826626</td><td>0.47276708619</td><td>25.0</td><td>0.0</td><td>8200.64521946</td><td>3.19746747319</td><td>-1.26201533809</td></tr>\n",
       "<tr><td>5</td><td>1362.4645143</td><td>26.9619775656</td><td>0.433803791444</td><td>0.112942501751</td><td>0.429731153435</td><td>25.0</td><td>0.0</td><td>21372.8993241</td><td>8.4265184611</td><td>-2.31412044147</td></tr>\n",
       "<tr><td>6</td><td>611.037214433</td><td>28.7876045811</td><td>0.401685842506</td><td>0.167921662456</td><td>0.419344081499</td><td>25.0</td><td>0.0</td><td>5478.05140726</td><td>1.98099565096</td><td>-0.742208805242</td></tr>\n",
       "<tr><td>7</td><td>358.189493857</td><td>34.4617814911</td><td>0.503866132614</td><td>-0.00948179156524</td><td>0.26618329465</td><td>25.0</td><td>0.0</td><td>90850.6289014</td><td>34.2481732847</td><td>-3.83659353055</td></tr>\n",
       "<tr><td>8</td><td>449.539598592</td><td>82.1069177762</td><td>0.379309325403</td><td>0.0980523701524</td><td>0.294225743108</td><td>25.0</td><td>0.0</td><td>7390.86311537</td><td>2.73814589578</td><td>-1.09364146198</td></tr>\n",
       "<tr><td>9</td><td>964.423830549</td><td>83.9894147141</td><td>0.385366389308</td><td>0.120651000841</td><td>0.347913771591</td><td>25.0</td><td>0.0</td><td>3227.59162244</td><td>1.11273974981</td><td>-0.115984005994</td></tr>\n",
       "<tr><td>10</td><td>1492.34117879</td><td>88.076552297</td><td>0.403327315803</td><td>0.0561076819001</td><td>0.486178359365</td><td>25.0</td><td>0.0</td><td>4414.12417045</td><td>1.67920723014</td><td>-0.562760738787</td></tr>\n",
       "<tr><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td></tr>\n",
       "<tr><td>106</td><td>827.079932634</td><td>1179.02025525</td><td>0.397659638526</td><td>-0.258297321853</td><td>0.70898599022</td><td>25.0</td><td>0.0</td><td>5463.571649</td><td>1.83315370736</td><td>-0.657997203678</td></tr>\n",
       "<tr><td>107</td><td>37.7068283536</td><td>1181.08561658</td><td>0.491163497203</td><td>-0.900016865609</td><td>0.596737124036</td><td>25.0</td><td>0.0</td><td>4461.48138946</td><td>1.29027419392</td><td>-0.276705028177</td></tr>\n",
       "<tr><td>108</td><td>36.6950532505</td><td>1181.68436243</td><td>0.438245603628</td><td>-0.469566088309</td><td>0.743915646567</td><td>25.0</td><td>0.0</td><td>4282.61376842</td><td>1.46004410038</td><td>-0.410914934433</td></tr>\n",
       "<tr><td>109</td><td>1614.10324382</td><td>1208.14146027</td><td>0.447498568029</td><td>0.299919428033</td><td>0.425186845443</td><td>25.0</td><td>0.0</td><td>3437.80892385</td><td>1.01760469367</td><td>-0.0189477536934</td></tr>\n",
       "<tr><td>110</td><td>1198.9960206</td><td>1224.78385165</td><td>0.504992145083</td><td>0.0904748634633</td><td>0.649853049223</td><td>25.0</td><td>0.0</td><td>4715.86564124</td><td>1.49606300697</td><td>-0.437374710765</td></tr>\n",
       "<tr><td>111</td><td>823.282215118</td><td>1226.31768657</td><td>0.472764055771</td><td>-0.260089223041</td><td>0.663208681197</td><td>25.0</td><td>0.0</td><td>11893.6554079</td><td>4.12646557823</td><td>-1.53894556676</td></tr>\n",
       "<tr><td>112</td><td>754.38068389</td><td>1230.61620576</td><td>0.459792598755</td><td>-0.404474536098</td><td>0.484569005181</td><td>25.0</td><td>0.0</td><td>4524.82913537</td><td>1.32212938817</td><td>-0.303184896983</td></tr>\n",
       "<tr><td>113</td><td>30.9364489476</td><td>1238.4237432</td><td>0.514502694717</td><td>-0.407484064484</td><td>0.772663088555</td><td>25.0</td><td>0.0</td><td>6368.25112169</td><td>1.93519483193</td><td>-0.716811738856</td></tr>\n",
       "<tr><td>114</td><td>152.107544222</td><td>1241.39331104</td><td>0.439041849586</td><td>-0.824367161386</td><td>0.626689872888</td><td>25.0</td><td>0.0</td><td>4610.07089741</td><td>1.35224097488</td><td>-0.32763522889</td></tr>\n",
       "<tr><td>115</td><td>630.434529709</td><td>1245.65972882</td><td>0.382036296985</td><td>-0.344707119657</td><td>0.696812392358</td><td>25.0</td><td>0.0</td><td>3931.92700694</td><td>1.17175061378</td><td>-0.172087974203</td></tr>\n",
       "</table>"
      ],
      "text/plain": [
       "<Table length=115>\n",
       "  id    xcentroid     ycentroid   ...      flux           mag       \n",
       "int64    float64       float64    ...    float64        float64     \n",
       "----- ------------- ------------- ... ------------- ----------------\n",
       "    1  645.55658029 7.69908327045 ...  28.269908551   -3.62831100903\n",
       "    2 1309.16598231 11.5995843996 ... 2.85692255025   -1.13974616759\n",
       "    3 1375.65792899 11.9995956879 ... 1.05632128534 -0.0594900777583\n",
       "    4 1122.94782151 16.7158994909 ... 3.19746747319   -1.26201533809\n",
       "    5  1362.4645143 26.9619775656 ...  8.4265184611   -2.31412044147\n",
       "    6 611.037214433 28.7876045811 ... 1.98099565096  -0.742208805242\n",
       "    7 358.189493857 34.4617814911 ... 34.2481732847   -3.83659353055\n",
       "    8 449.539598592 82.1069177762 ... 2.73814589578   -1.09364146198\n",
       "    9 964.423830549 83.9894147141 ... 1.11273974981  -0.115984005994\n",
       "   10 1492.34117879  88.076552297 ... 1.67920723014  -0.562760738787\n",
       "  ...           ...           ... ...           ...              ...\n",
       "  106 827.079932634 1179.02025525 ... 1.83315370736  -0.657997203678\n",
       "  107 37.7068283536 1181.08561658 ... 1.29027419392  -0.276705028177\n",
       "  108 36.6950532505 1181.68436243 ... 1.46004410038  -0.410914934433\n",
       "  109 1614.10324382 1208.14146027 ... 1.01760469367 -0.0189477536934\n",
       "  110  1198.9960206 1224.78385165 ... 1.49606300697  -0.437374710765\n",
       "  111 823.282215118 1226.31768657 ... 4.12646557823   -1.53894556676\n",
       "  112  754.38068389 1230.61620576 ... 1.32212938817  -0.303184896983\n",
       "  113 30.9364489476  1238.4237432 ... 1.93519483193  -0.716811738856\n",
       "  114 152.107544222 1241.39331104 ... 1.35224097488   -0.32763522889\n",
       "  115 630.434529709 1245.65972882 ... 1.17175061378  -0.172087974203"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# extract sources from image and add details to comp_stars\n",
    "fwhm = 3.0\n",
    "source_snr = 20\n",
    "hdulist = fits.open(FITS_FILE)\n",
    "data = hdulist[0].data.astype(float)\n",
    "header = hdulist[0].header\n",
    "wcs = WCS(header)\n",
    "bkg_sigma = mad_std(data)    \n",
    "daofind = DAOStarFinder(fwhm=fwhm, threshold=source_snr*bkg_sigma)    \n",
    "sources = daofind(data)\n",
    "\n",
    "sources"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 6 Find the sources that correspond to our target and comparison stars\n",
    "\n",
    "Any source within 4 arc seconds of a comparison star is considered a match. Here the `results` will be converted to a pandas dataframe and will now contain the x,y coords of the target and comparison stars."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>auid</th>\n",
       "      <th>dec</th>\n",
       "      <th>error</th>\n",
       "      <th>peak</th>\n",
       "      <th>ra</th>\n",
       "      <th>vmag</th>\n",
       "      <th>x</th>\n",
       "      <th>y</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>000-BLS-551</td>\n",
       "      <td>44:22:48.0</td>\n",
       "      <td>0.054</td>\n",
       "      <td>19038.106969</td>\n",
       "      <td>20:06:48.09</td>\n",
       "      <td>11.263</td>\n",
       "      <td>388.415326</td>\n",
       "      <td>817.537270</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>000-BML-045</td>\n",
       "      <td>44:27:05.0</td>\n",
       "      <td>0.036</td>\n",
       "      <td>8153.566262</td>\n",
       "      <td>20:06:36.68</td>\n",
       "      <td>12.415</td>\n",
       "      <td>491.879800</td>\n",
       "      <td>617.615307</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>000-BLS-555</td>\n",
       "      <td>44:30:54.0</td>\n",
       "      <td>0.050</td>\n",
       "      <td>5919.023955</td>\n",
       "      <td>20:06:21.25</td>\n",
       "      <td>12.789</td>\n",
       "      <td>628.526243</td>\n",
       "      <td>440.515839</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>000-BML-046</td>\n",
       "      <td>44:35:41.8</td>\n",
       "      <td>0.048</td>\n",
       "      <td>4777.952321</td>\n",
       "      <td>20:06:53.14</td>\n",
       "      <td>12.983</td>\n",
       "      <td>365.886101</td>\n",
       "      <td>203.792448</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>000-BML-047</td>\n",
       "      <td>44:25:53.7</td>\n",
       "      <td>0.038</td>\n",
       "      <td>4538.973603</td>\n",
       "      <td>20:06:00.39</td>\n",
       "      <td>13.327</td>\n",
       "      <td>797.954483</td>\n",
       "      <td>683.615847</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>target</td>\n",
       "      <td>44 27 24.793</td>\n",
       "      <td>NaN</td>\n",
       "      <td>13457.943342</td>\n",
       "      <td>20 06 15.4553</td>\n",
       "      <td>NaN</td>\n",
       "      <td>672.395700</td>\n",
       "      <td>607.642881</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "          auid           dec  error          peak             ra    vmag  \\\n",
       "0  000-BLS-551    44:22:48.0  0.054  19038.106969    20:06:48.09  11.263   \n",
       "1  000-BML-045    44:27:05.0  0.036   8153.566262    20:06:36.68  12.415   \n",
       "2  000-BLS-555    44:30:54.0  0.050   5919.023955    20:06:21.25  12.789   \n",
       "3  000-BML-046    44:35:41.8  0.048   4777.952321    20:06:53.14  12.983   \n",
       "4  000-BML-047    44:25:53.7  0.038   4538.973603    20:06:00.39  13.327   \n",
       "5       target  44 27 24.793    NaN  13457.943342  20 06 15.4553     NaN   \n",
       "\n",
       "            x           y  \n",
       "0  388.415326  817.537270  \n",
       "1  491.879800  617.615307  \n",
       "2  628.526243  440.515839  \n",
       "3  365.886101  203.792448  \n",
       "4  797.954483  683.615847  \n",
       "5  672.395700  607.642881  "
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "for star in results:\n",
    "    star_coord = SkyCoord(star['ra'],star['dec'], unit=(u.hourangle, u.deg))\n",
    "    xy = SkyCoord.to_pixel(star_coord, wcs=wcs, origin=1)\n",
    "    x = xy[0].item(0)\n",
    "    y = xy[1].item(0)\n",
    "    for source in sources:\n",
    "        if(source['xcentroid']-4 < x < source['xcentroid']+4) and source['ycentroid']-4 < y < source['ycentroid']+4:\n",
    "            star['x'] = x\n",
    "            star['y'] = y\n",
    "            star['peak'] = source['peak']\n",
    "results = pd.DataFrame(results)\n",
    "results"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 7 Plot the image annotated with the comparsion stars\n",
    "This is just to check everything looks ok"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {

I'm working through your code on my computer to understand all the bits and ran into a roadblock:

FITS_FILE = '/home/dokeeffe/Pictures/CalibratedLight/2017-10-02/KIC08462852-2017-10-02-20-56-22Light_005.fits'

is obviously a file on your local computer but you did not seem to include it in this archive. How can I get it so I can reproduce/finish studying your code?

Thank you,
Steven