Intermediate: Finding Flares and Variable Stars in TASOC Light Curves#
This notebook demostrates MAST’s programmatic tools for accessing TESS time series data while exploring a flaring star from the literature (Günther et al 2019) and a nearby variable star.
The following topics will be covered:
Using the MAST API to get mission pipeline and TASOC light curves
Plotting TESS light curves in Python
Using the MAST API to make an full frame image (FFI) cutout
Creating a movie of TPF frames in Python
Using the MAST API to get a list of TESS Input Catalog (TIC) sources
Over plotting TIC sources on TESS images
Exploring the period of a variable star
Table of contents#
Exploring a known stellar flare
Querying MAST
Downloading the light curves
Plotting the light curves
Making an animation
Exploring a variable star
Overlaying TESS Input Catalog sources
Getting the variable star light curve
Plotting the variable star light curve
Finding the period of the variable star
Additional Resources
About this Notebook
Terminology#
TESS: The Transiting Exoplanet Survey Satellite
TASOC: The TESS Asteroseismic Science Operations Center
Sector: TESS observed the sky in regions of 24x96 degrees along the southern, then northern, ecliptic hemispheres. Each of these regions is referred to as a “sector”, starting with Sector 1.
TIC: The TESS input catalog.
FFI: TESS periodically reads out the entire frame of all four cameras, nominally every 30 minutes, and stores them as full frame images (FFIs).
TPF: Target Pixel File, a fits file containing stacks of small images centered on a target, one image for every timestamp the telescope took data.
HDU: Header Data Unit. A FITS file is made up of HDUs that contain data and metadata relating to the file. The first HDU is called the primary HDU, and anything that follows is considered an “extension”, e.g., “the first FITS extension”, “the second FITS extension”, etc.
HDUList: A list of HDUs that comprise a fits file.
BJD: Barycentric Julian Date, the Julian Date that has been corrected for differences in the Earth’s position with respect to the Solar System center of mass.
BTJD: Barycentric TESS Julian Date, the timestamp measured in BJD, but offset by 2457000.0. I.e., BTJD = BJD - 2457000.0
WCS: World Coordinate System, the coordinates that locate an astronomical object on the sky.
Imports#
In this notbook we will use the MAST module of Astroquery (astroquery.mast) to query and download data, and various astropy and numpy functions to manipulate the data. We will use both the matplotlib and bokeh plotting packages to visualize our data as they have different strengths and weaknesses.
# For querying for data
from astroquery.mast import Tesscut, Observations, Catalogs
# For manipulating data
import numpy as np
from astropy.table import Table
from astropy.coordinates import SkyCoord
from astropy.io import fits
from astropy.wcs import WCS
from astropy.timeseries import LombScargle
from astropy.time import Time
import astropy.units as u
# For matplotlib plotting
import matplotlib
%matplotlib inline
import matplotlib.pyplot as plt
import matplotlib.animation as animation
# For animation display
from matplotlib import rc
from IPython.display import HTML
rc('animation', html='jshtml')
# For bokeh plotting
from bokeh import plotting
from bokeh.models import Span
plotting.output_notebook()
Exploring a stellar flare#
Selecting the flare#
We will start with a known flare from the literature, in this case from Günther, M. N., Zhan, Z., Seager, S., et al. 2019, arXiv e-prints, arXiv:1901.00443. We picked a particularly long flare to give us the best chance of finding it in the 30 minute cadence data as well as the 2 minute cadence data.
We’ve made note of the TIC ID and sector for our flare of interest, as well as its peak time in BJD format:
tic_id = 141914082
sector = 1
tpeak = 2458341.89227 # Julian Day
Querying MAST#
Mission light curves#
Using the query_criteria function in the astroquery.mast.Observations class, we will specify that we are looking for TESS mission data using the obs_collection argument, that we want a specific TIC ID using the target_name argument, and that we want observations from a particular sector using the sequence_number argument.
mission_res = Observations.query_criteria(obs_collection="TESS",
target_name=tic_id,
sequence_number=sector)
mission_res
| intentType | obs_collection | provenance_name | instrument_name | project | filters | wavelength_region | target_name | target_classification | obs_id | s_ra | s_dec | dataproduct_type | proposal_pi | calib_level | t_min | t_max | t_exptime | em_min | em_max | obs_title | t_obs_release | proposal_id | proposal_type | sequence_number | s_region | jpegURL | dataURL | dataRights | mtFlag | srcDen | obsid | objID |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| str7 | str4 | str4 | str10 | str4 | str4 | str7 | str9 | str1 | str47 | float64 | float64 | str10 | str14 | int64 | float64 | float64 | float64 | float64 | float64 | str1 | float64 | str23 | str1 | int64 | str41 | str1 | str73 | str6 | bool | float64 | str8 | str9 |
| science | TESS | SPOC | Photometer | TESS | TESS | Optical | 141914082 | -- | tess2018206045859-s0001-0000000141914082-0120-s | 94.6175354047675 | -72.0448462219924 | timeseries | Ricker, George | 3 | 58324.79323696759 | 58352.67632608796 | 120.0 | 600.0 | 1000.0 | -- | 58458.5833333 | G011175_G011180_G011176 | -- | 1 | CIRCLE 94.6175354 -72.04484622 0.00138889 | -- | mast:TESS/product/tess2018206045859-s0001-0000000141914082-0120-s_lc.fits | PUBLIC | False | nan | 60829534 | 110344410 |
TASOC light curves#
In addition to mission pipeline data, MAST also hosts a variety of community contributed High Level Science Products (HLSPs), all of which are given the mission (obs_collection) “HLSP”. In this case we will specifically search for HLSPs in the TESS project, which will return the light curves provided by the TASOC (note the provenance_name of “TASOC”). All other arguments remain the same.
tasoc_res = Observations.query_criteria(target_name=tic_id,
obs_collection="HLSP",
project="TESS",
sequence_number=sector)
tasoc_res['dataproduct_type',"obs_collection","target_name","t_exptime","filters",
"provenance_name","project","sequence_number","instrument_name"]
| dataproduct_type | obs_collection | target_name | t_exptime | filters | provenance_name | project | sequence_number | instrument_name |
|---|---|---|---|---|---|---|---|---|
| str10 | str4 | str9 | float64 | str4 | str17 | str4 | int64 | str10 |
| timeseries | HLSP | 141914082 | 1800.0 | TESS | TESS-SPOC | TESS | 1 | Photometer |
| timeseries | HLSP | 141914082 | 120.0 | TESS | TASOC | TESS | 1 | Photometer |
| timeseries | HLSP | 141914082 | 1800.0 | TESS | TASOC | TESS | 1 | Photometer |
| timeseries | HLSP | 141914082 | 1800.0 | TESS | GSFC-ELEANOR-LITE | TESS | 1 | Photometer |
| timeseries | HLSP | 141914082 | 1800.0 | TESS | TGLC | TESS | 1 | Photometer |
| timeseries | HLSP | 141914082 | 1800.0 | TESS | QLP | TESS | 1 | Photometer |
This query returns two light curves. To understand the difference between the two light curves we look at the t_exptime column, and note the different values. These exposure times correspond to 2 minutes (short cadence) and 30 minutes (long cadence). We will explore both light curves.
Downloading the light curves#
For the rest of this notebook we will work with the TASOC light curves only, although we could do the same analysis on the mission light curves.
Querying for the list of associated data products#
Each observation may have one or more associated data products. In the case of the TASOC light curves, there is simply a single light curve file for each observation.
tasoc_prod = Observations.get_product_list(tasoc_res)
tasoc_prod["dataproduct_type", "description", "dataURI", "size"]
| dataproduct_type | description | dataURI | size |
|---|---|---|---|
| str10 | str4 | str129 | int64 |
| timeseries | FITS | mast:HLSP/qlp/s0001/0000/0001/4191/4082/hlsp_qlp_tess_ffi_s0001-0000000141914082_tess_v01_llc.fits | 80640 |
| timeseries | Text | mast:HLSP/qlp/s0001/0000/0001/4191/4082/hlsp_qlp_tess_ffi_s0001-0000000141914082_tess_v01_llc.txt | 57710 |
| timeseries | FITS | mast:HLSP/tess-spoc/s0001/target/0000/0001/4191/4082/hlsp_tess-spoc_tess_phot_0000000141914082-s0001_tess_v1_lc.fits | 164160 |
| timeseries | FITS | mast:HLSP/tess-spoc/s0001/target/0000/0001/4191/4082/hlsp_tess-spoc_tess_phot_0000000141914082-s0001_tess_v1_tp.fits | 3925440 |
| timeseries | FITS | mast:HLSP/tasoc/s0001/c0120/0000/0001/4191/4082/hlsp_tasoc_tess_tpf_tic00141914082-s0001-cam4-ccd2-c0120_tess_v05_cbv-lc.fits | 1880640 |
| timeseries | FITS | mast:HLSP/tasoc/s0001/c1800/0000/0001/4191/4082/hlsp_tasoc_tess_ffi_tic00141914082-s0001-cam4-ccd2-c1800_tess_v05_cbv-lc.fits | 164160 |
| timeseries | FITS | mast:HLSP/tasoc/s0001/c1800/0000/0001/4191/4082/hlsp_tasoc_tess_ffi_tic00141914082-s0001-cam4-ccd2-c1800_tess_v05_ens-lc.fits | 167040 |
| timeseries | FITS | mast:HLSP/gsfc-eleanor-lite/s0001/0000/0001/4191/4082/hlsp_gsfc-eleanor-lite_tess_ffi_s0001-0000000141914082_tess_v1.0_lc.fits | 106560 |
| timeseries | FITS | mast:HLSP/tglc/s0001/cam4-ccd2/0052/6627/0443/4424/hlsp_tglc_tess_ffi_gaiaid-5266270443442455040-s0001-cam4-ccd2_tess_v1_llc.fits | 57600 |
Downloading files#
We can choose to download some or all of the associated data files, in this case since we just have the two light curves, we will download all of the products.
tasoc_manifest = Observations.download_products(tasoc_prod)
tasoc_manifest
Downloading URL https://mast.stsci.edu/api/v0.1/Download/file?uri=mast:HLSP/gsfc-eleanor-lite/s0001/0000/0001/4191/4082/hlsp_gsfc-eleanor-lite_tess_ffi_s0001-0000000141914082_tess_v1.0_lc.fits to ./mastDownload/HLSP/hlsp_gsfc-eleanor-lite_tess_ffi_s0001-0000000141914082_tess_v1.0_lc/hlsp_gsfc-eleanor-lite_tess_ffi_s0001-0000000141914082_tess_v1.0_lc.fits ...
[Done]
Downloading URL https://mast.stsci.edu/api/v0.1/Download/file?uri=mast:HLSP/qlp/s0001/0000/0001/4191/4082/hlsp_qlp_tess_ffi_s0001-0000000141914082_tess_v01_llc.fits to ./mastDownload/HLSP/hlsp_qlp_tess_ffi_s0001-0000000141914082_tess_v01_llc/hlsp_qlp_tess_ffi_s0001-0000000141914082_tess_v01_llc.fits ...
[Done]
Downloading URL https://mast.stsci.edu/api/v0.1/Download/file?uri=mast:HLSP/qlp/s0001/0000/0001/4191/4082/hlsp_qlp_tess_ffi_s0001-0000000141914082_tess_v01_llc.txt to ./mastDownload/HLSP/hlsp_qlp_tess_ffi_s0001-0000000141914082_tess_v01_llc/hlsp_qlp_tess_ffi_s0001-0000000141914082_tess_v01_llc.txt ...
[Done]
Downloading URL https://mast.stsci.edu/api/v0.1/Download/file?uri=mast:HLSP/tasoc/s0001/c0120/0000/0001/4191/4082/hlsp_tasoc_tess_tpf_tic00141914082-s0001-cam4-ccd2-c0120_tess_v05_cbv-lc.fits to ./mastDownload/HLSP/hlsp_tasoc_tess_tpf_tic00141914082-s0001-cam4-ccd2-c0120_tess_v05/hlsp_tasoc_tess_tpf_tic00141914082-s0001-cam4-ccd2-c0120_tess_v05_cbv-lc.fits ...
[Done]
Downloading URL https://mast.stsci.edu/api/v0.1/Download/file?uri=mast:HLSP/tasoc/s0001/c1800/0000/0001/4191/4082/hlsp_tasoc_tess_ffi_tic00141914082-s0001-cam4-ccd2-c1800_tess_v05_cbv-lc.fits to ./mastDownload/HLSP/hlsp_tasoc_tess_ffi_tic00141914082-s0001-cam4-ccd2-c1800_tess_v05/hlsp_tasoc_tess_ffi_tic00141914082-s0001-cam4-ccd2-c1800_tess_v05_cbv-lc.fits ...
[Done]
Downloading URL https://mast.stsci.edu/api/v0.1/Download/file?uri=mast:HLSP/tasoc/s0001/c1800/0000/0001/4191/4082/hlsp_tasoc_tess_ffi_tic00141914082-s0001-cam4-ccd2-c1800_tess_v05_ens-lc.fits to ./mastDownload/HLSP/hlsp_tasoc_tess_ffi_tic00141914082-s0001-cam4-ccd2-c1800_tess_v05/hlsp_tasoc_tess_ffi_tic00141914082-s0001-cam4-ccd2-c1800_tess_v05_ens-lc.fits ...
[Done]
Downloading URL https://mast.stsci.edu/api/v0.1/Download/file?uri=mast:HLSP/tess-spoc/s0001/target/0000/0001/4191/4082/hlsp_tess-spoc_tess_phot_0000000141914082-s0001_tess_v1_lc.fits to ./mastDownload/HLSP/hlsp_tess-spoc_tess_phot_0000000141914082-s0001_tess_v1_tp/hlsp_tess-spoc_tess_phot_0000000141914082-s0001_tess_v1_lc.fits ...
[Done]
Downloading URL https://mast.stsci.edu/api/v0.1/Download/file?uri=mast:HLSP/tess-spoc/s0001/target/0000/0001/4191/4082/hlsp_tess-spoc_tess_phot_0000000141914082-s0001_tess_v1_tp.fits to ./mastDownload/HLSP/hlsp_tess-spoc_tess_phot_0000000141914082-s0001_tess_v1_tp/hlsp_tess-spoc_tess_phot_0000000141914082-s0001_tess_v1_tp.fits ...
[Done]
Downloading URL https://mast.stsci.edu/api/v0.1/Download/file?uri=mast:HLSP/tglc/s0001/cam4-ccd2/0052/6627/0443/4424/hlsp_tglc_tess_ffi_gaiaid-5266270443442455040-s0001-cam4-ccd2_tess_v1_llc.fits to ./mastDownload/HLSP/hlsp_tglc_tess_ffi_gaiaid-5266270443442455040-s0001-cam4-ccd2_tess_v1_llc/hlsp_tglc_tess_ffi_gaiaid-5266270443442455040-s0001-cam4-ccd2_tess_v1_llc.fits ...
[Done]
| Local Path | Status | Message | URL |
|---|---|---|---|
| str172 | str8 | object | object |
| ./mastDownload/HLSP/hlsp_gsfc-eleanor-lite_tess_ffi_s0001-0000000141914082_tess_v1.0_lc/hlsp_gsfc-eleanor-lite_tess_ffi_s0001-0000000141914082_tess_v1.0_lc.fits | COMPLETE | None | None |
| ./mastDownload/HLSP/hlsp_qlp_tess_ffi_s0001-0000000141914082_tess_v01_llc/hlsp_qlp_tess_ffi_s0001-0000000141914082_tess_v01_llc.fits | COMPLETE | None | None |
| ./mastDownload/HLSP/hlsp_qlp_tess_ffi_s0001-0000000141914082_tess_v01_llc/hlsp_qlp_tess_ffi_s0001-0000000141914082_tess_v01_llc.txt | COMPLETE | None | None |
| ./mastDownload/HLSP/hlsp_tasoc_tess_tpf_tic00141914082-s0001-cam4-ccd2-c0120_tess_v05/hlsp_tasoc_tess_tpf_tic00141914082-s0001-cam4-ccd2-c0120_tess_v05_cbv-lc.fits | COMPLETE | None | None |
| ./mastDownload/HLSP/hlsp_tasoc_tess_ffi_tic00141914082-s0001-cam4-ccd2-c1800_tess_v05/hlsp_tasoc_tess_ffi_tic00141914082-s0001-cam4-ccd2-c1800_tess_v05_cbv-lc.fits | COMPLETE | None | None |
| ./mastDownload/HLSP/hlsp_tasoc_tess_ffi_tic00141914082-s0001-cam4-ccd2-c1800_tess_v05/hlsp_tasoc_tess_ffi_tic00141914082-s0001-cam4-ccd2-c1800_tess_v05_ens-lc.fits | COMPLETE | None | None |
| ./mastDownload/HLSP/hlsp_tess-spoc_tess_phot_0000000141914082-s0001_tess_v1_tp/hlsp_tess-spoc_tess_phot_0000000141914082-s0001_tess_v1_lc.fits | COMPLETE | None | None |
| ./mastDownload/HLSP/hlsp_tess-spoc_tess_phot_0000000141914082-s0001_tess_v1_tp/hlsp_tess-spoc_tess_phot_0000000141914082-s0001_tess_v1_tp.fits | COMPLETE | None | None |
| ./mastDownload/HLSP/hlsp_tglc_tess_ffi_gaiaid-5266270443442455040-s0001-cam4-ccd2_tess_v1_llc/hlsp_tglc_tess_ffi_gaiaid-5266270443442455040-s0001-cam4-ccd2_tess_v1_llc.fits | COMPLETE | None | None |
Plotting the light curves#
We will use the bokeh plotting library so that we can have interactivity, and will plot both the 2 minute and 30 minute cadence data together.
We can tell which file corresponds to which cadence length by examining the filenames and noting that one contains c0120 (2 minutes) and the other c1800 (30 minutes).
# Loading the short cadence light curve
hdu = fits.open(tasoc_manifest["Local Path"][0])
short_cad_lc = Table(hdu[1].data)
hdu.close()
# Loading the long cadence light curve
hdu = fits.open(tasoc_manifest["Local Path"][1])
long_cad_lc = Table(hdu[1].data)
hdu.close()
long_cad_lc
| TIME | CADENCENO | SAP_FLUX | KSPSAP_FLUX | KSPSAP_FLUX_ERR | QUALITY | ORBITID | SAP_X | SAP_Y | SAP_BKG | SAP_BKG_ERR | KSPSAP_FLUX_SML | KSPSAP_FLUX_LAG |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| float64 | int32 | float32 | float32 | float32 | int32 | int32 | float32 | float32 | float32 | float32 | float32 | float32 |
| 1325.3239205802092 | 4697 | 0.9907516 | 0.99937916 | 0.0004727787 | 4096 | 9 | 833.876 | 1609.046 | 323.96 | 545.07 | 0.9994108 | 0.99934494 |
| 1325.3447537297802 | 4698 | 0.99065953 | 0.99991655 | 0.0004727787 | 4096 | 9 | 833.877 | 1609.046 | 292.76 | 430.71 | 0.99993587 | 0.99988496 |
| 1325.3655868809078 | 4699 | 0.9909295 | 1.0005782 | 0.0004727787 | 4096 | 9 | 833.877 | 1609.046 | 180.94 | 480.8 | 1.0005702 | 1.0005482 |
| 1325.386420033564 | 4700 | 0.9906624 | 1.0004753 | 0.0004727787 | 4096 | 9 | 833.877 | 1609.045 | 18.43 | 442.45 | 1.0004779 | 1.0004911 |
| 1325.4072531877214 | 4701 | 0.99046475 | 1.0002398 | 0.0004727787 | 4096 | 9 | 833.878 | 1609.046 | 28.41 | 461.86 | 1.0002539 | 1.0002934 |
| 1325.4280863433216 | 4702 | 0.9909264 | 1.0004859 | 0.0004727787 | 4096 | 9 | 833.878 | 1609.044 | -89.99 | 462.61 | 1.0004859 | 1.0006546 |
| 1325.4489195003468 | 4703 | 0.99120694 | 1.000383 | 0.0004727787 | 4096 | 9 | 833.878 | 1609.044 | -118.24 | 472.78 | 1.0004015 | 1.0002135 |
| 1325.4697526587393 | 4704 | 0.9914977 | 1.000145 | 0.0004727787 | 4096 | 9 | 833.879 | 1609.044 | -307.69 | 311.05 | 1.000151 | 1.0001317 |
| 1325.4905858184716 | 4705 | 0.99213433 | 1.0001283 | 0.0004727787 | 4096 | 9 | 833.879 | 1609.044 | -228.2 | 359.26 | 1.0001414 | 1.0002065 |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 1352.9695117277513 | 6024 | 1.0127677 | 0.9998823 | 0.0004727787 | 0 | 10 | 833.884 | 1609.154 | 200.05 | 784.49 | 0.99989283 | 1.0001168 |
| 1352.990344774581 | 6025 | 1.0124162 | 0.9999271 | 0.0004727787 | 0 | 10 | 833.884 | 1609.155 | 181.46 | 665.6 | 0.9999403 | 1.0000422 |
| 1353.0111778210035 | 6026 | 1.0118057 | 0.999762 | 0.0004727787 | 0 | 10 | 833.882 | 1609.153 | 193.32 | 571.46 | 0.9997629 | 0.9999271 |
| 1353.0320108670521 | 6027 | 1.011702 | 1.000145 | 0.0004727787 | 0 | 10 | 833.884 | 1609.157 | 258.93 | 701.06 | 1.0001379 | 1.0001836 |
| 1353.0528439127938 | 6028 | 1.0111037 | 1.0000861 | 0.0004727787 | 0 | 10 | 833.883 | 1609.153 | 226.06 | 637.77 | 1.0001204 | 1.0000738 |
| 1353.0736769582584 | 6029 | 1.0105393 | 1.0001072 | 0.0004727787 | 0 | 10 | 833.883 | 1609.156 | 267.1 | 719.3 | 1.0000957 | 1.0001739 |
| 1353.0945100035185 | 6030 | 1.011731 | 1.0019132 | 0.0004727787 | 0 | 10 | 833.883 | 1609.156 | 237.25 | 576.67 | 1.00187 | 1.0020751 |
| 1353.1153430486104 | 6031 | 1.0092119 | 1.0000905 | 0.0004727787 | 0 | 10 | 833.882 | 1609.156 | 220.36 | 776.46 | 1.0000939 | 1.0001388 |
| 1353.136176093613 | 6032 | 1.0085595 | 1.0001625 | 0.0004727787 | 4096 | 10 | 833.883 | 1609.156 | 229.04 | 772.09 | 1.0001705 | 1.0002179 |
| 1353.1570091385684 | 6033 | 1.0075595 | 0.99993503 | 0.0004727787 | 4096 | 10 | 833.882 | 1609.157 | 417.53 | 742.14 | 0.99995434 | 0.99971807 |
long_cad_lc.columns
<TableColumns names=('TIME','CADENCENO','SAP_FLUX','KSPSAP_FLUX','KSPSAP_FLUX_ERR','QUALITY','ORBITID','SAP_X','SAP_Y','SAP_BKG','SAP_BKG_ERR','KSPSAP_FLUX_SML','KSPSAP_FLUX_LAG')>
bfig = plotting.figure(width=850, height=300, title=f"Detrended Lightcurve (TIC{tic_id})")
# Short cadence
bfig.circle(short_cad_lc["TIME"], short_cad_lc["RAW_FLUX"]/np.median(short_cad_lc["RAW_FLUX"]), fill_color="black",size=2, line_color=None)
bfig.line(short_cad_lc["TIME"], short_cad_lc["RAW_FLUX"]/np.median(short_cad_lc["RAW_FLUX"]), line_color='black')
# Long cadence
bfig.circle(long_cad_lc["TIME"],long_cad_lc["SAP_FLUX"], fill_color="#0384f7",size=6, line_color=None)
bfig.line(long_cad_lc["TIME"],long_cad_lc["SAP_FLUX"], line_color='#0384f7')
# Marking the flare (tpeak is in BJD, while the time column in the light curve is BTJD, so we must convert)
vline = Span(location=(tpeak - 2457000), dimension='height', line_color='#bf006e', line_width=3)
bfig.renderers.extend([vline])
# Labeling the axes
bfig.xaxis.axis_label = "Time (BTJD)"
bfig.yaxis.axis_label = "Flux"
plotting.show(bfig)
BokehDeprecationWarning: 'circle() method with size value' was deprecated in Bokeh 3.4.0 and will be removed, use 'scatter(size=...) instead' instead.
BokehDeprecationWarning: 'circle() method with size value' was deprecated in Bokeh 3.4.0 and will be removed, use 'scatter(size=...) instead' instead.
Experiment with the controls on the right to zoom in and and pan around on this light curve to look at the marked flare and other features (some sort of periodic feature maybe?)
Making an animation#
Looking at the above plot we can see the flare event in both the long and short cadence light curves. Since we can see it even in the 30 minute data, we should be able to make an animation of the area around the flaring star and see the flare happen.
We will use TESScut, the MAST cutout tool for full-frame images to cutout the area around the flaring star across the entire sector, and then make a movie that shows how it changes over time.
We will use the astroquery.mast Tesscut class to make this cutout.
We will use two functions:
Find the sky coordinate of our flare star:
Observations._resolve_object*Query for cutouts and get the result as a list of HDUList objects:
Tesscut.get_cutouts**
* Observations._resolve_object is a private (not documented) function which is being removed in favor of the public function Observations.resolve_object in the next Astroquery release.
** We must start by finding the sky coordinate of our star, however starting with the next Astroquery release, Tesscut functions will be able to take an object name such as a TIC ID as well.
coord = Observations.resolve_object(f"TIC {tic_id}")
**Requesting a cutout target pixel file. **
This query will return a list of HDUList objects, each of which is the cutout target pixel file (TPF) for a single sector. In this case, because we specified a single sector we know that the resulting list will only have one element and can pull it out directly.
cutout_hdu = Tesscut.get_cutouts(coordinates=coord, size=40, sector=1)[0]
cutout_hdu.info()
Filename: <class '_io.BytesIO'>
No. Name Ver Type Cards Dimensions Format
0 PRIMARY 1 PrimaryHDU 57 ()
1 PIXELS 1 BinTableHDU 281 1282R x 12C [D, E, J, 1600J, 1600E, 1600E, 1600E, 1600E, J, E, E, 38A]
2 APERTURE 1 ImageHDU 82 (40, 40) int32
A TESScut TPF has three extensions:
No. 0 (Primary): This HDU contains meta-data related to the entire file.
No. 1 (Pixels): This HDU contains a binary table that holds data like cutout image arrays and times. We will extract information from here to get our set of cutout images.
No. 2 (Aperture): This HDU contains the image extension with data collected from the aperture. We will use this extension to get the WCS associated with out cutout.
cutout_table = Table(cutout_hdu[1].data)
cutout_table.columns
<TableColumns names=('TIME','TIMECORR','CADENCENO','RAW_CNTS','FLUX','FLUX_ERR','FLUX_BKG','FLUX_BKG_ERR','QUALITY','POS_CORR1','POS_CORR2','FFI_FILE')>
The Pixels extension contains the cutout data table, which has a number of columns. For our puposes we care about the “TIME” column which has the observation times in BTJD, and the “FLUX” column which contains the cutout images (they are calibrated, but not background subtracted).
Exploring the cutout time series#
We want to explore what is happening within our cutout area over the time that the flare occurs, so we will make an animated plot of the cutout frames.
We can’t make a movie of the whole sector (it would take too long), so we will choose only the time range around the flare.
Use the light curve plot to figure out what time range you want to animate, or use our selections.
start_btjd = 1341.5
end_btjd = 1342.5
start_index = (np.abs(cutout_table['TIME'] - start_btjd)).argmin()
end_index = (np.abs(cutout_table['TIME'] - end_btjd)).argmin()
print(f"Frames {start_index}-{end_index} ({end_index-start_index} frames)")
Frames 721-769 (48 frames)
Looking at the animated cutout#
def make_animation(data_array, start_frame=0, end_frame=None, vmin=None, vmax=None, delay=50):
"""
Function that takes an array where each frame is a 2D image array and make an animated plot
that runs through the frames.
Note: This can take a long time to run if you have a lot of frames.
Parameters
----------
data_array : array
Array of 2D images.
start_frame : int
The index of the initial frame to show. Default is the first frame.
end_frame : int
The index of the final frame to show. Default is the last frame.
vmin : float
Data range min for the colormap. Defaults to data minimum value.
vmax : float
Data range max for the colormap. Defaults to data maximum value.
delay:
Delay before the next frame is shown in milliseconds.
Returns
-------
response : `animation.FuncAnimation`
"""
if not vmin:
vmin = np.min(data_array)
if not vmax:
vmax = np.max(data_array)
if not end_frame:
end_frame = len(data_array) - 1 # set to the end of the array
num_frames = end_frame - start_frame + 1 # include the end frame
def animate(i, fig, ax, binarytab, start=0):
"""Function used to update the animation"""
ax.set_title("Epoch #" + str(i+start))
im = ax.imshow(binarytab[i+start], cmap=plt.cm.YlGnBu_r, vmin=vmin, vmax=vmax)
return im,
# Create initial plot.
fig, ax = plt.subplots(figsize=(10,10))
ax.imshow(data_array[start_frame], cmap=plt.cm.YlGnBu_r, vmin=vmin, vmax=vmax)
ani = animation.FuncAnimation(fig, animate, fargs=(fig, ax, data_array, start_frame), frames=num_frames,
interval=delay, repeat_delay=1000)
plt.close()
return ani
make_animation(cutout_table['FLUX'], start_index, end_index, vmax=700, delay=150)
We can see three things in this plot:
The flare that occures in frames 740-743
An aberration that appears in frame 754
A variable star pulsing in the lower right corner
Exploring a variable star#
Now we will look more closely at the variable star we can see in the animation. We will use the TESS Input Catalog (TIC) to figure out the TESS ID of the variable star, and then using that ID pull down and explore the star’s light curve from MAST.
Overlaying TIC sources#
First we will use the astroquery.mast.Catalog class to query the TESS Input Catalog (TIC) and get a list of sources that appear in our cutout.
Here we do a simple cone search and the apply a magnitude limit. Play around with the magnitude limit to see how it changes the number of sources.
sources = Catalogs.query_object(catalog="TIC", objectname=f"TIC {tic_id}", radius=10*u.arcmin)
sources = sources[sources["Tmag"] < 12]
print(f"Number of sources: {len(sources)}")
print(sources)
Number of sources: 9
ID ra dec ... wdflag dstArcSec
--------- ---------------- ----------------- ... ------ ------------------
141914082 94.6175354047675 -72.0448462219924 ... 0 0.0
141914038 94.7353652417206 -72.0837957851839 ... 0 191.637404585262
141914103 94.7749592275334 -72.0227841314024 ... 0 192.00656069294254
141914130 94.5074642175418 -71.998930582098 ... 0 205.6248622381558
141914317 94.9032150266819 -72.0334742503164 ... 0 319.770063872556
141913994 94.5588022852622 -72.1321979465057 ... 0 321.11922083466345
166975135 94.6067354129356 -71.9255734581109 ... 0 429.55027273236567
141869504 94.2183733431791 -72.0232095513896 ... 0 450.0317520158321
141913929 94.7099836100023 -72.1924762405753 ... 0 541.2030879067767
Plotting sources on an individual cutout#
We will get the WCS infomation associated with our cutout from the Aperture extension, and use it to make a WCS-aware plot of a single cutout image. Then we display the image and sources together, and label the sources with their row number in the catalog table.
cutout_wcs = WCS(cutout_hdu[2].header)
cutout_img = cutout_table["FLUX"][start_index]
fig, ax = plt.subplots(subplot_kw={'projection':cutout_wcs})
fig.set_size_inches(10,10)
plt.grid(color='white', ls='solid')
# Setup WCS axes.
xcoords = ax.coords[0]
ycoords = ax.coords[1]
xcoords.set_major_formatter('d.ddd')
ycoords.set_major_formatter('d.ddd')
xcoords.set_axislabel("RA (deg)")
ycoords.set_axislabel("Dec (deg)")
ax.imshow(cutout_img, cmap=plt.cm.YlGnBu_r,vmin=0,vmax=700)
ax.plot(sources['ra'],sources['dec'],'x',transform=ax.get_transform('icrs'),color="red")
# Annotating the sources with their row nnumber in the sources table
for i,star in enumerate(sources):
ax.text(star['ra']+0.01,star['dec'],i,transform=ax.get_transform('icrs'))
ax.set_xlim(0,cutout_img.shape[1]-1)
ax.set_ylim(cutout_img.shape[0]-1,0)
plt.show()
The variable star is row 4 in the catalog sources table (note that if you changed the magnitude threshold the variale star may be in a different row).
sources["ID","ra","dec"][4]
| ID | ra | dec |
|---|---|---|
| str11 | float64 | float64 |
| 141914317 | 94.9032150266819 | -72.0334742503164 |
Getting the variable star light curve#
Again, we will look specifically for the TASOC light curve(s) associated with this star, rather than the mission pipeline ones. Below we go through the same process as in the Downloading the light curves section to search for the observation, then find the associated data products, and download them.
variable_tic_id = sources["ID"][4]
variable_res = Observations.query_criteria(target_name=str(variable_tic_id),
obs_collection="HLSP",
filters="TESS")
print(f"Number of tasoc light curves for {variable_tic_id}: {len(variable_res)}\n")
variable_prod = Observations.get_product_list(variable_res[0])
variable_manifest = Observations.download_products(variable_prod)
Number of tasoc light curves for 141914317: 83
Downloading URL https://mast.stsci.edu/api/v0.1/Download/file?uri=mast:HLSP/qlp/s0013/0000/0001/4191/4317/hlsp_qlp_tess_ffi_s0013-0000000141914317_tess_v01_llc.fits to ./mastDownload/HLSP/hlsp_qlp_tess_ffi_s0013-0000000141914317_tess_v01_llc/hlsp_qlp_tess_ffi_s0013-0000000141914317_tess_v01_llc.fits ...
[Done]
Downloading URL https://mast.stsci.edu/api/v0.1/Download/file?uri=mast:HLSP/qlp/s0013/0000/0001/4191/4317/hlsp_qlp_tess_ffi_s0013-0000000141914317_tess_v01_llc.txt to ./mastDownload/HLSP/hlsp_qlp_tess_ffi_s0013-0000000141914317_tess_v01_llc/hlsp_qlp_tess_ffi_s0013-0000000141914317_tess_v01_llc.txt ...
[Done]
Note that this time there is only one (30 minute cadence) TASOC light curve. This is because this star was not one of the targets that TESS observed at the shorter cadence.
hdu = fits.open(variable_manifest["Local Path"][0])
variable_lc = Table(hdu[1].data)
hdu.close()
Plotting the variable star light curve#
We will again plot the light curve using bokeh so that we have interactivity.
variable_lc
| TIME | CADENCENO | SAP_FLUX | KSPSAP_FLUX | KSPSAP_FLUX_ERR | QUALITY | ORBITID | SAP_X | SAP_Y | SAP_BKG | SAP_BKG_ERR | KSPSAP_FLUX_SML | KSPSAP_FLUX_LAG |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| float64 | int32 | float32 | float32 | float32 | int32 | int32 | float32 | float32 | float32 | float32 | float32 | float32 |
| 1653.9281264088759 | 20470 | 1.0268497 | 1.0558715 | 0.13665526 | 4128 | 33 | 421.72015 | 1216.5792 | 55.87 | 303.0 | 1.0491947 | 1.0580579 |
| 1653.948959649118 | 20471 | 0.9360396 | 0.9612099 | 0.13665526 | 4096 | 33 | 421.72424 | 1216.5784 | 63.58 | 413.66 | 0.96619904 | 0.95966774 |
| 1653.969792891784 | 20472 | 0.8209844 | 0.8419681 | 0.13665526 | 4096 | 33 | 421.73187 | 1216.5752 | -218.77 | 387.51 | 0.86041814 | 0.83654326 |
| 1653.9906261366607 | 20473 | 0.7616569 | 0.7801408 | 0.13665526 | 4096 | 33 | 421.735 | 1216.5725 | -0.87 | 454.05 | 0.80635834 | 0.7720926 |
| 1654.011459383893 | 20474 | 0.8281789 | 0.8472413 | 0.13665526 | 4096 | 33 | 421.73022 | 1216.5757 | -157.54 | 447.66 | 0.8650209 | 0.84183574 |
| 1654.0322926332728 | 20475 | 0.94315326 | 0.96371984 | 0.13665526 | 4096 | 33 | 421.7236 | 1216.5797 | 84.09 | 448.03 | 0.96824825 | 0.96247375 |
| 1654.053125884903 | 20476 | 1.0280553 | 1.0492675 | 0.13665526 | 4096 | 33 | 421.71915 | 1216.5809 | -48.59 | 350.11 | 1.0434113 | 1.0512749 |
| 1654.0739591385864 | 20477 | 1.0823321 | 1.1034349 | 0.13665526 | 4096 | 33 | 421.71863 | 1216.5817 | -22.17 | 435.48 | 1.0914406 | 1.1074061 |
| 1654.0947923943963 | 20478 | 1.105596 | 1.1259379 | 0.13665526 | 4096 | 33 | 421.7166 | 1216.583 | 42.58 | 475.55 | 1.1116827 | 1.1304736 |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 1682.1571617084755 | 21825 | 0.8778642 | 0.8501837 | 0.13665526 | 4096 | 34 | 421.72787 | 1216.628 | -97.16 | 516.83 | 0.8661454 | 0.84522027 |
| 1682.1779948309704 | 21826 | 0.7011507 | 0.6742332 | 0.13665526 | 4132 | 34 | 421.72543 | 1216.6162 | -2064.34 | 1480.94 | 0.70495605 | 0.66895413 |
| 1682.1988279531406 | 21827 | 0.7898172 | 0.76078343 | 0.13665526 | 0 | 34 | 421.73325 | 1216.6271 | -102.71 | 397.85 | 0.7869868 | 0.7528833 |
| 1682.2196610751407 | 21828 | 0.89272857 | 0.8574418 | 0.13665526 | 0 | 34 | 421.7258 | 1216.6309 | 20.56 | 437.69 | 0.8730871 | 0.8527917 |
| 1682.2404941969792 | 21829 | 0.9951512 | 0.95296043 | 0.13665526 | 0 | 34 | 421.7207 | 1216.6326 | 223.22 | 546.28 | 0.95803046 | 0.95133626 |
| 1682.2613273188597 | 21830 | 1.0591404 | 1.0110862 | 0.13665526 | 0 | 34 | 421.7172 | 1216.6337 | 74.16 | 340.84 | 1.0100105 | 1.0116068 |
| 1682.2821604407818 | 21831 | 1.1034971 | 1.0500327 | 0.13665526 | 0 | 34 | 421.71664 | 1216.6342 | 250.13 | 395.14 | 1.0449181 | 1.0516453 |
| 1682.3029935630043 | 21832 | 1.1116096 | 1.0542127 | 0.13665526 | 0 | 34 | 421.71613 | 1216.6346 | 306.66 | 411.17 | 1.0488179 | 1.0560246 |
| 1682.3238266855037 | 21833 | 1.0830014 | 1.0235174 | 0.13665526 | 0 | 34 | 421.71698 | 1216.6335 | 295.34 | 449.06 | 1.0217227 | 1.0241388 |
| 1682.3446598085977 | 21834 | 1.0243448 | 0.96460146 | 0.13665526 | 4096 | 34 | 421.719 | 1216.6324 | 146.27 | 408.53 | 0.9687561 | 0.96340925 |
bfig = plotting.figure(width=850, height=300, title=f"Detrended Lightcurve (TIC{variable_tic_id})")
bfig.circle(variable_lc["TIME"],variable_lc["SAP_FLUX"], fill_color="black",size=4, line_color=None)
bfig.line(variable_lc["TIME"],variable_lc["SAP_FLUX"], line_color='black')
# Labeling the axes
bfig.xaxis.axis_label = "Time (BTJD)"
bfig.yaxis.axis_label = "SAP Flux"
plotting.show(bfig)
BokehDeprecationWarning: 'circle() method with size value' was deprecated in Bokeh 3.4.0 and will be removed, use 'scatter(size=...) instead' instead.
That looks variable all right!
Finding the period#
We’ll run a Lomb Scargle priodogram on this light curve to see if we can quantify the periodic behavior. To do this we will use the astropy.timeseries class LombScargle.
lomb = LombScargle(variable_lc["TIME"], variable_lc["SAP_FLUX"])
frequency, power = lomb.autopower(maximum_frequency=25)
Plotting the periodogram#
bfig = plotting.figure(width=850, height=300, x_range=(0,25),
title=f"Periodogram (TIC{variable_tic_id})")
bfig.line(frequency, power, line_color='black')
# Labeling the axes
bfig.xaxis.axis_label = "Frequency (1/day)"
bfig.yaxis.axis_label = "Power"
plotting.show(bfig)
Phasing on the highest power period/frequency#
There is a clear dominant frequency in the above plot, with a small harmonic also visible. We will phase the stellar light curve on the period corresponding to the dominant frequency and plot both it and the corresponding sinusoidal fit.
dominant_freq = frequency[np.argmax(power)].value
print(f"The dominant period: {1/dominant_freq*24:.3} hours")
The dominant period: 5.24 hours
bfig = plotting.figure(width=850, height=300, title=f"Phased Lightcurve (TIC{variable_tic_id})")
# Plotting the phased light curve
bfig.circle(variable_lc["TIME"]%(1/dominant_freq),variable_lc["SAP_FLUX"], fill_color="black",size=4, line_color=None)
# Plotting the periodic fit
t_fit = np.linspace(0,1/dominant_freq,100)
bfig.line(t_fit, lomb.model(t_fit, dominant_freq), color='#1b9f00', line_width=2)
# Labeling the axes
bfig.xaxis.axis_label = "Phase (days)"
bfig.yaxis.axis_label = "Flux"
plotting.show(bfig)
BokehDeprecationWarning: 'circle() method with size value' was deprecated in Bokeh 3.4.0 and will be removed, use 'scatter(size=...) instead' instead.
Additional Resources#
About this Notebook#
Author: C. E. Brasseur, STScI Software Engineer
Last Updated: May 2023
Top of Page
