In [ ]:
import wobble
import numpy as np
import matplotlib.pyplot as plt
import glob
import os
from tqdm import tqdm
np.random.seed(0)
In [ ]:
data = wobble.Data()
In [ ]:
filenames = glob.glob('/Users/mbedell/python/wobble/data/toi/TOI-*_CCF_A.fits')
for filename in tqdm(filenames):
try:
sp = wobble.Spectrum()
sp.from_ESPRESSO(filename, process=True)
data.append(sp)
except Exception as e:
print("File {0} failed; error: {1}".format(filename, e))
In [ ]:
data.write('../data/toi.hdf5')
In [ ]:
data = wobble.Data(filename='../data/toi.hdf5')
R = np.copy(data.R) # we'll need this later
In [ ]:
data
In [ ]:
data.drop_bad_orders(min_snr=3)
data.drop_bad_epochs(min_snr=3)
In [ ]:
data.orders
In [ ]:
r = 0
good = data.ivars[r] > 0.
for e in [0,10,20]:
plt.errorbar(data.xs[r][e][good[e]], data.ys[r][e][good[e]],
1./np.sqrt(data.ivars[r][e][good[e]]), ls='', fmt='o', ms=2, alpha=0.5)
plt.title('Echelle order #{0}'.format(data.orders[r]), fontsize=14);
Since we don't have any existing regularization parameter files for ESPRESSO, we have to make some new ones.
This is needed because the default wobble regularization is tuned to HARPS, which has a different number of spectral orders and different wavelength coverage - if we try to run with those files, the optimization will (a) be non-optimal and (b) eventually crash when we try to access an order than does not exist for HARPS.
In [ ]:
star_filename = '../wobble/regularization/toi_star.hdf5'
tellurics_filename = '../wobble/regularization/toi_tellurics.hdf5'
wobble.generate_regularization_file(star_filename, R, type='star')
wobble.generate_regularization_file(tellurics_filename, R, type='telluric')
In [ ]:
plot_dir = '../regularization/toi/'
if not os.path.exists(plot_dir):
os.makedirs(plot_dir)
We'll tune the regularization using a train-and-validate approach, so let's set aside some epochs to be the validation set:
In [ ]:
validation_epochs = np.random.choice(data.N, data.N//6, replace=False) # 3 epochs for validation set
In [ ]:
r = 100
for e in [validation_epochs[0]]:
plt.errorbar(data.xs[r][e][good[e]], data.ys[r][e][good[e]],
1./np.sqrt(data.ivars[r][e][good[e]]), ls='', fmt='o', ms=2, alpha=0.5)
Here's an example of how this regularization tuning will go for one order:
In [ ]:
r = 100
o = data.orders[r]
objs = wobble.setup_for_order(r, data, validation_epochs)
In [ ]:
wobble.improve_order_regularization(o, star_filename, tellurics_filename,
*objs,
verbose=False, plot=False,
basename='{0}o{1}'.format(plot_dir, o),
K_t=0, L1=True, L2=True)
This only does one order. To save time & print statements, we'll do the full loop over all orders from a script. See wobble/scripts/tune_regularization.py for an example.
OK, let's assume that we have run the regularization tuning script and we're satisfied with the settings there. Now we can finally get some solutions!