lightkurve.lightcurve.
TessLightCurve
(time=None, flux=None, flux_err=None, flux_unit=Unit("electron / s"), time_format='btjd', time_scale='tdb', centroid_col=None, centroid_row=None, quality=None, quality_bitmask=None, cadenceno=None, sector=None, camera=None, ccd=None, targetid=None, ra=None, dec=None, label=None, meta=None)¶Bases: lightkurve.lightcurve.LightCurve
Subclass of LightCurve
which holds extra data specific to the TESS mission.
Time measurements
Data flux for every time point
Uncertainty on each flux data point
Unit
or strUnit of the flux values. If a string is passed, it will be passed
on the the constructor of Unit
.
String specifying how an instant of time is represented, e.g. ‘bkjd’ or ‘jd’.
String which specifies how the time is measured, e.g. tdb’, ‘tt’, ‘ut1’, or ‘utc’.
Centroid column and row coordinates as a function of time
Array indicating the quality of each data point
Bitmask specifying quality flags of cadences that should be ignored
Cadence numbers corresponding to every time measurement
Tess Input Catalog ID number
Attributes Summary
Returns the time values as an Astropy |
|
Returns the flux as an astropy.units.Quantity object. |
|
Methods Summary
|
Append LightCurve objects. |
|
Bins a lightcurve in chunks defined by |
|
Returns a copy of the LightCurve object. |
|
Plots the light curve using Matplotlib’s |
|
Estimate the CDPP noise metric using the Savitzky-Golay (SG) method. |
|
Fill in gaps in time. |
|
Removes the low frequency trend using scipy’s Savitzky-Golay filter. |
|
Folds the lightcurve at a specified |
|
Create a new |
|
Creates a new |
|
Display an interactive Jupyter Notebook widget to find planets. |
|
Returns a normalized version of the light curve. |
|
Plot the light curve using Matplotlib’s |
|
Removes cadences where the flux is NaN. |
|
Removes outlier data points using sigma-clipping. |
|
Plots the light curve using Matplotlib’s |
|
Prints a description of all non-callable attributes. |
|
Returns a corrector object to remove instrument systematics. |
|
Writes the light curve to a CSV file. |
|
Writes the KeplerLightCurve to a FITS file. |
|
Converts the light curve to a Pandas |
|
Converts the light curve to a |
|
Returns a |
|
Returns a |
|
Converts the light curve to an Astropy |
|
Converts the light curve to an AstroPy |
Attributes Documentation
astropy_time
¶Returns the time values as an Astropy Time
object.
The Time object will be created based on the values of the light curve’s
time
, time_format
, and time_scale
attributes.
If the time_format
attribute is not set or not one of the formats
allowed by AstroPy.
Examples
The section below demonstrates working with time values using the TESS light curve of Pi Mensae as an example, which we obtained as follows:
>>> import lightkurve as lk
>>> lc = lk.search_lightcurvefile("Pi Mensae", mission="TESS", sector=1).download().PDCSAP_FLUX
>>> lc
TessLightCurve(TICID: 261136679)
Every LightCurve
object has a time
attribute, which provides access
to the original array of time values given in the native format and
scale used by the data product from which the light curve was obtained:
>>> lc.time
array([1325.29698328, 1325.29837215, 1325.29976102, ..., 1353.17431099,
1353.17569985, 1353.17708871])
>>> lc.time_format
'btjd'
>>> lc.time_scale
'tdb'
To enable users to convert these time values to different formats or scales, Lightkurve provides an easy way to access the time values as an AstroPy Time object:
>>> lc.astropy_time
<Time object: scale='tdb' format='jd' value=[2458325.29698328 2458325.29837215 2458325.29976102 ... 2458353.17431099
2458353.17569985 2458353.17708871]>
This is convenient because AstroPy Time objects provide a lot of useful features. For example, we can now obtain the Julian Day or ISO values that correspond to the raw time values:
>>> lc.astropy_time.iso
array(['2018-07-25 19:07:39.356', '2018-07-25 19:09:39.354',
'2018-07-25 19:11:39.352', ..., '2018-08-22 16:11:00.470',
'2018-08-22 16:13:00.467', '2018-08-22 16:15:00.464'], dtype='<U23')
>>> lc.astropy_time.jd
array([2458325.29698328, 2458325.29837215, 2458325.29976102, ...,
2458353.17431099, 2458353.17569985, 2458353.17708871])
flux_quantity
¶Returns the flux as an astropy.units.Quantity object.
flux_unit
¶Methods Documentation
append
(self, others, inplace=False)¶Append LightCurve objects.
Light curves to be appended to the current one.
Concatenated light curve.
bin
(self, binsize=None, bins=None, method='mean')¶Bins a lightcurve in chunks defined by binsize
or bins
.
The flux value of the bins will be computed by taking the mean
(method='mean'
) or the median (method='median'
) of the flux.
The default is mean.
Number of cadences to include in every bin. The default
is 13 if neither bins
nor binsize
is assigned.
Requires Astropy version >3.1 and >2.10
Instruction for how to assign bin locations grouping by the time of
samples rather than index; overrides the binsize=
if given.
If bins
is an int, it is the number of bins. If it is a list
it is taken to be the bin edges. If it is a string, it must be one
of ‘blocks’, ‘knuth’, ‘scott’ or ‘freedman’.
See histogram
for description of these algorithms.
The summary statistic to return for each bin. Default: ‘mean’.
LightCurve
A new light curve which has been binned.
Notes
If the ratio between the lightcurve length and the binsize is not a whole number, then the remainder of the data points will be ignored.
If the original light curve contains flux uncertainties (flux_err
),
the binned lightcurve will report the root-mean-square error.
If no uncertainties are included, the binned curve will return the
standard deviation of the data.
If the original lightcurve contains a quality attribute, then the bitwise OR of the quality flags will be returned per bin.
copy
(self)¶Returns a copy of the LightCurve object.
This method uses the copy.deepcopy
function to ensure that all
objects stored within the LightCurve are copied (e.g. time and flux).
A new LightCurve
object which is a copy of the original.
errorbar
(self, linestyle='', **kwargs)¶Plots the light curve using Matplotlib’s errorbar
method.
Axes
A matplotlib axes object to plot into. If no axes is provided, a new one will be generated.
Normalize the lightcurve before plotting?
Plot x axis label
Plot y axis label
Plot set_title
Path or URL to a matplotlib style file, or name of one of matplotlib’s built-in stylesheets (e.g. ‘ggplot’). Lightkurve’s custom stylesheet is used by default.
Connect the error bars using a line?
Dictionary of arguments to be passed to matplotlib.pyplot.scatter
.
Axes
The matplotlib axes object.
estimate_cdpp
(self, transit_duration=13, savgol_window=101, savgol_polyorder=2, sigma=5.0)¶Estimate the CDPP noise metric using the Savitzky-Golay (SG) method.
A common estimate of the noise in a lightcurve is the scatter that remains after all long term trends have been removed. This is the idea behind the Combined Differential Photometric Precision (CDPP) metric. The official Kepler Pipeline computes this metric using a wavelet-based algorithm to calculate the signal-to-noise of the specific waveform of transits of various durations. In this implementation, we use the simpler “sgCDPP proxy algorithm” discussed by Gilliland et al (2011ApJS..197….6G) and Van Cleve et al (2016PASP..128g5002V).
Remove low frequency signals using a Savitzky-Golay filter with
window length savgol_window
and polynomial order savgol_polyorder
.
Remove outliers by rejecting data points which are separated from
the mean by sigma
times the standard deviation.
Compute the standard deviation of a running mean with
a configurable window length equal to transit_duration
.
We use a running mean (as opposed to block averaging) to strongly attenuate the signal above 1/transit_duration whilst retaining the original frequency sampling. Block averaging would set the Nyquist limit to 1/transit_duration.
The transit duration in units of number of cadences. This is the length of the window used to compute the running mean. The default is 13, which corresponds to a 6.5 hour transit in data sampled at 30-min cadence.
Width of Savitsky-Golay filter in cadences (odd number). Default value 101 (2.0 days in Kepler Long Cadence mode).
Polynomial order of the Savitsky-Golay filter. The recommended value is 2.
The number of standard deviations to use for clipping outliers. The default is 5.
Savitzky-Golay CDPP noise metric in units parts-per-million (ppm).
Notes
This implementation is adapted from the Matlab version used by Jeff van Cleve but lacks the normalization factor used there: svn+ssh://murzim/repo/so/trunk/Develop/jvc/common/compute_SG_noise.m
fill_gaps
(self, method='gaussian_noise')¶Fill in gaps in time.
By default, the gaps will be filled with random white Gaussian noise distributed according to \(\mathcal{N} (\mu=\overline{\mathrm{flux}}, \sigma=\mathrm{CDPP})\). No other methods are supported at this time.
Method to use for gap filling. Fills with Gaussian noise by default.
LightCurve
A new light curve object in which all NaN values and gaps in time have been filled.
flatten
(self, window_length=101, polyorder=2, return_trend=False, break_tolerance=5, niters=3, sigma=3, mask=None, **kwargs)¶Removes the low frequency trend using scipy’s Savitzky-Golay filter.
This method wraps scipy.signal.savgol_filter
.
The length of the filter window (i.e. the number of coefficients).
window_length
must be a positive odd integer.
The order of the polynomial used to fit the samples. polyorder
must be less than window_length.
If True
, the method will return a tuple of two elements
(flattened_lc, trend_lc) where trend_lc is the removed trend.
If there are large gaps in time, flatten will split the flux into
several sub-lightcurves and apply savgol_filter
to each
individually. A gap is defined as a period in time larger than
break_tolerance
times the median gap. To disable this feature,
set break_tolerance
to None.
Number of iterations to iteratively sigma clip and flatten. If more than one, will perform the flatten several times, removing outliers each time.
Number of sigma above which to remove outliers from the flatten
Boolean array to mask data with before flattening. Flux values where mask is True will not be used to flatten the data. An interpolated result will be provided for these points. Use this mask to remove data you want to preserve, e.g. transits.
Dictionary of arguments to be passed to scipy.signal.savgol_filter
.
Flattened lightcurve.
return_trend
is True
, the method will also return:Trend in the lightcurve data
fold
(self, period, t0=None, transit_midpoint=None)¶Folds the lightcurve at a specified period
and reference time t0
.
This method returns a FoldedLightCurve
object
in which the time values range between -0.5 to +0.5 (i.e. the phase).
Data points which occur exactly at t0
or an integer multiple of
t0 + n*period
will have phase value 0.0.
The period upon which to fold, in the same units as this
LightCurve’s time
attribute.
Time corresponding to zero phase, in the same units as this
LightCurve’s time
attribute. Defaults to 0 if not set.
Deprecated. Use t0
instead.
FoldedLightCurve
A new light curve object in which the data are folded and sorted by
phase. The object contains an extra phase
attribute.
Examples
The example below shows a light curve with a period dip which occurs near
time value 1001 and has a period of 5 days. Calling the fold
method
will transform the light curve into a
FoldedLightCurve
object:
>>> import lightkurve as lk
>>> lc = lk.LightCurve(time=range(1001, 1012), flux=[0.5, 1.0, 1.0, 1.0, 1.0, 0.5, 1.0, 1.0, 1.0, 1.0, 0.5])
>>> folded_lc = lc.fold(period=5., t0=1006.)
>>> folded_lc
<lightkurve.lightcurve.FoldedLightCurve>
An object of type FoldedLightCurve
is useful
because it provides convenient access to the phase values and the
phase-folded fluxes:
>>> folded_lc.phase
array([-0.4, -0.4, -0.2, -0.2, 0. , 0. , 0. , 0.2, 0.2, 0.4, 0.4])
>>> folded_lc.flux
array([1. , 1. , 1. , 1. , 0.5, 0.5, 0.5, 1. , 1. , 1. , 1. ])
We can still access the original time values as well:
>>> folded_lc.time_original
array([1004, 1009, 1005, 1010, 1001, 1006, 1011, 1002, 1007, 1003, 1008])
A FoldedLightCurve
inherits all the features
of a standard LightCurve
object. For example, we can very quickly
obtain a phase-folded plot using:
>>> folded_lc.plot()
from_stingray
(lc)¶Create a new LightCurve
from a stingray.Lightcurve
.
stingray.Lightcurve
A stingray Lightcurve object.
from_timeseries
(ts)¶Creates a new LightCurve
from an AstroPy
TimeSeries
object.
TimeSeries
The AstroPy TimeSeries object. The object must contain columns named ‘time’, ‘flux’, and ‘flux_err’.
interact_bls
(self, notebook_url='localhost:8888', minimum_period=None, maximum_period=None, resolution=2000)¶Display an interactive Jupyter Notebook widget to find planets.
The Box Least Squares (BLS) periodogram is a statistical tool used for detecting transiting exoplanets and eclipsing binaries in light curves. This method will display a Jupyter Notebook Widget which enables the BLS algorithm to be used interactively. Behind the scenes, the widget uses the AstroPy implementation of BLS [1].
This feature only works inside an active Jupyter Notebook. It requires Bokeh v1.0 (or later) and AstroPy v3.1 (or later), which are optional dependencies. An error message will be shown if these dependencies are not available.
Location of the Jupyter notebook page (default: “localhost:8888”) When showing Bokeh applications, the Bokeh server must be explicitly configured to allow connections originating from different URLs. This parameter defaults to the standard notebook host and port. If you are running on a different location, you will need to supply this value for the application to display properly. If no protocol is supplied in the URL, e.g. if it is of the form “localhost:8888”, then “http” will be used.
Minimum period to assess the BLS to. If None, default value of 0.3 days will be used.
Maximum period to evaluate the BLS to. If None, the time coverage of the lightcurve / 4 will be used.
Number of points to use in the BLS panel. Lower this value for faster but less accurate performance. You can also vary this value using the widget’s Resolution Slider.
References
Examples
Load the light curve for Kepler-10, remove long-term trends, and display the BLS tool as follows:
>>> import lightkurve as lk
>>> lc = lk.search_lightcurvefile('kepler-10', quarter=3).download()
>>> lc = lc.PDCSAP_FLUX.normalize().flatten()
>>> lc.interact_bls()
normalize
(self, unit='unscaled')¶Returns a normalized version of the light curve.
The normalized light curve is obtained by dividing the flux
and
flux_err
object attributes by the by the median flux.
Optionally, the result will be multiplied by 1e2 (if unit='percent'
),
1e3 (unit='ppt'
), or 1e6 (unit='ppm'
).
The desired relative units of the normalized light curve; ‘ppt’ means ‘parts per thousand’, ‘ppm’ means ‘parts per million’.
LightCurve
A new light curve object in which flux
and flux_err
have
been divided by the median flux.
If the median flux is negative or within half a standard deviation from zero.
Examples
>>> import lightkurve as lk
>>> lc = lk.LightCurve(time=[1, 2, 3], flux=[25945.7, 25901.5, 25931.2], flux_err=[6.8, 4.6, 6.2])
>>> normalized_lc = lc.normalize()
>>> normalized_lc.flux
array([1.00055917, 0.99885466, 1. ])
>>> normalized_lc.flux_err
array([0.00026223, 0.00017739, 0.00023909])
plot
(self, **kwargs)¶Plot the light curve using Matplotlib’s plot
method.
Axes
A matplotlib axes object to plot into. If no axes is provided, a new one will be created.
Normalize the lightcurve before plotting?
Plot x axis label
Plot y axis label
Plot set_title
Path or URL to a matplotlib style file, or name of one of matplotlib’s built-in stylesheets (e.g. ‘ggplot’). Lightkurve’s custom stylesheet is used by default.
Dictionary of arguments to be passed to matplotlib.pyplot.plot
.
Axes
The matplotlib axes object.
remove_nans
(self)¶Removes cadences where the flux is NaN.
LightCurve
A new light curve object from which NaNs fluxes have been removed.
remove_outliers
(self, sigma=5.0, sigma_lower=None, sigma_upper=None, return_mask=False, **kwargs)¶Removes outlier data points using sigma-clipping.
This method returns a new LightCurve
object from which data points
are removed if their flux values are greater or smaller than the median
flux by at least sigma
times the standard deviation.
Sigma-clipping works by iterating over data points, each time rejecting values that are discrepant by more than a specified number of standard deviations from a center value. If the data contains invalid values (NaNs or infs), they are automatically masked before performing the sigma clipping.
Note
This function is a convenience wrapper around
astropy.stats.sigma_clip()
and provides the same functionality.
Any extra arguments passed to this method will be passed on to
sigma_clip
.
The number of standard deviations to use for both the lower and
upper clipping limit. These limits are overridden by
sigma_lower
and sigma_upper
, if input. Defaults to 5.
The number of standard deviations to use as the lower bound for
the clipping limit. Can be set to float(‘inf’) in order to avoid
clipping outliers below the median at all. If None
then the
value of sigma
is used. Defaults to None
.
The number of standard deviations to use as the upper bound for
the clipping limit. Can be set to float(‘inf’) in order to avoid
clipping outliers above the median at all. If None
then the
value of sigma
is used. Defaults to None
.
Whether or not to return a mask (i.e. a boolean array) indicating
which data points were removed. Entries marked as True
in the
mask are considered outliers. This mask is not returned by default.
Dictionary of arguments to be passed to astropy.stats.sigma_clip
.
LightCurve
A new light curve object from which outlier data points have been removed.
Boolean array flagging which cadences were removed.
Only returned if return_mask=True
.
Examples
This example generates a new light curve in which all points that are more than 1 standard deviation from the median are removed:
>>> lc = LightCurve(time=[1, 2, 3, 4, 5], flux=[1, 1000, 1, -1000, 1])
>>> lc_clean = lc.remove_outliers(sigma=1)
>>> lc_clean.time
array([1, 3, 5])
>>> lc_clean.flux
array([1, 1, 1])
Instead of specifying sigma
, you may specify separate sigma_lower
and sigma_upper
parameters to remove only outliers above or below
the median. For example:
>>> lc = LightCurve(time=[1, 2, 3, 4, 5], flux=[1, 1000, 1, -1000, 1])
>>> lc_clean = lc.remove_outliers(sigma_lower=float('inf'), sigma_upper=1)
>>> lc_clean.time
array([1, 3, 4, 5])
>>> lc_clean.flux
array([ 1, 1, -1000, 1])
Optionally, you may use the return_mask
parameter to return a boolean
array which flags the outliers identified by the method. For example:
>>> lc_clean, mask = lc.remove_outliers(sigma=1, return_mask=True)
>>> mask
array([False, True, False, True, False])
scatter
(self, colorbar_label='', show_colorbar=True, **kwargs)¶Plots the light curve using Matplotlib’s scatter
method.
Axes
A matplotlib axes object to plot into. If no axes is provided, a new one will be generated.
Normalize the lightcurve before plotting?
Plot x axis label
Plot y axis label
Plot set_title
Path or URL to a matplotlib style file, or name of one of matplotlib’s built-in stylesheets (e.g. ‘ggplot’). Lightkurve’s custom stylesheet is used by default.
Label to show next to the colorbar (if c
is given).
Show the colorbar if colors are given using the c
argument?
Dictionary of arguments to be passed to matplotlib.pyplot.scatter
.
Axes
The matplotlib axes object.
show_properties
(self)¶Prints a description of all non-callable attributes.
Prints in order of type (ints, strings, lists, arrays, others).
to_corrector
(self, method='sff')¶Returns a corrector object to remove instrument systematics.
Currently, only “sff” is supported. This will return a
SFFCorrector
class instance.
lightkurve.Correcter
Instance of a Corrector class, which typically provides correct()
and diagnose()
methods.
to_csv
(self, path_or_buf=None, **kwargs)¶Writes the light curve to a CSV file.
This method will convert the light curve into the Comma-Separated Values
(CSV) text format. By default this method will return the result as a
string, but you can also write the string directly to disk by providing
a file name or handle via the path_or_buf
parameter.
File path or object. By default, the result is returned as a string.
Dictionary of arguments to be passed to pandas.DataFrame.to_csv()
.
Returns a csv-formatted string if path_or_buf=None
.
Returns None
otherwise.
to_fits
(self, path=None, overwrite=False, flux_column_name='FLUX', aperture_mask=None, **extra_data)¶Writes the KeplerLightCurve to a FITS file.
File path, if None
returns an astropy.io.fits.HDUList object.
Whether or not to overwrite the file
The name of the label for the FITS extension, e.g. SAP_FLUX or FLUX
Optional 2D aperture mask to save with this lightcurve object, if defined. The mask can be either a boolean mask or an integer mask mimicking the Kepler/TESS convention; boolean masks are automatically converted to the Kepler/TESS conventions
Extra keywords or columns to include in the FITS file. Arguments of type str, int, float, or bool will be stored as keywords in the primary header. Arguments of type np.array or list will be stored as columns in the first extension.
Returns an astropy.io.fits object if path is None
to_pandas
(self, columns=('time', 'flux', 'flux_err'))¶Converts the light curve to a Pandas DataFrame
object.
By default, the object returned will contain the columns ‘time’, ‘flux’,
and ‘flux_err’. This can be changed using the columns
parameter.
List of columns to include in the DataFrame. The names must match
attributes of the LightCurve
object (e.g. time
, flux
).
pandas.DataFrame
A data frame indexed by time
and containing the columns flux
and flux_err
.
to_periodogram
(self, method='lombscargle', **kwargs)¶Converts the light curve to a Periodogram
power spectrum object.
This method will call either
lightkurve.periodogram.LombScarglePeriodogram.from_lightcurve()
or
lightkurve.periodogram.BoxLeastSquaresPeriodogram.from_lightcurve()
,
which in turn wrap astropy.stats.LombScargle
and astropy.stats.BoxLeastSquares
.
Optional keywords accepted if method='lombscargle'
are:
minimum_frequency
, maximum_frequency
, mininum_period
,
maximum_period
, frequency
, period
, nterms
,
nyquist_factor
, oversample_factor
, freq_unit
,
normalization
, ls_method
.
Optional keywords accepted if method='bls'
are
minimum_period
, maximum_period
, period
,
frequency_factor
, duration
.
Use the Lomb Scargle or Box Least Squares (BLS) method to
extract the power spectrum. Defaults to 'lombscargle'
.
'ls'
and 'bls'
are shorthands for 'lombscargle'
and 'boxleastsquares'
.
Keyword arguments passed to either
LombScarglePeriodogram
or
BoxLeastSquaresPeriodogram
.
Periodogram
objectThe power spectrum object extracted from the light curve.
to_seismology
(self, **kwargs)¶Returns a Seismology
object for estimating
quick-look asteroseismic quantities.
All **kwargs
will be passed to the to_periodogram()
method.
Seismology
objectObject which can be used to estimate quick-look asteroseismic quantities.
to_stingray
(self)¶Returns a stingray.Lightcurve
object.
This feature requires Stingray
to be installed (e.g. pip install stingray
). An ImportError
will
be raised if this package is not available.
stingray.Lightcurve
An stingray Lightcurve object.
to_table
(self)¶Converts the light curve to an Astropy Table
object.
astropy.table.Table
An AstroPy Table with columns ‘time’, ‘flux’, and ‘flux_err’.
to_timeseries
(self)¶Converts the light curve to an AstroPy
TimeSeries
object.
This feature requires AstroPy v3.2 or later (released in 2019).
An ImportError
will be raised if this version is not available.
TimeSeries
An AstroPy TimeSeries object.