<jupyter_start><jupyter_text>This tutorial continues part one, and recovers planets *c* and *d* in the K2-3 system. Planet *d* is a Super-Earth in the habitable zone. First, we recover planet *b* again:<jupyter_code>import numpy
import scipy
import everest
from astropy.stats import sigma_clip
import matplotlib.pyplot as plt
import logging
logging.basicConfig(level=logging.DEBUG)
logger = logging.getLogger('matplotlib')
logger.setLevel(logging.CRITICAL)
from matplotlib import rcParams; rcParams["figure.dpi"] = 150
EPIC_id = 201367065        
star = everest.Everest(EPIC_id)
t = numpy.delete(star.time, star.badmask)
y = numpy.delete(star.fcor, star.badmask)
t = numpy.array(t[~numpy.isnan(y)], dtype='float32')
y = numpy.array(y[~numpy.isnan(y)], dtype='float32')
trend = scipy.signal.medfilt(y, 25)
y_filt = y /trend
y_filt = sigma_clip(y_filt, sigma_upper=2, sigma_lower=float('inf'))
from transitleastsquares import transitleastsquares
model = transitleastsquares(t, y_filt)
results = model.power()
plt.figure()
plt.plot(
    results.model_folded_phase,
    results.model_folded_model,
    color='red')
plt.scatter(
    results.folded_phase,
    results.folded_y,
    color='blue',
    s=10,
    alpha=0.5,
    zorder=2)
plt.xlim(0.49, 0.51)
plt.xlabel('Phase')
plt.ylabel('Relative flux');<jupyter_output>/home/michael/.local/lib/python3.6/site-packages/pyfits/__init__.py:22: PyFITSDeprecationWarning: PyFITS is deprecated, please use astropy.io.fits
  PyFITSDeprecationWarning)  # noqa
/usr/lib/python3.6/importlib/_bootstrap.py:219: RuntimeWarning: numpy.dtype size changed, may indicate binary incompatibility. Expected 96, got 88
  return f(*args, **kwds)<jupyter_text>Now, we use the ``transit_mask`` function to mask the in-transit points of planet *b*. There are apparently more transits of another planet:<jupyter_code>from transitleastsquares import transit_mask, cleaned_array
intransit = transit_mask(t, results.period, 2*results.duration, results.T0)
y_second_run = y_filt[~intransit]
t_second_run = t[~intransit]
t_second_run, y_second_run = cleaned_array(t_second_run, y_second_run)
from matplotlib import rcParams; rcParams["figure.dpi"] = 150
plt.figure()
ax = plt.gca()
ax.plot(t_second_run, y_second_run, "k")
ax.set_xlim(t.min(), t.max())
ax.set_xlabel("Time (days)")
ax.set_ylabel("Normalized flux");<jupyter_output><empty_output><jupyter_text>The second run and resulting power spectrum:<jupyter_code>model_second_run = transitleastsquares(t_second_run, y_second_run)
results_second_run = model_second_run.power()
plt.figure()
ax = plt.gca()
ax.axvline(results_second_run.period, alpha=0.4, lw=3)
plt.xlim(numpy.min(results_second_run.periods), numpy.max(results_second_run.periods))
for n in range(2, 10):
    ax.axvline(n*results_second_run.period, alpha=0.4, lw=1, linestyle="dashed")
    ax.axvline(results_second_run.period / n, alpha=0.4, lw=1, linestyle="dashed")
plt.ylabel(r'SDE')
plt.xlabel('Period (days)')
plt.plot(results_second_run.periods, results_second_run.power, color='black', lw=0.5)
plt.xlim(0, max(results_second_run.periods));
plt.figure()
plt.plot(
    results_second_run.model_folded_phase,
    results_second_run.model_folded_model,
    color='red')
plt.scatter(
    results_second_run.folded_phase,
    results_second_run.folded_y,
    color='blue',
    s=10,
    alpha=0.5,
    zorder=2)
plt.xlim(0.49, 0.51)
plt.ylim(0.9985, 1.0003)
plt.xlabel('Phase')
plt.ylabel('Relative flux');

