Plot the spectrum of a transiting planet

This example runs VSPEC with a transiting planet scenario.

from pathlib import Path
import numpy as np
import matplotlib.pyplot as plt
import pypsg

from VSPEC import ObservationModel,PhaseAnalyzer

try:
    CFG_PATH = Path(__file__).parent / 'transit_spectrum.yaml'
except NameError:
    CFG_PATH = Path('transit_spectrum.yaml')

pypsg.docker.set_url_and_run()

Saved settings to /home/ted/.pypsg/settings.json
Reloading settings...

Initialize the VSPEC run

We read in the config file and run the model.

model = ObservationModel.from_yaml(CFG_PATH)
model.build_planet()
model.build_spectra()

/home/ted/github/VSPEC/VSPEC/gcm/heat_transfer.py:538: RuntimeWarning: Energy balance off by -9.6 %. eps = 6.0, T0 = 0.7
  warnings.warn(msg, RuntimeWarning)
Starting at phase 178.8 deg, observe for 1.5 h in 18 steps
Phases = [178.8  178.92 179.05 179.18 179.3  179.42 179.55 179.68 179.8  179.92
 180.05 180.17 180.3  180.42 180.55 180.67 180.8  180.92 181.05] deg

Build Planet:   0%|          | 0/18 [00:00<?, ?it/s]
Build Planet:   6%|▌         | 1/18 [00:00<00:06,  2.43it/s]
Build Planet:  11%|█         | 2/18 [00:00<00:05,  2.99it/s]
Build Planet:  17%|█▋        | 3/18 [00:00<00:04,  3.11it/s]
Build Planet:  22%|██▏       | 4/18 [00:01<00:04,  2.81it/s]
Build Planet:  28%|██▊       | 5/18 [00:01<00:04,  2.68it/s]
Build Planet:  33%|███▎      | 6/18 [00:02<00:04,  2.58it/s]
Build Planet:  39%|███▉      | 7/18 [00:02<00:04,  2.50it/s]
Build Planet:  44%|████▍     | 8/18 [00:03<00:04,  2.43it/s]
Build Planet:  50%|█████     | 9/18 [00:03<00:03,  2.44it/s]
Build Planet:  56%|█████▌    | 10/18 [00:03<00:03,  2.47it/s]
Build Planet:  61%|██████    | 11/18 [00:04<00:02,  2.47it/s]
Build Planet:  67%|██████▋   | 12/18 [00:04<00:02,  2.40it/s]
Build Planet:  72%|███████▏  | 13/18 [00:05<00:02,  2.41it/s]
Build Planet:  78%|███████▊  | 14/18 [00:05<00:01,  2.43it/s]
Build Planet:  83%|████████▎ | 15/18 [00:05<00:01,  2.41it/s]
Build Planet:  89%|████████▉ | 16/18 [00:06<00:00,  2.42it/s]
Build Planet:  94%|█████████▍| 17/18 [00:06<00:00,  2.39it/s]
Build Planet: 100%|██████████| 18/18 [00:07<00:00,  2.62it/s]
Build Planet: 19it [00:07,  2.75it/s]
Build Planet: 19it [00:07,  2.56it/s]

Build Spectra:   0%|          | 0/18 [00:00<?, ?it/s]

Loading Spectra:   0%|          | 0/17 [00:00<?, ?it/s]

Loading Spectra:   6%|▌         | 1/17 [00:00<00:04,  3.24it/s]

Loading Spectra:  12%|█▏        | 2/17 [00:00<00:04,  3.08it/s]

Loading Spectra:  18%|█▊        | 3/17 [00:00<00:04,  3.06it/s]

Loading Spectra:  24%|██▎       | 4/17 [00:01<00:04,  3.07it/s]

Loading Spectra:  29%|██▉       | 5/17 [00:01<00:03,  3.03it/s]

Loading Spectra:  35%|███▌      | 6/17 [00:01<00:03,  3.06it/s]

Loading Spectra:  41%|████      | 7/17 [00:02<00:03,  3.09it/s]

Loading Spectra:  47%|████▋     | 8/17 [00:02<00:02,  3.06it/s]

Loading Spectra:  53%|█████▎    | 9/17 [00:02<00:02,  3.11it/s]

Loading Spectra:  59%|█████▉    | 10/17 [00:03<00:02,  3.17it/s]

Loading Spectra:  65%|██████▍   | 11/17 [00:03<00:01,  3.22it/s]

Loading Spectra:  71%|███████   | 12/17 [00:03<00:01,  3.26it/s]

Loading Spectra:  76%|███████▋  | 13/17 [00:04<00:01,  3.20it/s]

Loading Spectra:  82%|████████▏ | 14/17 [00:04<00:00,  3.18it/s]

Loading Spectra:  88%|████████▊ | 15/17 [00:04<00:00,  3.18it/s]

Loading Spectra:  94%|█████████▍| 16/17 [00:05<00:00,  3.14it/s]

Loading Spectra: 100%|██████████| 17/17 [00:05<00:00,  3.10it/s]
Loading Spectra: 100%|██████████| 17/17 [00:05<00:00,  3.13it/s]

Build Spectra:   6%|▌         | 1/18 [02:30<42:43, 150.77s/it]
Build Spectra:  11%|█         | 2/18 [02:31<16:35, 62.23s/it]
Build Spectra:  17%|█▋        | 3/18 [02:31<08:28, 33.93s/it]
Build Spectra:  22%|██▏       | 4/18 [02:31<04:49, 20.66s/it]
Build Spectra:  28%|██▊       | 5/18 [02:31<02:53, 13.34s/it]
Build Spectra:  33%|███▎      | 6/18 [02:32<01:46,  8.90s/it]
Build Spectra:  39%|███▉      | 7/18 [02:32<01:06,  6.09s/it]
Build Spectra:  44%|████▍     | 8/18 [02:32<00:42,  4.25s/it]
Build Spectra:  50%|█████     | 9/18 [02:33<00:27,  3.01s/it]
Build Spectra:  56%|█████▌    | 10/18 [02:33<00:17,  2.17s/it]
Build Spectra:  61%|██████    | 11/18 [02:33<00:11,  1.59s/it]
Build Spectra:  67%|██████▋   | 12/18 [02:33<00:07,  1.20s/it]
Build Spectra:  72%|███████▏  | 13/18 [02:34<00:04,  1.08it/s]
Build Spectra:  78%|███████▊  | 14/18 [02:34<00:02,  1.36it/s]
Build Spectra:  83%|████████▎ | 15/18 [02:34<00:01,  1.63it/s]
Build Spectra:  89%|████████▉ | 16/18 [02:35<00:01,  1.88it/s]
Build Spectra:  94%|█████████▍| 17/18 [02:35<00:00,  2.17it/s]
Build Spectra: 100%|██████████| 18/18 [02:35<00:00,  2.60it/s]
Build Spectra: 100%|██████████| 18/18 [02:35<00:00,  8.65s/it]

Load in the data

We can use VSPEC to read in the synthetic data we just created.

data = PhaseAnalyzer(model.directories['all_model'])

Calculate the transit depth

Since this star model has no limb darkeing, no spots, and no noise, we don’t need to fit a model to our data to extract the transit depth.

cmb_data = data.total
continuum = cmb_data[:,0] # the first epoch
data_min = np.min(cmb_data,axis=1)
transit_depth = (continuum-data_min)/continuum
rp_rs = np.sqrt(transit_depth)
wavelength = data.wavelength

plt.plot(wavelength,rp_rs*100)
plt.xlabel(f'Wavelength {wavelength.unit}')
plt.ylabel(r'$\frac{R_p}{R_*}$ (%)')

Text(24.972222222222214, 0.5, '$\\frac{R_p}{R_*}$ (%)')