A look at the shadow of HD 332231 b

Here we will be using tracit to look at the shadow of HD 332231 b


First we import tracit. The two run commands are just to setup text rendering in the plots to LaTeX.

[1]:
import tracit

We want to create two dictionaries:

  • par in which we specify the values for the parameters, priors, boundaries, etc.

  • dat in which we give the filenames for our data, and we can also specify the de-trending/noise model

[2]:
par = tracit.par_struct(n_planets=1,n_phot=0,n_spec=1)
dat = tracit.dat_struct(n_phot=0,n_rvs=0,n_ls=1,n_sl=0)

For the parameter structure we specified that we have 1 planet in the system, 1 photometric system, and 1 spectroscopic system (RVs, shadow, or slope).

For the data we also define the number of photometric systems, while we further specify the type of spectroscopic measurements\(-\)n_rvs for RVs, n_ls for the shadow (lineshape), or n_sl for the subplanetary velocity (slope).

Here we want to look at the shadow, so we set n_ls=1.

Below we set the path to our shadow data, and give the indices (counting from 0) for the out-of-transit timestamps.

We also specify a range in the CCF, where there’s no bump. In this case, for instance, the bump of the CCF does not extend out to beyond 12 km/s. It’s of course important to ensure that the value you give is also inside your velocity grid. This is used to de-trend and estimate the scatter of the CCF.

[3]:
dat['LS filename_1'] = 'shadow.hdf5' # filename
dat['idxs_1'] = [39,40,41] # The out-of-transit indices
dat['No_bump_1'] = 12 # No bump in the CCF, i.e., from (vel < - 12 km/s) or (vel > 12 km/s)
dat['Velocity_range_1'] = 20

We then initialize the data, and the business module.

[4]:
tracit.ini_data(dat)
tracit.run_bus(par,dat)

Here we import some results from a previous run to make a nice plot.

[5]:
saved_results = 1 # set to 1 if you want to read in the results instead
if saved_results:
    import pandas as pd
    rdf = pd.read_csv('results_from_old_fit.csv')
    tracit.update_pars(rdf,par,best_fit=False)

If you don’t have results from a previous run, you’ll have to set them individually as shown below.

[6]:
par['T0_b']['Value'] = 2458729.6940
par['lam_b']['Value'] = 0.0

Let’s take a look at the shadow of HD 332231 b. We do this by calling plot_shadow.

[7]:
tracit.run_exp(1)
tracit.plot_shadow(par,dat)
/home/emil/Desktop/PhD/tracit/src/tracit/shady.py:482: RuntimeWarning: divide by zero encountered in true_divide
  tan = np.exp(-1*np.power(vel_1d/(zeta*y),2))/y
/home/emil/Desktop/PhD/tracit/src/tracit/shady.py:482: RuntimeWarning: invalid value encountered in true_divide
  tan = np.exp(-1*np.power(vel_1d/(zeta*y),2))/y
Using indices [39, 40, 41] as out-of-transit spectra
../../_images/examples_shadow_hd332231_fit_14_2.png

We can also look at the out-of-transit CCFs. Here we’ll also look at all the in-transit CCFs after subtracting our average out-of-transit CCF, i.e., the 1D shadow.

[8]:
tracit.plot_oot_ccf(par,dat)
/home/emil/Desktop/PhD/tracit/src/tracit/shady.py:482: RuntimeWarning: divide by zero encountered in true_divide
  tan = np.exp(-1*np.power(vel_1d/(zeta*y),2))/y
/home/emil/Desktop/PhD/tracit/src/tracit/shady.py:482: RuntimeWarning: invalid value encountered in true_divide
  tan = np.exp(-1*np.power(vel_1d/(zeta*y),2))/y
Using indices [39, 40, 41] as out-of-transit spectra
## Spectroscopic system 1/Spectrograph\ 1 ##:

Reduced chi-squared for the oot CCF is:
         12.089
Factor to apply to get a reduced chi-squared around 1.0 is:
         3.477

Number of data points: 31
Number of fitting parameters: 0
#########################
/home/emil/Desktop/PhD/tracit/src/tracit/shady.py:482: RuntimeWarning: divide by zero encountered in true_divide
  tan = np.exp(-1*np.power(vel_1d/(zeta*y),2))/y
/home/emil/Desktop/PhD/tracit/src/tracit/shady.py:482: RuntimeWarning: invalid value encountered in true_divide
  tan = np.exp(-1*np.power(vel_1d/(zeta*y),2))/y
../../_images/examples_shadow_hd332231_fit_16_3.png
../../_images/examples_shadow_hd332231_fit_16_4.png

If you want to fit the shadow you can do as shown below.

[9]:
dat['Fit LS_1'] = 1 # whether to fit the shadow or not
par['lam_b']['Prior'] = 'uni'
par['lam_b']['Prior_vals'] = [90,2,-180,180]
par['FPs'] = ['lam_b']
mc = 0
if mc:
    ndraws = 1500
    nwalkers = 10
    nproc = 1
    ndf = tracit.mcmc(par,dat,ndraws,nwalkers,nproc=1,corner=True,chains=True,path_to_tracit=path)

You can also take a look at the slope. For that you also need to set the LC to true (for now at least).

[10]:
slope = 0
if slope:
    par = tracit.par_struct(n_planets=1,n_phot=1,n_spec=1)
    dat = tracit.dat_struct(n_phot=1,n_rvs=0,n_ls=0,n_sl=1)

    dat['LC filename_1'] = 'lc_HD332231.txt'

    import pandas as pd
    rdf = pd.read_csv('results_from_old_fit.csv')
    tracit.update_pars(rdf,par,best_fit=False)

    par['T0_b']['Value'] = 2458729.6940
    par['lam_b']['Value'] = 0.0

    dat['SLs'] = 1
    dat['SL filename_1'] = 'shadow.hdf5' # filename
    dat['idxs_1'] = [39,40,41] # The out-of-transit indices
    dat['No_bump_1'] = 12 # No bump in the CCF, i.e., from (vel < - 12 km/s) or (vel > 12 km/s)
    tracit.ini_data(dat)
    tracit.run_bus(par,dat)

    tracit.plot_slope(par,dat)