In [1]:

    

%matplotlib inline
import numpy as np
from scipy import stats
import matplotlib 
import matplotlib.pyplot as plt
from matplotlib.colors import LogNorm
import matplotlib.cm as cm
import pandas as pd
from astropy.io import fits
from astropy.stats import LombScargle
from gatspy import datasets, periodic


    

In [2]:

    

matplotlib.rcParams.update({'font.size':18}) 
matplotlib.rcParams.update({'font.family':'serif'})


    

In [3]:

    

target = fits.open('kplr010536761-2009166043257_llc.fits')


    

In [4]:

    

target.info()


    

    




Filename: kplr010536761-2009166043257_llc.fits
No.    Name         Type      Cards   Dimensions   Format
  0  PRIMARY     PrimaryHDU      58   ()      
  1  LIGHTCURVE  BinTableHDU    155   1639R x 20C   [D, E, J, E, E, E, E, E, E, J, D, E, D, E, D, E, D, E, E, E]   
  2  APERTURE    ImageHDU        48   (5, 4)   int32   

In [5]:

    

target[0].header


    

    
Out[5]:





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    = '2015-09-01'         / file creation date.                            
CREATOR = '796031 FluxExporter2PipelineModule' / pipeline job and program used t
PROCVER = 'svn+ssh://murzim/repo/soc/tags/release/9.3.22 r60269' / 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  = 'KIC 10536761'       / string version of target id                    
KEPLERID=             10536761 / unique Kepler target identifier                
CHANNEL =                   58 / CCD channel                                    
SKYGROUP=                   18 / roll-independent location of channel           
MODULE  =                   17 / CCD module                                     
OUTPUT  =                    2 / CCD output                                     
QUARTER =                    1 / Observing quarter                              
SEASON  =                    3 / mission season during which data was collected 
DATA_REL=                   25 / data release version number                    
OBSMODE = 'long cadence'       / observing mode                                 
MISSION = 'Kepler  '           / Mission name                                   
TTABLEID=                   20 / target table id                                
RADESYS = 'ICRS    '           / reference frame of celestial coordinates       
RA_OBJ  =           292.547199 / [deg] right ascension                          
DEC_OBJ =            47.733840 / [deg] declination                              
EQUINOX =               2000.0 / equinox of celestial coordinate system         
PMRA    =               0.0000 / [arcsec/yr] RA proper motion                   
PMDEC   =               0.0000 / [arcsec/yr] Dec proper motion                  
PMTOTAL =               0.0000 / [arcsec/yr] total proper motion                
PARALLAX=                      / [arcsec] parallax                              
GLON    =            79.902195 / [deg] galactic longitude                       
GLAT    =            13.692583 / [deg] galactic latitude                        
GMAG    =               16.295 / [mag] SDSS g band magnitude                    
RMAG    =               14.930 / [mag] SDSS r band magnitude                    
IMAG    =               13.881 / [mag] SDSS i band magnitude                    
ZMAG    =               13.311 / [mag] SDSS z band magnitude                    
D51MAG  =               15.992 / [mag] D51 magnitude,                           
JMAG    =               11.990 / [mag] J band magnitude from 2MASS              
HMAG    =               11.391 / [mag] H band magnitude from 2MASS              
KMAG    =               11.163 / [mag] K band magnitude from 2MASS              
KEPMAG  =               14.605 / [mag] Kepler magnitude (Kp)                    
GRCOLOR =                1.365 / [mag] (g-r) color, SDSS bands                  
JKCOLOR =                0.827 / [mag] (J-K) color, 2MASS bands                 
GKCOLOR =                5.132 / [mag] (g-K) color, SDSS g - 2MASS K            
TEFF    =                 3461 / [K] Effective temperature                      
LOGG    =                4.870 / [cm/s2] log10 surface gravity                  
FEH     =                0.000 / [log10([Fe/H])]  metallicity                   
EBMINUSV=                0.038 / [mag] E(B-V) reddening                         
AV      =                0.119 / [mag] A_v extinction                           
RADIUS  =                0.380 / [solar radii] stellar radius                   
TMINDEX =            328893856 / unique 2MASS catalog ID                        
SCPID   =                      / unique SCP processing ID                       
CHECKSUM= 'B7QdC4ObB4ObB4Ob'   / HDU checksum updated 2015-09-01T14:43:53Z      

In [6]:

    

data = target[1].data


    

In [7]:

    

time = data['TIME'] 
flux = data['PDCSAP_FLUX'] 
rel_flux = (flux-np.nanmedian(flux))/np.nanmedian(flux)
err = data['PDCSAP_FLUX_ERR']
rel_err = (err-np.nanmedian(err))/np.nanmedian(err)


    

In [8]:

    

fig = plt.figure(figsize = (11,5))

plt.plot(time,rel_flux,color='b')

plt.xlim(150,160)
plt.ylim(-0.02,0.05)
plt.xlabel('Time [BJD-2454833]')
plt.ylabel('Relative Flux')
plt.title('KIC 10536761')


    

    
Out[8]:





<matplotlib.text.Text at 0x119552048>






    

In [9]:

    

cut = np.where(np.isnan(flux) == False)


    

In [10]:

    

np.size(time[cut])


    

    
Out[10]:





1624

In [11]:

    

fig = plt.figure(figsize = (11,5))

plt.plot(time[cut],rel_flux[cut],color='b')

#plt.xlim(150,156)
plt.ylim(-0.02,0.05)
plt.xlabel('Time [BJD-2454833]')
plt.ylabel('Relative Flux')
plt.title('KIC 10536761')

plt.savefig('KIC10536761_lc.pdf', dpi = 600)


    

In [12]:

    

model = periodic.LombScargleFast(fit_period = True)

model.optimizer.period_range = (0.2, 6.0)

model.fit(time[cut], rel_flux[cut], err[cut])


    

    




Finding optimal frequency:
 - Estimated peak width = 0.188
 - Using 5 steps per peak; omega_step = 0.0375
 - User-specified period range:  0.2 to 6
 - Computing periods at 810 steps
Zooming-in on 5 candidate peaks:
 - Computing periods at 1000 steps






    
Out[12]:





<gatspy.periodic.lomb_scargle_fast.LombScargleFast at 0x10b27e710>

In [13]:

    

model.best_period


    

    
Out[13]:





5.6495940473148138

In [20]:

    

# Compute the scores on a grid of periods
periods = np.linspace(0.2, 10.0, 10000)

import warnings
with warnings.catch_warnings():
    warnings.simplefilter("ignore")
    scores = model.score(periods)

# Plot the results
fig, ax = plt.subplots(figsize=(11, 5))
fig.subplots_adjust(bottom=0.2)
ax.plot(periods, scores, 'red')
ax.set(xlabel='period (days)', ylabel='Lomb Scargle Power', title='L-S Periodogram',
       xlim=(0.0, 10.0), ylim=(0, 1), xscale = 'log')

plt.savefig('KIC10536761_periodogram.pdf', dpi = 600)


    

In [ ]:

    




    

In [ ]: