diff --git a/changelog.md b/changelog.md index ee2ba62..4e4e442 100644 --- a/changelog.md +++ b/changelog.md @@ -1,6 +1,11 @@ # Changelog All notable changes to this project will be documented in this file. +### [1.3.1] -- 2024-01-19 +#### Added +- Fixed implementation of STScI up-the-ramp JumpStep which would not run consistently. +- Made it easier to select stellar models to use for limb-darkening calculation in Stage 4. + ### [1.3.0] -- 2023-12-21 #### Added - Compatibility with v1.12.5 of the STScI jwst package. diff --git a/setup.py b/setup.py index f23d9ea..a2ff2a0 100644 --- a/setup.py +++ b/setup.py @@ -4,7 +4,7 @@ from setuptools import setup setup(name='supreme_spoon', - version='1.3.0', + version='1.3.1', license='MIT', author='Michael Radica', author_email='michael.radica@umontreal.ca', diff --git a/supreme_spoon/fit_lightcurves.py b/supreme_spoon/fit_lightcurves.py index 12a8acd..40468d5 100644 --- a/supreme_spoon/fit_lightcurves.py +++ b/supreme_spoon/fit_lightcurves.py @@ -207,7 +207,8 @@ assert np.all(np.array([m_h, logg, teff]) != None), msg c1, c2 = stage4.gen_ld_coefs(config['spectrace_ref'], wave_low, wave_up, order, m_h, logg, teff, - config['ld_data_path']) + config['ld_data_path'], + model_type=config['ld_model_type']) q1, q2 = juliet.reverse_q_coeffs('quadratic', c1, c2) # Save calculated coefficients. target = fits.getheader(config['infile'], 0)['TARGET'] diff --git a/supreme_spoon/fit_lightcurves.yaml b/supreme_spoon/fit_lightcurves.yaml index a2bd41a..44637d1 100644 --- a/supreme_spoon/fit_lightcurves.yaml +++ b/supreme_spoon/fit_lightcurves.yaml @@ -63,6 +63,9 @@ teff : None # Path to ExoTiC-LD Reference Files # See ExoTiC-LD documentation for more info: https://exotic-ld.readthedocs.io/en/latest/views/installation.html. ld_data_path : +# Stellar models to use for limb-darkening calculation. +# See ExoTiC-LD documentation for more info. +ld_model_type : 'stagger' # Path to JWST spectrace Reference File # Will be in crds_cache with file name like jwst_niriss_spectrace_XXXX.fits spectrace_ref : './crds_cache/references/jwst/niriss/jwst_niriss_spectrace_0023.fits' diff --git a/supreme_spoon/stage1.py b/supreme_spoon/stage1.py index 965d27f..d20525f 100644 --- a/supreme_spoon/stage1.py +++ b/supreme_spoon/stage1.py @@ -355,8 +355,15 @@ def run(self, save_results=True, force_redo=False, do_plot=False, '.fits', '_2.pdf') else: plot_file1, plot_file2 = None, None + + # For scale-chromatic correction, collapse the 2D timeseries to + # 1D for plotting purposes. + if self.method == 'scale-chromatic': + this_ts = np.nanmedian(self.timeseries, axis=1) + else: + this_ts = self.timeseries plotting.make_oneoverf_plot(results, - timeseries=self.timeseries, + timeseries=this_ts, baseline_ints=self.baseline_ints, outfile=plot_file1, show_plot=show_plot) @@ -366,7 +373,7 @@ def run(self, save_results=True, force_redo=False, do_plot=False, else: window = False plotting.make_oneoverf_psd(results, self.datafiles, - timeseries=self.timeseries, + timeseries=this_ts, baseline_ints=self.baseline_ints, pixel_masks=self.pixel_masks, outfile=plot_file2, @@ -465,8 +472,14 @@ def run(self, save_results=True, force_redo=False, flag_up_ramp=False, # Get number of groups in the observation - ngroup=2 must be # treated in a special way as the default pipeline JumpStep # will fail. + # Also need to set minimum_sigclip_groups to something >nints, + # else the up-the-ramp jump detection will be replaced by a + # time-domain sigma clipping. testfile = datamodels.open(self.datafiles[0]) ngroups = testfile.meta.exposure.ngroups + i_start = testfile.meta.exposure.integration_start + i_end = testfile.meta.exposure.integration_end + nints = i_end - i_start testfile.close() # For ngroup > 2, default JumpStep can be used. if ngroups > 2 and flag_up_ramp is True: @@ -474,7 +487,9 @@ def run(self, save_results=True, force_redo=False, flag_up_ramp=False, res = step.call(segment, output_dir=self.output_dir, save_results=save_results, rejection_threshold=rejection_threshold, - maximum_cores='quarter', **kwargs) + maximum_cores='quarter', + minimum_sigclip_groups=nints+100, + **kwargs) # Time domain jump step must be run for ngroup=2. else: res = segment diff --git a/supreme_spoon/stage4.py b/supreme_spoon/stage4.py index 82fa9f9..a16a36f 100644 --- a/supreme_spoon/stage4.py +++ b/supreme_spoon/stage4.py @@ -344,8 +344,8 @@ def fit_lightcurves(data_dict, prior_dict, order, output_dir, fit_suffix, return results -def gen_ld_coefs(spectrace_ref, wavebin_low, wavebin_up, order, m_h, logg, teff, - ld_data_path, model_type='stagger'): +def gen_ld_coefs(spectrace_ref, wavebin_low, wavebin_up, order, m_h, logg, + teff, ld_data_path, model_type='stagger'): """Generate estimates of quadratic limb-darkening coefficients using the ExoTiC-LD package. diff --git a/supreme_spoon/utils.py b/supreme_spoon/utils.py index bb3556c..c465c5b 100644 --- a/supreme_spoon/utils.py +++ b/supreme_spoon/utils.py @@ -310,10 +310,10 @@ def get_filename_root(datafiles): # Now assuming everything is in chronological order, just increment the # segment number. split = fileroot.split('seg') - for segment in range(seg_start+1, len(datafiles)+1): + for segment in range(seg_start+1, seg_start+len(datafiles)): if segment < 10: seg_no = 'seg00{}'.format(segment) - elif 10 <= segment < 99: + elif 10 <= segment <= 99: seg_no = 'seg0{}'.format(segment) else: seg_no = 'seg{}'.format(segment) @@ -709,6 +709,7 @@ def outlier_resistant_variance(data): return var +# TODO: Include something to check that inputs are correct format/not nonsensical. def parse_config(config_file): """Parse a yaml config file. @@ -793,11 +794,11 @@ def save_extracted_spectra(filename, wl1, wu1, f1, e1, wl2, wu2, f2, e2, t, hdu3 = fits.ImageHDU(wu1, header=hdr) hdr = fits.Header() hdr['EXTNAME'] = "Flux O1" - hdr['UNITS'] = "Electrons" + hdr['UNITS'] = "DN/s" hdu4 = fits.ImageHDU(f1, header=hdr) hdr = fits.Header() hdr['EXTNAME'] = "Flux Err O1" - hdr['UNITS'] = "Electrons" + hdr['UNITS'] = "DN/s" hdu5 = fits.ImageHDU(e1, header=hdr) # Pack order 2 values. @@ -811,17 +812,17 @@ def save_extracted_spectra(filename, wl1, wu1, f1, e1, wl2, wu2, f2, e2, t, hdu7 = fits.ImageHDU(wu2, header=hdr) hdr = fits.Header() hdr['EXTNAME'] = "Flux O2" - hdr['UNITS'] = "Electrons" + hdr['UNITS'] = "DN/s" hdu8 = fits.ImageHDU(f2, header=hdr) hdr = fits.Header() hdr['EXTNAME'] = "Flux Err O2" - hdr['UNITS'] = "Electrons" + hdr['UNITS'] = "DN/s" hdu9 = fits.ImageHDU(e2, header=hdr) # Pack time axis. hdr = fits.Header() hdr['EXTNAME'] = "Time" - hdr['UNITS'] = "BJD" + hdr['UNITS'] = "MJD" hdu10 = fits.ImageHDU(t, header=hdr) if save_results is True: