What’s new in PyKE 3?

Developed since 2012, PyKE offers a user-friendly way to inspect and analyze the pixels and lightcurves obtained by NASA’s Kepler and K2.

The latest version of PyKE, v3.1, was released in January 2018 and adds a new object-oriented Python API which is intended to aid the development of custom pipelines and tools by the community.

In [47]:
%load_ext autoreload
%autoreload 2
%matplotlib inline
from matplotlib import rcParams
rcParams["figure.figsize"] = (14, 5)
The autoreload extension is already loaded. To reload it, use:
  %reload_ext autoreload

Introducing a generic LightCurve class

The most notable change is the introduction of a generic LightCurve class which provides operations that are intended to suit time series data from any astronomical survey. A light curve is simply instantiated as follows:

In [45]:
from pyke import LightCurve
lc = LightCurve(time=[1, 2, 3], flux=[78.4, 79.6, 76.5])

A LightCurve object provides easy access to a range of common operations, such as fold(), flatten(), remove_outliers(), cdpp(), plot(), and more. To demonstrate these operations, let’s create a LightCurve object from a KeplerLightCurveFile we obtain from the data archive at MAST:

In [46]:
from pyke import KeplerLightCurveFile
lcfile = KeplerLightCurveFile("https://archive.stsci.edu/missions/kepler/lightcurves/0119/011904151/kplr011904151-2010009091648_llc.fits")
lc = lcfile.SAP_FLUX

Now lc is a LightCurve object on which you can run operations. For example, we can plot it:

In [43]:
<matplotlib.axes._subplots.AxesSubplot at 0x7f4aadfd9b00>

We can access several of the metadata properties:

In [94]:
In [95]:
In [96]:

We can access the time and flux as arrays:

In [97]:
array([ 352.37710511,  352.39753829,  352.43840464,  352.45883781,
        352.47927099,  352.49970417,  352.52013735,  352.54057042,
        352.5610037 ,  352.58143688])
In [98]:
array([ 562954.5625,  562910.3125,  562950.375 ,  562939.1875,
        562964.5625,  562967.    ,  562952.1875,  562950.375 ,
        562979.4375,  562973.4375], dtype=float32)

We don’t particularly care about the long-term trends, so let’s use a Savitzky-Golay filter to flatten the lightcurve:

In [86]:
detrended_lc, _ = lc.flatten(polyorder=1)
<matplotlib.axes._subplots.AxesSubplot at 0x7f4aac1e1d30>
In [92]:
folded_lc = detrended_lc.fold(period=0.837495, phase=0.92)

We can also compute the CDPP noise metric:

In [88]:

Target Pixel File (TPF)

PyKE 3.1 includes class called KeplerTargetPixelFile which is used to handle target pixel files:

In [10]:
from pyke import KeplerTargetPixelFile

A KeplerTargetPixelFile can be instantiated either from a local file or a url:

