Stage 30
Stage 30 fits either:
the white light curve generated in Stage 20
or the spectroscopic light curves generated in Stage 21
PACMAN automatically determines which previous stage products to use based on the settings in obs_par.pcf.
The relevant settings are:
s30_fit_white True/False
s30_fit_spec True/False
If s30_fit_white is True, PACMAN automatically uses the most recent:
stage20/s20_run_*
If s30_fit_spec is True, PACMAN automatically uses the most recent:
stage21/s21_run_*
A new Stage 30 run directory is then created.
The Stage 30 directory structure is:
stage30/
├── white_lc/
│ └── s30_run_YYYY-MM-DD_HH-MM-SS/
│ ├── extracted_lc/
│ └── fit_white/
│ ├── raw_lc/
│ ├── fit_lc/
│ ├── lsq_res/
│ ├── mcmc_res/
│ └── nested_res/
└── spec_lc/
└── s30_run_YYYY-MM-DD_HH-MM-SS/
├── extracted_sp/
└── fit_spec/
├── raw_lc/
├── fit_lc/
├── lsq_res/
├── mcmc_res/
└── nested_res/
At the beginning of the stage, the current obs_par.pcf and fit_par.txt
from pacman_run_files are copied into the new Stage 30 run directory.
The copied obs_par.pcf is then used to update the metadata for the run.
This ensures every Stage 30 run preserves the exact fitting settings used.
1) Preparation
Here we can fit the broadband (“white”) light curve (which was created in S20) or spectroscopic light curves (which were created in S21).
Let’s remove the first orbit from every visit and the first exposure from every orbit as they are typically strongly affected by instrument systematics:
We can choose if we also want to run an MCMC using the emcee package after running the least squares routine. Alternatively, the user can run a nested sampling routine using the dynesty package. Here, we run the least squares routine and then dynesty:
For the nested sampling, let’s use the dynamic approach:
Note: The user can also use an MCMC approach with emcee.
Let’s use the following model:
‘constant’: a normalization constant
‘upstream_downstream’: accounts for the forward and reverse scanning effect which creates an offset in the measured flux
‘model_ramp’: a ramp for every orbit
‘polynomial1’: a slope over the visit
‘transit’: a BATMAN transit model
‘uncmulti’: a free parameter which rescales the uncertainties at every step of the sampler
Additional functions are listed on the models page.
Let’s have the following free parameters:
t0_0, rp_0, u1_0, c_0, c_1, v_0, v_1, r1_0, r1_1, r2_0, r2_1, scale_0, scale_1, uncmulti_val_0, uncmulti_val_1
Other important parameters (per, ars, inc) are fixed to the literature values.
The user can set in the pcf whether the uncertainties should be rescaled to achieve a reduced chi2 of unity, using ‘rescale_uncert’. An alternative method which we are using here is using an additional free parameter which rescales the uncertainties at every step of the sampler. This model is called uncmulti.
2) Run PACMAN: White light curve fit
In the obs_par.pcf file, we set:
s30_fit_white True
s30_fit_spec False
Here’s the fit_par.txt file which was used in this example to fit the white light curve:
#parameter fixed tied value prior p1 p2
per True -1 1.5804046 X 2.25314 2e-05
t0 False -1 0.18 U 0.18 0.19
t_secondary True -1 0.0 X 0.0 0.0
w True -1 90.0 X 0.0 0.0
a True -1 15.23 X 4.98 0.05
inc True -1 89.1 X 85.3 0.2
rp False -1 0.116 U 0.01 0.2
fp True -1 0.0 X 0.0 1.0
u1 False -1 0.29 U 0.0 1.0
u2 True -1 0.0 X 0.0 0.0
ecc True -1 0.0 X 0.0 0.0
c False 0 8.45 U 8.2 8.5
c False 1 8.45 U 8.2 8.5
v False 0 0 U -1e-05 1e-05
v False 1 0 U -1e-05 1e-05
v2 True -1 0.0 X 0.0 0.0
r1 False 0 0.1 U 0 1.0
r1 False 1 0.1 U 0 1.0
r2 False 0 0.1 U 0.0 0.1
r2 False 1 0.1 U 0.0 0.1
r3 True -1 0.0 U -10.0 10.0
scale False 0 0.0 U -0.1 0.1
scale False 1 0.0 U -0.1 0.1
uncmulti_val False 0 2.0 U 0.1 8.0
uncmulti_val False 1 2.0 U 0.1 8.0
If a parameter is free and not jointly shared across visits, the user has to repeat the parameter in the fit_par.txt file as shown above (this might be fixed in a later version as it is not ideal, sepacially when the user wants to analyze many visits). Furthermore, the visit number must be added in the “tied” column. tied = -1 means that the parameters is shared across all visits, which makes sense for orbital parameters but less for systematics like the normalization constant.
Let’s run the white light curve fit now:
Starting s30
Using Stage 20 input directory: /Users/sebastianzieba/Desktop/Projects/Observations/Hubble/GJ1214_13021_2026/stage20/s20_run_2026-06-15_10-28-16
Location of the new Stage 30 run directory: /Users/sebastianzieba/Desktop/Projects/Observations/Hubble/GJ1214_13021_2026/stage30/white_lc/s30_run_2026-06-15_10-42-38
White light curve fit will be performed
Identified file(s) for fitting: [PosixPath('/Users/sebastianzieba/Desktop/Projects/Observations/Hubble/GJ1214_13021_2026/stage30/white_lc/s30_run_2026-06-15_10-42-38/extracted_lc/lc_white.txt')]
****** File: 1/1
Removed 8 exposures because they were the first exposures in the orbit.
Removed 34 exposures because they were the first orbit in the visit.
Current wavelength: 1.4 microns
Read in fit_par.pcf file:
parameter fixed tied value prior p1 p2
------------ ----- ---- --------- ----- ------- -----
per True -1 1.5804046 X 2.25314 2e-05
t0 False -1 0.18 U 0.18 0.19
t_secondary True -1 0.0 X 0.0 0.0
w True -1 90.0 X 0.0 0.0
a True -1 15.23 X 4.98 0.05
inc True -1 89.1 X 85.3 0.2
rp False -1 0.116 U 0.01 0.2
fp True -1 0.0 X 0.0 1.0
u1 False -1 0.29 U 0.0 1.0
u2 True -1 0.0 X 0.0 0.0
ecc True -1 0.0 X 0.0 0.0
c False 0 8.45 U 8.2 8.5
c False 1 8.45 U 8.2 8.5
v False 0 0.0 U -1e-05 1e-05
v False 1 0.0 U -1e-05 1e-05
v2 True -1 0.0 X 0.0 0.0
r1 False 0 0.1 U 0.0 1.0
r1 False 1 0.1 U 0.0 1.0
r2 False 0 0.1 U 0.0 0.1
r2 False 1 0.1 U 0.0 0.1
r3 True -1 0.0 U -10.0 10.0
scale False 0 0.0 U -0.1 0.1
scale False 1 0.0 U -0.1 0.1
uncmulti_val False 0 2.0 U 0.1 8.0
uncmulti_val False 1 2.0 U 0.1 8.0
Median log10 raw flux of full light curve: 8.431066421116899
The highest amount of exposures in an orbit is 18
Number of free parameters: 15
Names of free parameters: ['t0', 'rp', 'u1', 'c', 'c', 'v', 'v', 'r1', 'r1', 'r2', 'r2', 'scale', 'scale', 'uncmulti_val', 'uncmulti_val']
c sanity check passed:
Visit 0: log10(flux) mean=8.430215, std=0.002850, sanity range=[8.139674, 8.718145], c_0=8.450000
Visit 1: log10(flux) mean=8.430278, std=0.002774, sanity range=[8.146780, 8.710389], c_1=8.450000
The predicted rms is 63.69 ppm
*STARTS LEAST SQUARED*
Runs MPFIT...
t0_0 1.82781881e-01 9.61920459e-06
rp_0 1.15840824e-01 8.24474371e-05
u1_0 2.70979110e-01 5.00735301e-03
c_0 8.43146869e+00 4.38106609e-05
c_1 8.43126990e+00 3.82298955e-05
v_0 -1.23012909e-06 1.10621165e-07
v_1 -1.58628863e-07 1.10685866e-07
r1_0 3.75542137e-02 4.21618627e-03
r1_1 4.04501573e-02 3.95761217e-03
r2_0 1.71143133e-03 8.87915185e-05
r2_1 1.78757967e-03 7.64847900e-05
scale_0 4.16736148e-03 1.75043899e-05
scale_1 4.18389841e-03 1.75006456e-05
rms, chi2red = 117.32411228761931 3.9407785958955266
Saved white_systematics.txt file
*STARTS NESTED SAMPLING*
Using multiprocessing...
Run dynesty...
1it [00:01, 1.54s/it, batch: 0 | bound: 0 | nc: 1 | ncall: 1 | eff(%): 0.498 | loglstar: -inf < -inf < inf | lo
4457it [00:11, 345.85it/s, batch: 0 | bound: 44 | nc: 35 | ncall: 137556 | eff(%): 3.235 | loglstar: -inf < -5009.859 < inf | logz: -5037.519 +/- 0.365 | dlogz: 26
25487it [01:10, 359.00it/s, batch: 6 | bound: 2 | nc: 1 | ncall: 857511 | eff(%): 2.926 | loglstar: -1284.668 < -1277.056 < -1283.932 | logz: -1353.017 +/- 0.344 | stop: 0.967]
Saved white_systematics.txt file for nested sampling run
Saved fit_results.txt file
Saving Metadata
Finished s30
There are several plots created then:
The raw light curve:
** From the least squares routine **
The fitted light curve without the systematics:
The time averaging plot (formally known as Allan deviation plot):
** Using dynesty **
The fitted light curve without the systematics:
Corner plot from the nested sampling:
** Using emcee **
MCMC chains with burn-in:
MCMC chains without burn-in
3) Run PACMAN: Spectroscopic light curve fit
In the obs_par.pcf file, we set:
s30_fit_white False
s30_fit_spec True
Here’s the fit_par.txt file which was used in this example to fit the spectroscopic light curves:
#parameter fixed tied value prior p1 p2
per True -1 1.5804046 X 2.25314 2e-05
t0 False -1 0.18 U 0.18 0.19
t_secondary True -1 0.0 X 0.0 0.0
w True -1 90.0 X 0.0 0.0
a True -1 15.23 X 4.98 0.05
inc True -1 89.1 X 85.3 0.2
rp False -1 0.116 U 0.01 0.2
fp True -1 0.0 X 0.0 1.0
u1 False -1 0.29 U 0.0 1.0
u2 True -1 0.0 X 0.0 0.0
ecc True -1 0.0 X 0.0 0.0
c False 0 6.25 U 6.1 6.4
c False 1 6.25 U 6.1 6.4
v False 0 0 U -1e-05 1e-05
v False 1 0 U -1e-05 1e-05
v2 True -1 0.0 X 0.0 0.0
r1 False 0 0.1 U 0 1.0
r1 False 1 0.1 U 0 1.0
r2 False 0 0.1 U 0.0 0.1
r2 False 1 0.1 U 0.0 0.1
r3 True -1 0.0 U -10.0 10.0
scale False 0 0.0 U -0.1 0.1
scale False 1 0.0 U -0.1 0.1
uncmulti_val False 0 2.0 U 0.1 8.0
uncmulti_val False 1 2.0 U 0.1 8.0
It is worth noting that PACMAN returns and error if it detects that the constants in the fit_par.txt file are not within a sanity check range which is determined from the data. This is to prevent the user from running a fit with unphysical parameters which can lead to very long runtimes and/or non-convergence of the fit. This might lead a raised error message like this one:
Suspicious c values in fit_par.txt:
- Visit 0, c_0: initial value 7.450000 is outside the log10(flux) sanity range [5.954690, 6.535117]. Flux stats: mean=1.763025e+06, min=1.740295e+06, max=1.774943e+06; log10(flux) stats: mean=6.246249, std=0.002859.
- Visit 0, c_0: lower prior bound 7.200000 is outside the log10(flux) sanity range [5.954690, 6.535117]. Flux stats: mean=1.763025e+06, min=1.740295e+06, max=1.774943e+06; log10(flux) stats: mean=6.246249, std=0.002859.
- Visit 0, c_0: upper prior bound 7.500000 is outside the log10(flux) sanity range [5.954690, 6.535117]. Flux stats: mean=1.763025e+06, min=1.740295e+06, max=1.774943e+06; log10(flux) stats: mean=6.246249, std=0.002859.
If you have set the fit_par.txt correctly and run the fit, you should get an output similar to this one:
Most plots which are created during the white light curve fit will be also created during and after running the spectroscopic fits. Let’s look at some examples:
The first raw spectroscopic light curve:
** Using least squared **
The fitted spectroscopic light curve without the systematics:
All fitted parameters as a function of wavelength:
** Using dynesty **
The spectrum (rprs vs wavelength):
There is also a fit that stores the transmission spectrum:
#wavelength wave_width wave_low wave_high p50 p16 p84 minus plus maxL
1.15804545e+00 2.30454545e-02 1.13500000e+00 1.18109091e+00 1.16251726e-01 1.15861551e-01 1.16636190e-01 3.90175124e-04 3.84464096e-04 1.16215265e-01
1.20413636e+00 2.30454545e-02 1.18109091e+00 1.22718182e+00 1.16163257e-01 1.15805905e-01 1.16541892e-01 3.57351848e-04 3.78635596e-04 1.15974567e-01
1.25022727e+00 2.30454545e-02 1.22718182e+00 1.27327273e+00 1.15685002e-01 1.15398621e-01 1.15974707e-01 2.86381423e-04 2.89705342e-04 1.15759608e-01
1.29631818e+00 2.30454545e-02 1.27327273e+00 1.31936364e+00 1.15931077e-01 1.15628220e-01 1.16242732e-01 3.02857838e-04 3.11654603e-04 1.15983434e-01
1.34240909e+00 2.30454545e-02 1.31936364e+00 1.36545455e+00 1.15806511e-01 1.15494810e-01 1.16120327e-01 3.11701242e-04 3.13816180e-04 1.15835558e-01
1.38850000e+00 2.30454545e-02 1.36545455e+00 1.41154545e+00 1.16415949e-01 1.16092146e-01 1.16748365e-01 3.23802236e-04 3.32416003e-04 1.16524251e-01
1.43459091e+00 2.30454545e-02 1.41154545e+00 1.45763636e+00 1.15854876e-01 1.15514348e-01 1.16189964e-01 3.40527938e-04 3.35087435e-04 1.15893762e-01
1.48068182e+00 2.30454545e-02 1.45763636e+00 1.50372727e+00 1.15878373e-01 1.15566288e-01 1.16205192e-01 3.12085400e-04 3.26819316e-04 1.15864558e-01
1.52677273e+00 2.30454545e-02 1.50372727e+00 1.54981818e+00 1.15953617e-01 1.15664661e-01 1.16235760e-01 2.88956634e-04 2.82142311e-04 1.15863863e-01
1.57286364e+00 2.30454545e-02 1.54981818e+00 1.59590909e+00 1.16102938e-01 1.15794136e-01 1.16409748e-01 3.08802494e-04 3.06809688e-04 1.16073204e-01
1.61895455e+00 2.30454545e-02 1.59590909e+00 1.64200000e+00 1.16158327e-01 1.15821426e-01 1.16505991e-01 3.36901548e-04 3.47663982e-04 1.15937821e-01
And finally a file that summarizes the fit results for the light curves: .. include:: media/s30/spectroscopic/fit_results.txt
- literal: