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 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]:
lc.plot()
Out[43]:
We can access several of the metadata properties:
In [94]:
lc.keplerid
Out[94]:
In [95]:
lc.channel
Out[95]:
In [96]:
lc.quarter
Out[96]:
We can access the time and flux as arrays:
In [97]:
lc.time[:10]
Out[97]:
In [98]:
lc.flux[:10]
Out[98]:
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)
detrended_lc.plot()
Out[86]:
In [92]:
folded_lc = detrended_lc.fold(period=0.837495, phase=0.92)
folded_lc.plot();
We can also compute the CDPP noise metric:
In [88]:
lc.cdpp()
Out[88]:
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/'
'200100000/82000/ktwo200182949-c14_lpd-targ.fits.gz')
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/'
'200100000/82000/ktwo200182949-c14_lpd-targ.fits.gz',
quality_bitmask=KeplerQualityFlags.CONSERVATIVE_BITMASK)
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/'
'200100000/82000/ktwo200182949-c14_lpd-targ.fits.gz',
aperture_mask='kepler-pipeline',
quality_bitmask=KeplerQualityFlags.CONSERVATIVE_BITMASK)
In [14]:
tpf.aperture_mask
Out[14]:
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]:
tpf.header(ext=0)
Out[15]:
The pixel fluxes time series can be accessed using the flux property:
In [16]:
tpf.flux.shape
Out[16]:
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]:
tpf.plot(frame=1)
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)
Out[19]:
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/'
'c14/200100000/82000/ktwo200182949-c14_llc.fits',
quality_bitmask=KeplerQualityFlags.CONSERVATIVE_BITMASK)
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)')
Out[23]:
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)')
Out[25]:
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)')
Out[27]:
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
In [3]:
module_output_to_channel(module=19, output=3)
Out[3]:
In [4]:
channel_to_module_output(67)
Out[4]:
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]:
KeplerQualityFlags.decode(1)
Out[6]:
It also can handle multiple flags:
In [7]:
KeplerQualityFlags.decode(1 + 1024 + 1048576)
Out[7]:
A few quality flags are already computed:
In [8]:
KeplerQualityFlags.decode(KeplerQualityFlags.DEFAULT_BITMASK)
Out[8]:
In [9]:
KeplerQualityFlags.decode(KeplerQualityFlags.CONSERVATIVE_BITMASK)
Out[9]:
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/'
'201500000/43000/ktwo201543306-c14_lpd-targ.fits.gz',
quality_bitmask=KeplerQualityFlags.CONSERVATIVE_BITMASK)
In [30]:
tpf.plot(frame=100)
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]:
scene.plot(*unif_prior.mean)
In [35]:
prf_phot = PRFPhotometry(scene_model=scene, prior=unif_prior)
In [36]:
results = prf_phot.fit(tpf.flux + tpf.flux_bkg)
In [37]:
plt.imshow(prf_phot.residuals[1], origin='lower')
plt.colorbar()
Out[37]:
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)')
Out[39]:
In [40]:
plt.figure(figsize=[17, 4])
plt.plot(tpf.time, xcenter)
plt.ylabel('Column position')
plt.xlabel('Time (BJD - 2454833)')
Out[40]:
In [41]:
plt.figure(figsize=[17, 4])
plt.plot(tpf.time, ycenter)
plt.ylabel('Row position')
plt.xlabel('Time (BJD - 2454833)')
Out[41]:
In [42]:
plt.figure(figsize=[17, 4])
plt.plot(tpf.time, bkg_density)
plt.ylabel('Background density')
plt.xlabel('Time (BJD - 2454833)')
Out[42]:
For more examples on PRF photometry, please see our tutorials page: http://pyke.keplerscience.org/tutorials/index.html