In [11]:
tpf = KeplerTargetPixelFile('https://archive.stsci.edu/missions/k2/target_pixel_files/c14/'

Additionally, we can mask out cadences that are flagged using the quality_bitmask argument in the constructor:

In [12]:
tpf = KeplerTargetPixelFile('https://archive.stsci.edu/missions/k2/target_pixel_files/c14/'

Furthermore, we can mask out pixel values using the aperture_mask argument. The default behaviour is to use all pixels that have real values. This argument can also get a string value 'kepler-pipeline', in which case the default aperture used by Kepler’s pipeline is applied.

In [13]:
tpf = KeplerTargetPixelFile('https://archive.stsci.edu/missions/k2/target_pixel_files/c14/'
In [14]:
array([[False, False, False, ..., False, False, False],
       [False, False, False, ..., False, False, False],
       [False, False, False, ..., False, False, False],
       [False, False, False, ..., False, False, False],
       [False, False, False, ..., False, False, False],
       [False, False, False, ..., False, False, False]], dtype=bool)

The TPF objects stores both data and a few metadata information, e.g., channel number, EPIC number, reference column and row, module, and shape. The whole header is also available:

In [15]:
SIMPLE  =                    T / conforms to FITS standards
BITPIX  =                    8 / array data type
NAXIS   =                    0 / number of array dimensions
EXTEND  =                    T / file contains extensions
NEXTEND =                    2 / number of standard extensions
EXTNAME = 'PRIMARY '           / name of extension
EXTVER  =                    1 / extension version number (not format version)
ORIGIN  = 'NASA/Ames'          / institution responsible for creating this file
DATE    = '2017-10-30'         / file creation date.
CREATOR = '472097 TargetPixelExporterPipelineModule' / pipeline job and program
PROCVER = 'svn+ssh://murzim/repo/soc/tags/release/9.3.75 r62933' / SW version
FILEVER = '6.1     '           / file format version
TIMVERSN= 'OGIP/93-003'        / OGIP memo number for file format
TELESCOP= 'Kepler  '           / telescope
INSTRUME= 'Kepler Photometer'  / detector type
OBJECT  = 'EPIC 200182949'     / string version of target id
KEPLERID=            200182949 / unique Kepler target identifier
CHANNEL =                    1 / CCD channel
MODULE  =                    2 / CCD module
OUTPUT  =                    1 / CCD output
CAMPAIGN=                   14 / Observing campaign number
DATA_REL=                   20 / data release version number
OBSMODE = 'long cadence'       / observing mode
MISSION = 'K2      '           / Mission name
TTABLEID=                   92 / target table id
RADESYS = 'ICRS    '           / reference frame of celestial coordinates
RA_OBJ  =           165.706242 / [deg] right ascension
DEC_OBJ =             2.769444 / [deg] declination
EQUINOX =               2000.0 / equinox of celestial coordinate system
PMRA    =                      / [arcsec/yr] RA proper motion
PMDEC   =                      / [arcsec/yr] Dec proper motion
PMTOTAL =                      / [arcsec/yr] total proper motion
PARALLAX=                      / [arcsec] parallax
GLON    =                      / [deg] galactic longitude
GLAT    =                      / [deg] galactic latitude
GMAG    =                      / [mag] SDSS g band magnitude
RMAG    =                      / [mag] SDSS r band magnitude
IMAG    =                      / [mag] SDSS i band magnitude
ZMAG    =                      / [mag] SDSS z band magnitude
JMAG    =                      / [mag] J band magnitude from 2MASS
HMAG    =                      / [mag] H band magnitude from 2MASS
KMAG    =                      / [mag] K band magnitude from 2MASS
KEPMAG  =                      / [mag] Kepler magnitude (Kp)
GRCOLOR =                      / [mag] (g-r) color, SDSS bands
JKCOLOR =                      / [mag] (J-K) color, 2MASS bands
GKCOLOR =                      / [mag] (g-K) color, SDSS g - 2MASS K
TEFF    =                      / [K] Effective temperature
LOGG    =                      / [cm/s2] log10 surface gravity
FEH     =                      / [log10([Fe/H])]  metallicity
EBMINUSV=                      / [mag] E(B-V) reddening
AV      =                      / [mag] A_v extinction
RADIUS  =                      / [solar radii] stellar radius
TMINDEX =                      / unique 2MASS catalog ID
CHECKSUM= 'QkmmQijlQijlQijl'   / HDU checksum updated 2017-10-30T20:42:50Z

The pixel fluxes time series can be accessed using the flux property:

In [16]:
(3209, 35, 35)

This shows that our TPF is a 35 x 35 image recorded over 3209 cadences.

One can visualize the pixel data at a given cadence using the plot method:

In [17]:

We can perform aperture photometry using the method to_lightcurve:

In [18]:
lc = tpf.to_lightcurve()
In [19]:
plt.figure(figsize=[17, 4])
plt.plot(lc.time, lc.flux)
[<matplotlib.lines.Line2D at 0x129c079e8>]

Let’s see how the previous light curve compares against the 'SAP_FLUX' produced by Kepler’s pipeline. For that, we are going to explore the KeplerLightCurveFile class:

In [20]:
from pyke.lightcurve import KeplerLightCurveFile
In [21]:
klc = KeplerLightCurveFile('https://archive.stsci.edu/missions/k2/lightcurves/'
In [22]:
sap_lc = klc.SAP_FLUX
In [23]:
plt.figure(figsize=[17, 4])
plt.plot(lc.time, lc.flux)
plt.plot(sap_lc.time, sap_lc.flux)
plt.ylabel('Flux (e-/s)')
plt.xlabel('Time (BJD - 2454833)')
<matplotlib.text.Text at 0x12a0af400>

Now, let’s correct this light curve using by fitting cotrending basis vectors. That can be achived either with the KeplerCBVCorrector class or the compute_cotrended_lightcurve in KeplerLightCurveFile. Let’s try the latter:

In [24]:
klc_corrected = klc.compute_cotrended_lightcurve(cbvs=range(1, 17))
In [25]:
plt.figure(figsize=[17, 4])
plt.plot(klc_corrected.time, klc_corrected.flux)
plt.ylabel('Flux (e-/s)')
plt.xlabel('Time (BJD - 2454833)')
<matplotlib.text.Text at 0x12a1c1c50>
In [26]:
pdcsap_lc = klc.PDCSAP_FLUX
In [27]:
plt.figure(figsize=[17, 4])
plt.plot(klc_corrected.time, klc_corrected.flux)
plt.plot(pdcsap_lc.time, pdcsap_lc.flux)
plt.ylabel('Flux (e-/s)')
plt.xlabel('Time (BJD - 2454833)')
<matplotlib.text.Text at 0x12a0bc0f0>

Utility functions

PyKE has included two convinience functions to convert between module.output to channel and vice-versa:

In [1]:
%matplotlib inline
import matplotlib.pyplot as plt
In [2]:
from pyke.utils import module_output_to_channel, channel_to_module_output
/Users/jvmirca/anaconda3/lib/python3.6/site-packages/matplotlib/__init__.py:1405: UserWarning:
This call to matplotlib.use() has no effect because the backend has already
been chosen; matplotlib.use() must be called *before* pylab, matplotlib.pyplot,
or matplotlib.backends is imported for the first time.

In [3]:
module_output_to_channel(module=19, output=3)
In [4]:
(19, 3)

PyKE 3.1 includes KeplerQualityFlags class which encodes the meaning of the Kepler QUALITY bitmask flags as documented in the Kepler Archive Manual (Table 2.3):

In [5]:
from pyke.utils import KeplerQualityFlags
In [6]:
['Attitude tweak']

It also can handle multiple flags:

In [7]:
KeplerQualityFlags.decode(1 + 1024 + 1048576)
['Attitude tweak', 'Sudden sensitivity dropout', 'Thruster firing']

A few quality flags are already computed:

In [8]:
['Attitude tweak',
 'Safe mode',
 'Coarse point',
 'Earth point',
 'Desaturation event',
 'Cosmic ray in optimal aperture',
 'Manual exclude',
 'No data',
 'Thruster firing']
In [9]:
['Attitude tweak',
 'Safe mode',
 'Coarse point',
 'Earth point',
 'Desaturation event',
 'Cosmic ray in optimal aperture',
 'Manual exclude',
 'Sudden sensitivity dropout',
 'Cosmic ray in collateral data',
 'No data',
 'Possible thruster firing',
 'Thruster firing']

Pixel Response Function (PRF) Photometry

PyKE 3.1 also includes tools to perform PRF Photometry:

In [28]:
from pyke.prf import PRFPhotometry, SceneModel, SimpleKeplerPRF

For that, let’s create a SceneModel which will be fitted to the object of the following TPF:

In [29]:
tpf = KeplerTargetPixelFile('https://archive.stsci.edu/missions/k2/target_pixel_files/c14/'
In [30]:
In [31]:
scene = SceneModel(prfs=[SimpleKeplerPRF(channel=tpf.channel, shape=tpf.shape[1:],
                                         column=tpf.column, row=tpf.row)])

We also need to define prior distributions on the parameters of our SceneModel model. Those parameters are the flux, center positions of the target, and a constant background level. We can do that with oktopus:

In [32]:
from oktopus.prior import UniformPrior
In [33]:
unif_prior = UniformPrior(lb=[0, 1090., 706., 0.],
                          ub=[1e5, 1096., 712., 1e5])
In [34]:
In [35]:
prf_phot = PRFPhotometry(scene_model=scene, prior=unif_prior)
In [36]:
results = prf_phot.fit(tpf.flux + tpf.flux_bkg)
  0%|          | 0/3633 [00:00<?, ?it/s]/Users/jvmirca/anaconda3/lib/python3.6/site-packages/autograd/core.py:82: RuntimeWarning: invalid value encountered in log
  result_value = self.fun(*argvals, **kwargs)
/Users/jvmirca/anaconda3/lib/python3.6/site-packages/scipy/optimize/optimize.py:1927: RuntimeWarning: invalid value encountered in double_scalars
  tmp2 = (x - v) * (fx - fw)
/Users/jvmirca/anaconda3/lib/python3.6/site-packages/scipy/optimize/optimize.py:1928: RuntimeWarning: invalid value encountered in double_scalars
  p = (x - v) * tmp2 - (x - w) * tmp1
/Users/jvmirca/anaconda3/lib/python3.6/site-packages/scipy/optimize/optimize.py:1929: RuntimeWarning: invalid value encountered in double_scalars
  tmp2 = 2.0 * (tmp2 - tmp1)
100%|██████████| 3633/3633 [02:49<00:00, 21.46it/s]
In [37]:
plt.imshow(prf_phot.residuals[1], origin='lower')
<matplotlib.colorbar.Colorbar at 0x12aed6be0>
In [38]:
flux = results[:, 0]
xcenter = results[:, 1]
ycenter = results[:, 2]
bkg_density = results[:, 3]
In [39]:
plt.figure(figsize=[17, 4])
plt.plot(tpf.time, flux)
plt.ylabel('Flux (e-/s)')
plt.xlabel('Time (BJD - 2454833)')
<matplotlib.text.Text at 0x12b2960f0>
In [40]:
plt.figure(figsize=[17, 4])
plt.plot(tpf.time, xcenter)
plt.ylabel('Column position')
plt.xlabel('Time (BJD - 2454833)')
<matplotlib.text.Text at 0x128a809b0>
In [41]:
plt.figure(figsize=[17, 4])
plt.plot(tpf.time, ycenter)
plt.ylabel('Row position')
plt.xlabel('Time (BJD - 2454833)')
<matplotlib.text.Text at 0x1299983c8>
In [42]:
plt.figure(figsize=[17, 4])
plt.plot(tpf.time, bkg_density)
plt.ylabel('Background density')
plt.xlabel('Time (BJD - 2454833)')
<matplotlib.text.Text at 0x12b2bb278>

For more examples on PRF photometry, please see our tutorials page: http://pyke.keplerscience.org/tutorials/index.html