plt.figure()
in_transit = transit_mask(
    t_second_run,
    results_second_run.period,
    results_second_run.duration,
    results_second_run.T0)
plt.scatter(
    t_second_run[in_transit],
    y_second_run[in_transit],
    color='red',
    s=2,
    zorder=0)
plt.scatter(
    t_second_run[~in_transit],
    y_second_run[~in_transit],
    color='blue',
    alpha=0.5,
    s=2,
    zorder=0)
plt.plot(
    results_second_run.model_lightcurve_time,
    results_second_run.model_lightcurve_model, alpha=0.5, color='red', zorder=1)
plt.xlim(min(t), max(t))
plt.ylim(min(y_filt), max(y_filt))
plt.xlabel('Time (days)')
plt.ylabel('Relative flux');<jupyter_output><empty_output><jupyter_text>Looks like another planet (*c*). Again, we mask *c* and inspect the time series:<jupyter_code>intransit = transit_mask(t_second_run, results_second_run.period, 2*results_second_run.duration, results_second_run.T0)

y_third_run = y_second_run[~intransit]
t_third_run = t_second_run[~intransit]
t_third_run, y_third_run = cleaned_array(t_third_run, y_third_run)


from matplotlib import rcParams; rcParams["figure.dpi"] = 150
plt.figure()
ax = plt.gca()
ax.plot(t_third_run, y_third_run, "k")
ax.set_xlim(t.min(), t.max())
ax.set_xlabel("Time (days)")
ax.set_ylabel("Normalized flux");

model_third_run = transitleastsquares(t_third_run, y_third_run)
results_third_run = model_third_run.power()
plt.figure()
ax = plt.gca()
ax.axvline(results_third_run.period, alpha=0.4, lw=3)
plt.xlim(numpy.min(results_third_run.periods), numpy.max(results_third_run.periods))
for n in range(2, 10):
    ax.axvline(n*results_third_run.period, alpha=0.4, lw=1, linestyle="dashed")
    ax.axvline(results_third_run.period / n, alpha=0.4, lw=1, linestyle="dashed")
plt.ylabel(r'SDE')
plt.xlabel('Period (days)')
plt.plot(results_third_run.periods, results_third_run.power, color='black', lw=0.5)
plt.xlim(0, max(results_third_run.periods));
plt.figure()
plt.plot(
    results_third_run.model_folded_phase,
    results_third_run.model_folded_model,
    color='red')
plt.scatter(
    results_third_run.folded_phase,
    results_third_run.folded_y,
    color='blue',
    s=10,
    alpha=0.5,
    zorder=2)
plt.xlim(0.49, 0.51)
plt.ylim(0.9985, 1.0003)
plt.xlabel('Phase')
plt.ylabel('Relative flux');

plt.figure()
in_transit = transit_mask(
    t_third_run,
    results_third_run.period,
    results_third_run.duration,
    results_third_run.T0)
plt.scatter(
    t_third_run[in_transit],
    y_third_run[in_transit],
    color='red',
    s=2,
    zorder=0)
plt.scatter(
    t_third_run[~in_transit],
    y_third_run[~in_transit],
    color='blue',
    alpha=0.5,
    s=2,
    zorder=0)
plt.plot(
    results_third_run.model_lightcurve_time,
    results_third_run.model_lightcurve_model, alpha=0.5, color='red', zorder=1)
plt.xlim(min(t), max(t))
plt.ylim(min(y_filt), max(y_filt))
plt.xlabel('Time (days)')
plt.ylabel('Relative flux');<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.13821775964857233, 'fraction_numerical': 0.02306244116724192, 'mean_word_length': 5.832797427652733, 'pattern_counts': {'":': 0, '<': 16, '<?xml version=': 0, '>': 16, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '4', 'dir': 'jupyter-structured-clean-dedup', 'id': '514970', 'license': 'MIT', 'n_tokens_mistral': 2647, 'n_tokens_neox': 2523, 'n_words': 414, 'path': 'tutorials/03 Multiple planets - K2-3, a red dwarf with 3 Super-Earths.ipynb', 'repo_name': 'lgbouma/tls', 'size': '6374'}