Source code for GUIBRUSHR.Retrieval.ModelCalculation.ModelSetup
"""
ModelSetup Module
This module contains the ModelSetup class, extracted from ModelData.py as part of
the codebase decomposition. It handles initialization and setup of core state objects
for atmospheric retrieval, including reading configuration files, extracting night
observation data, and setting up the atmosphere/retrieval/bestpars objects.
"""
import os
from pathlib import Path
import pandas as pd
from GUIBRUSHR.General_Constants.FunctionsAndConstants.Constant_Variables import ConstantVariables
from GUIBRUSHR.General_Constants.Classes.Instrument import Instrument
from GUIBRUSHR.General_Constants.Classes.Night import Night
from GUIBRUSHR.General_Constants.Classes.Target import Target
from GUIBRUSHR.Retrieval.ModelCalculation.Classes.Atmosphere import Atmosphere
from GUIBRUSHR.Retrieval.ModelCalculation.Classes.Retrieval import Retrieval
from GUIBRUSHR.Retrieval.ModelCalculation.Classes.Bestpars import Bestpars
from GUIBRUSHR.Retrieval.ModelCalculation.Classes.Random import Random
from GUIBRUSHR.General_Constants.FunctionsAndConstants.general_functions import get_csv_value
from GUIBRUSHR.core.io.retrieval_io import read_general_info_df
from GUIBRUSHR.core.io._coerce import coerce_bool as _coerce_bool, is_nullish as _is_nullish
def _gi_bool(df_general_info, key, default=False):
"""Read a flag from df_general_info, coercing legacy CSV strings and
typed YAML values to the same Python ``bool``.
Old CSV cells were strings (``"0"``, ``"False"``, ``"true"``); the YAML
migration in commit d5c01f8 made them typed (``int``/``bool``). Plain
``bool(get_csv_value(...))`` produced inconsistent results across the
two formats — e.g. ``bool("0")`` is ``True`` but ``bool(0)`` is ``False``.
This helper accepts both encodings and returns the semantically correct
boolean. Missing values fall back to ``default``.
"""
raw = get_csv_value(df_general_info, key)
if raw is None or _is_nullish(raw):
return default
try:
return _coerce_bool(raw)
except ValueError:
return default
def _gi_float(df_general_info, key, default=None):
"""Read a numeric value from df_general_info as ``float``, returning
``default`` when the cell is missing/null. Avoids ``float(None)`` crashes
on YAML ``null`` entries."""
raw = get_csv_value(df_general_info, key)
if raw is None or _is_nullish(raw):
return default
return float(raw)
def _gi_int(df_general_info, key, default=None):
"""Read a numeric value from df_general_info as ``int``, returning
``default`` when the cell is missing/null."""
raw = get_csv_value(df_general_info, key)
if raw is None or _is_nullish(raw):
return default
return int(raw)
def _gi_str(df_general_info, key, default=None):
"""Read a string value from df_general_info, returning ``default`` when
the cell is missing/null. Strips surrounding whitespace and lower-cases
so users typing ``"Isotropic"`` or ``" diagonal "`` round-trip cleanly."""
raw = get_csv_value(df_general_info, key)
if raw is None or _is_nullish(raw):
return default
return str(raw).strip().lower()
[docs]
class ModelSetup:
"""
Handles initialization and setup of core state objects for atmospheric retrieval.
This class extracts the setup methods from ModelData, responsible for reading
configuration files, parsing instrument and target parameters, and initializing
the Atmosphere, Retrieval, Bestpars, and Random objects that drive the retrieval.
Attributes:
param_handler: ParameterHandler instance for parameter management
lbl_sampling_hr: Line-by-line sampling override for high resolution
lbl_sampling_lr: Line-by-line sampling override for low resolution
path_default: Default path for file operations
model_type: String identifying the run mode (e.g. 'Retrieval', 'Model')
initial_param_array: Initial parameter values array
start_molec: Starting molecule configuration
atmosphere: Atmosphere object (populated by setup methods)
retrieval_data: Retrieval object (populated by setup methods)
bestpars_data: Bestpars object (populated by setup methods)
random_obj: Random object (populated by setup methods)
wlen_list_overplot: Wavelength lists for overplotting (populated by night_data_extraction)
table_output_file: Path to output log file
"""
[docs]
def __init__(
self,
param_handler,
lbl_sampling_hr,
lbl_sampling_lr,
path_default,
model_type,
initial_param_array,
start_molec,
):
"""
Initialize ModelSetup with configuration parameters.
Args:
param_handler: ParameterHandler instance for parameter management
lbl_sampling_hr: Line-by-line sampling override for high resolution
lbl_sampling_lr: Line-by-line sampling override for low resolution
path_default: Default path for file operations
model_type: String identifying the run mode (e.g. 'Retrieval', 'Model')
initial_param_array: Initial parameter values array
start_molec: Starting molecule configuration
"""
self.param_handler = param_handler
self.lbl_sampling_hr = lbl_sampling_hr
self.lbl_sampling_lr = lbl_sampling_lr
self.path_default = path_default
self.model_type = model_type
self.initial_param_array = initial_param_array
self.start_molec = start_molec
# These will be populated by the setup methods
self.atmosphere = None
self.retrieval_data = None
self.bestpars_data = None
self.random_obj = None
self.wlen_list_overplot = None
self.table_output_file = None
[docs]
def populate_from_manual_model(
self,
manual_model_obj,
load_new_opacities,
minwlen_lr,
maxwlen_lr,
minwlen_hr,
maxwlen_hr,
range_min_press=None,
range_max_press=None,
nlayers=None,
):
"""
Populate model data from a manual model object.
This method initializes the model using parameters from a manual model object
instead of reading from parameter files. It processes both parameter arrays
and molecule lists to set up the atmospheric model.
Args:
manual_model_obj: Manual model object containing parameters and molecules
load_new_opacities: Whether to load new opacity data
minwlen_lr: Minimum wavelength for low resolution
maxwlen_lr: Maximum wavelength for low resolution
minwlen_hr: Minimum wavelength for high resolution
maxwlen_hr: Maximum wavelength for high resolution
range_min_press: Minimum pressure range
range_max_press: Maximum pressure range
nlayers: Number of atmospheric layers
"""
# Initialize parameter lists
scale = []
rayleigh_species = []
line_species = []
line_species_isotope = []
line_species_complete_name_hr = []
line_species_complete_name_lr = []
list_bestpars = []
list_bestpars_initial_value = []
list_fixed = []
list_condensed_molecules = []
isotope = None
# Background gases (H2, He, H-, H, e-) listed in
# ``manual_model_composition`` are required by petitRADTRANS as CIA
# species (H2--H2, H2--He) regardless of whether the user explicitly
# checked their presence boxes. Treat their inclusion in
# ``manual_model_composition`` as the authoritative signal.
composition = getattr(manual_model_obj, "manual_model_composition", None) or {}
forced_inclusion = set(composition.keys()) if composition else set()
# Process parameters from manual model object
for elem in manual_model_obj.param_array:
if elem is None:
continue
considered = elem.is_considered()
elem_name_check = getattr(elem, "name_for_retrieval", None)
if not considered and elem_name_check not in forced_inclusion:
continue
# Extract parameter information
name, name_for_list_molec, molec_formula = elem.get_name()
mass_elem = elem.mass_molec
value_during_retrieval = elem.get_value()
constant_VMR = elem.VMR
rayleigh = elem.get_rayleigh_value()
# Add parameter to model
parameters = self.param_handler.add_param_manual_model(
name, value_during_retrieval, constant_VMR, mass_elem,
molec_formula, name_for_list_molec, scale, rayleigh,
rayleigh_species, line_species, line_species_isotope,
line_species_complete_name_hr, line_species_complete_name_lr, list_bestpars,
list_bestpars_initial_value, isotope, list_fixed,
list_condensed_molecules
)
# Update parameter lists with returned values
(
scale, rayleigh_species, line_species, line_species_isotope,
line_species_complete_name_hr, line_species_complete_name_lr, list_bestpars,
list_bestpars_initial_value, list_fixed,
list_condensed_molecules
) = parameters
# Process molecules from manual model object
molec_list = manual_model_obj.molecs
for molec in molec_list.split("\n"):
string_comp = molec.split(",")
# Parse molecule information
name = molec_formula = string_comp[0]
isotope = string_comp[1]
name_for_list_molec = string_comp[2]
mass_elem = float(string_comp[4])
value_during_retrieval = float(string_comp[5])
constant_VMR = int(string_comp[6])
rayleigh = int(string_comp[8])
# Add molecule parameter to model
parameters = self.param_handler.add_param_manual_model(
name, value_during_retrieval, constant_VMR, mass_elem,
molec_formula, name_for_list_molec, scale, rayleigh,
rayleigh_species, line_species, line_species_isotope,
line_species_complete_name_hr, line_species_complete_name_lr, list_bestpars,
list_bestpars_initial_value, isotope, list_fixed,
list_condensed_molecules
)
# Update parameter lists with returned values
(
scale, rayleigh_species, line_species, line_species_isotope,
line_species_complete_name_hr, line_species_complete_name_lr, list_bestpars,
list_bestpars_initial_value, list_fixed,
list_condensed_molecules
) = parameters
# Process molecules from manual model object
condensed_list = manual_model_obj.condensed
if condensed_list != "":
for condensed in condensed_list.split("\n"):
string_comp = condensed.split(",")
# Parse molecule information
name = name_for_list_molec = string_comp[0]
molec_formula = string_comp[-1]
isotope = None
mass_elem = float(string_comp[2])
value_during_retrieval = float(string_comp[3])
constant_VMR = True
rayleigh = False
# Add molecule parameter to model
parameters = self.param_handler.add_param_manual_model(
name, value_during_retrieval, constant_VMR, mass_elem,
molec_formula, name_for_list_molec, scale, rayleigh,
rayleigh_species, line_species, line_species_isotope,
line_species_complete_name_hr, line_species_complete_name_lr, list_bestpars,
list_bestpars_initial_value, isotope, list_fixed,
list_condensed_molecules
)
# Update parameter lists with returned values
(
scale, rayleigh_species, line_species, line_species_isotope,
line_species_complete_name_hr, line_species_complete_name_lr, list_bestpars,
list_bestpars_initial_value, list_fixed,
list_condensed_molecules
) = parameters
# Process molecules from manual model object
hybrid_list = manual_model_obj.hybrid
if hybrid_list != "":
for hybrid in hybrid_list.split("\n"):
string_hybrid = hybrid.split(",")
value_during_retrieval = float(string_hybrid[2])
name = name_for_list_molec = string_hybrid[0]
molec_formula = None
isotope = None
mass_elem = None
constant_VMR = True
rayleigh = False
# Add molecule parameter to model
parameters = self.param_handler.add_param_manual_model(
name, value_during_retrieval, constant_VMR, mass_elem,
molec_formula, name_for_list_molec, scale, rayleigh,
rayleigh_species, line_species, line_species_isotope,
line_species_complete_name_hr, line_species_complete_name_lr, list_bestpars,
list_bestpars_initial_value, isotope, list_fixed,
list_condensed_molecules
)
# Update parameter lists with returned values
(
scale, rayleigh_species, line_species, line_species_isotope,
line_species_complete_name_hr, line_species_complete_name_lr, list_bestpars,
list_bestpars_initial_value, list_fixed,
list_condensed_molecules
) = parameters
# Initialize core objects
self.random_obj = Random(1990) # not needed in forward model
self.bestpars_data = Bestpars(list_bestpars, list_bestpars_initial_value)
# Initialize atmosphere if needed
if load_new_opacities:
self.atmosphere = Atmosphere()
# Update atmosphere parameters
self.atmosphere.update_params(
line_species,
line_species_isotope,
line_species_complete_name_hr,
line_species_complete_name_lr,
list_condensed_molecules,
rayleigh_species,
range_min_press,
range_max_press,
nlayers,
lbl_high_res=self.lbl_sampling_hr,
lbl_low_res=1,
manual_model_obj=manual_model_obj,
chemistry=manual_model_obj.chemistry,
start_molecs=self.start_molec
)
# Initialize retrieval data
self.retrieval_data = Retrieval(scale, manual_model_obj=manual_model_obj)
# Set up opacities and wavelength ranges if needed
if load_new_opacities:
self.atmosphere.set_wl_range(
minwlen_lr=minwlen_lr, maxwlen_lr=maxwlen_lr,
minwlen_hr=minwlen_hr, maxwlen_hr=maxwlen_hr,
)
self.atmosphere.read_opacities(
table_output_file=manual_model_obj.table_output_file,
path_target=manual_model_obj.path_targets,
stellar_spectrum_type_hr=manual_model_obj.path_stellar_spectrum
)
[docs]
def read_df_information(
self,
path_df_information,
df_parameters,
id_process,
table_output_file,
minwlen_lr,
maxwlen_lr,
minwlen_hr,
maxwlen_hr,
range_min, range_max, nlayers
):
"""
Read and process general information from configuration CSV file.
This method reads the main configuration file containing instrument setup,
atmospheric parameters, target properties, and retrieval settings. It
initializes all core objects needed for the atmospheric retrieval including
instruments, atmosphere, target, and retrieval configuration.
Args:
path_df_information: Path to the general information CSV file
df_parameters: Processed parameter data from read_df_parameters
id_process: Process ID for parallel operations
table_output_file: Output file for logging
minwlen_lr: Minimum wavelength for low resolution
maxwlen_lr: Maximum wavelength for low resolution
minwlen_hr: Minimum wavelength for high resolution
maxwlen_hr: Maximum wavelength for high resolution
range_min: Minimum pressure range override (overrides CSV value if provided)
range_max: Maximum pressure range override (overrides CSV value if provided)
nlayers: Number of atmospheric layers for retrieval
"""
# Extract processed parameter data via NAMED attribute access
# (df_parameters is a ReadDfParametersResult NamedTuple). Avoids the
# silent-reorder failure mode of the previous positional indexing.
mass_vector = df_parameters.mass_vector
scale = df_parameters.scale
rayleigh_species = df_parameters.rayleigh_species
line_species = df_parameters.line_species
line_species_isotope = df_parameters.line_species_isotope
line_species_complete_name_hr = df_parameters.line_species_complete_name_hr
line_species_complete_name_lr = df_parameters.line_species_complete_name_lr
list_bestpars = df_parameters.list_bestpars
list_bestpars_initial_value = df_parameters.list_bestpars_initial_value
list_fixed = df_parameters.list_fixed
list_condensed_molecules = df_parameters.list_condensed_molecules
# Per-parameter prior bounds for ``in_bestpars=1`` rows, aligned 1:1
# with ``scale``. Forwarded to Retrieval and consumed by the MAP
# Hessian extractor (compute_hessian_at_map) for the σ_p cap.
range_min_vector = df_parameters.range_min_vector
range_max_vector = df_parameters.range_max_vector
# Read general configuration information.
# The folder containing the file is what the YAML/CSV-fallback reader
# needs; the filename component of ``path_df_information`` is ignored.
df_general_info = read_general_info_df(Path(path_df_information).parent)
# Set up petitRADTRANS data path environment
input_data_folder = get_csv_value(df_general_info, "path_petitRadTrans")
if input_data_folder[-1] != os.sep:
input_data_folder += os.sep
os.environ["pRT_input_data_path"] = input_data_folder
# Initialize high-resolution instruments
instrument_HR = []
list_instruments_hr = get_csv_value(df_general_info, "instruments_hr_list") or []
total_nights = 0
total_hr_instruments = 0
min_hwhm_hr = 999
#
use_hr_linelist_for_lr = _gi_bool(df_general_info, "use_hr_linelist_for_lr")
if list_instruments_hr and list_instruments_hr != ["Nothing"]:
# Parse instrument string format: name&path&pixel_dv&min_wl&max_wl&nights
list_instruments = list_instruments_hr
for index_inst, curr_instrument in enumerate(list_instruments):
curr_instrument = curr_instrument.split("&")
pixel_dv = float(curr_instrument[2])
min_hwhm_hr = min(min_hwhm_hr, pixel_dv)
instrument_HR.append(
Instrument(
curr_instrument[0], curr_instrument[1], "HR", pixel_dv,
float(curr_instrument[3]), float(curr_instrument[4]),
nights=curr_instrument[5]
)
)
total_hr_instruments += 1
total_nights += instrument_HR[index_inst].ndate
# Initialize low-resolution instruments
instrument_LR = []
list_instruments_lr = get_csv_value(df_general_info, "instruments_lr_list") or []
if list_instruments_lr and list_instruments_lr != ["Nothing"]:
list_instruments = list_instruments_lr
for curr_instrument in list_instruments:
curr_instrument = curr_instrument.split("&")
# Handle optional pixel_dv parameter
if curr_instrument[2] != "":
try:
pixel_dv = float(curr_instrument[2])
except Exception as e:
pixel_dv = None
# min_hwhm_lr = min(min_hwhm_lr, pixel_dv)
else:
pixel_dv = None
# Handle optional wavelength ranges
min_wl = float(curr_instrument[3]) if curr_instrument[3] != "" else None
max_wl = float(curr_instrument[4]) if curr_instrument[4] != "" else None
instrument_LR.append(
Instrument(
curr_instrument[0], curr_instrument[1], "LR",
pixel_dv, min_wl, max_wl
)
)
#
chemistry = get_csv_value(df_general_info, "chemistry")
# chemcat_bool = get_csv_value(df_general_info, "chemcat_bool")
# chemcat_bool = chemcat_bool == "True"
#
format_temperature = ConstantVariables.normalize_pt_profile(
get_csv_value(df_general_info, "t_p_profile")
)
#
continum_opacity_arr = get_csv_value(
df_general_info, "continum_opacity"
)
continum_opacity = continum_opacity_arr.split(",")
#
maxsteps = _gi_int(df_general_info, "maxsteps")
#
seme = _gi_int(df_general_info, "seed")
#
multiplier_chains = get_csv_value(df_general_info, "multiplier_chains")
multiplier_chains = 2 if multiplier_chains is None else int(multiplier_chains)
#
multiplier_cores = get_csv_value(df_general_info, "multiplier_cores")
multiplier_cores = 1 if multiplier_cores is None else int(multiplier_cores)
#
use_pool = _gi_bool(df_general_info, "use_pool", default=False)
#
use_parallel_init = _gi_bool(df_general_info, "use_parallel_init", default=False)
#
order_sel = get_csv_value(df_general_info, "order_sel")
#
# Transmission or emission
rad_mode = get_csv_value(df_general_info, "rad_mode")
#
target = Target(path_df_information, rad_mode=rad_mode)
path_results = get_csv_value(df_general_info, "path_results")
path_target = Path(path_results).parent.parent
#
ecc_opi = _gi_bool(df_general_info, "ecc_opi")
#
tell_rm_method = get_csv_value(df_general_info, "tell_rm_method")
#
resolution = get_csv_value(df_general_info, "resolution")
#
stellar_spectrum_type_hr = get_csv_value(df_general_info, "stellar_spectrum_type_hr")
stellar_spectrum_type_lr = get_csv_value(df_general_info, "stellar_spectrum_type_lr")
stellar_spectrum_type_lr_mode = get_csv_value(df_general_info, "stellar_spectrum_type_lr_mode")
if stellar_spectrum_type_hr is None or stellar_spectrum_type_hr == "None":
stellar_spectrum_type_lr, = stellar_spectrum_type_hr = get_csv_value(df_general_info, "stellar_spectrum")
print("stellar_spectrum_type_lr:", stellar_spectrum_type_lr, " Mode:", stellar_spectrum_type_lr_mode)
print("stellar_spectrum_type_hr:", stellar_spectrum_type_hr)
if rad_mode == "Transmission":
stellar_spectrum_type_hr = None
stellar_spectrum_type_lr = None
#
range_min_pressures = (
float(get_csv_value(df_general_info, "range_min_pressures")) if range_min is None
else range_min
)
#
range_max_pressures = (
float(get_csv_value(df_general_info, "range_max_pressures")) if range_max is None
else range_max
)
#
layers = (
int(get_csv_value(df_general_info, "layers")) if nlayers is None
else nlayers
)
#
lbl_hr = get_csv_value(df_general_info, "lbl_hr")
if lbl_hr is None:
lbl_hr = int(get_csv_value(df_general_info, "lbl")) # RETROCOMPATIBILITY
else:
lbl_hr = int(lbl_hr)
lbl_hr = (
lbl_hr if self.lbl_sampling_hr is None else self.lbl_sampling_hr
)
#
lbl_lr = get_csv_value(df_general_info, "lbl_lr")
lbl_lr = (
lbl_lr if self.lbl_sampling_lr is None else self.lbl_sampling_lr
)
if lbl_lr is None:
lbl_lr = 4
else:
lbl_lr = int(lbl_lr)
#
model_reprocessing = get_csv_value(df_general_info, "model_reprocessing")
#
standardize_PCA = get_csv_value(df_general_info, "standardize_PCA")
if standardize_PCA is None:
standardize_PCA = "False"
#
rv_sampling = get_csv_value(df_general_info, "rv_sampling")
rv_sampling = float(rv_sampling) if rv_sampling is not None else 0.1
#
self.random_obj = Random(seme)
#
retrieval_model = get_csv_value(df_general_info, "retrieval_model")
#
# Mask phase settings. ``mask_phase_enabled`` is read with the bool
# coercer so that a YAML ``0`` and a legacy CSV ``"0"`` both disable
# masking (previously the YAML int branch entered the else and
# crashed on ``float(None)`` when ``min_phase: null``).
mask_phase_enabled = _gi_bool(df_general_info, "mask_phase_enabled")
if not mask_phase_enabled:
min_phase = None
max_phase = None
mask_night_index = None
mask_inside_range = None
else:
mask_night = get_csv_value(df_general_info, "mask_night")
mask_night_index = (
[int(value) for value in str(mask_night).split(",")]
if mask_night is not None else None
)
min_phase = _gi_float(df_general_info, "min_phase")
max_phase = _gi_float(df_general_info, "max_phase")
mask_inside_range = get_csv_value(df_general_info, "which_mask")
if mask_inside_range is None or mask_inside_range == "Inside range":
mask_inside_range = True
else:
mask_inside_range = False
transit_limit_mode = get_csv_value(df_general_info, "transit_limit_mode")
if transit_limit_mode is None:
transit_limit_mode = "T14"
#
nthin = get_csv_value(df_general_info, "nthin")
nthin = 1 if nthin is None else int(nthin)
#
mode_jump_threshold = get_csv_value(df_general_info, "mode_jump_threshold")
mode_jump_threshold = 0.05 if mode_jump_threshold is None else float(mode_jump_threshold)
#
init_scatter_factor = get_csv_value(df_general_info, "init_scatter_factor")
init_scatter_factor = 2.0 if init_scatter_factor is None else float(init_scatter_factor)
#
epsilon_scale_divisor = get_csv_value(df_general_info, "epsilon_scale_divisor")
epsilon_scale_divisor = 100.0 if epsilon_scale_divisor is None else float(epsilon_scale_divisor)
if epsilon_scale_divisor <= 0:
raise ValueError(
f"epsilon_scale_divisor must be > 0 (got {epsilon_scale_divisor}). "
f"This parameter divides the scale in the uniform epsilon perturbation "
f"(U-0.5)*2*scale/epsilon_scale_divisor; zero or negative values are undefined."
)
#
use_map_optimizer = _gi_bool(df_general_info, "use_map_optimizer", default=False)
# MAP optimizer knobs (defaults mirror find_map_start signature).
# No separate seed: run_map_warm_start derives one from the MCMC
# SeedSequence via spawn() so a single user-facing seed controls both.
map_maxiter_global = _gi_int(df_general_info, "map_maxiter_global", default=40)
map_popsize = _gi_int(df_general_info, "map_popsize", default=10)
map_tol_global = _gi_float(df_general_info, "map_tol_global", default=1.0e-2)
map_maxiter_local = _gi_int(df_general_info, "map_maxiter_local", default=150)
# Chain-init scatter strategy after MAP. Default ``"isotropic"`` is
# byte-identical to legacy behaviour; ``"diagonal"`` and
# ``"correlated"`` opt into Hessian-based scatter (see map_optimizer).
# Validated inside Retrieval.__init__ — invalid values raise loudly.
init_mode = _gi_str(df_general_info, "init_mode", default="isotropic")
# MCMC sampler choice. Default ``"DE-MC"`` keeps the legacy byte-identical
# path; ``"Snooker"`` selects DE-MCzs. Validated inside Retrieval.__init__.
mcmc_sampler = _gi_str(df_general_info, "mcmc_sampler", default="DE-MC")
# Gelman-Rubin R-hat threshold for declaring MCMC convergence (canonical
# default 1.01). Validated inside Retrieval.__init__ (must be > 1.0).
gelman_rubin_threshold = _gi_float(
df_general_info, "gelman_rubin_threshold", default=1.01
)
# Number of consecutive G-R checks required to declare convergence and stop
# (canonical default 6). Validated inside Retrieval.__init__ (must be >= 1).
maxnpass = _gi_int(df_general_info, "maxnpass", default=6)
# Chain debug JSONL stream toggle. Default True for backward compatibility:
# legacy runs (CSV without the field) keep producing the debug log so the
# GUI "Analyse Chains" tool works unchanged. Set to False from the GUI
# checkbox to skip the per-step writes on long production runs.
save_chain_debug = _gi_bool(df_general_info, "save_chain_debug", default=True)
#
# fmt: off
self.bestpars_data = Bestpars(
list_bestpars, list_bestpars_initial_value, multiplier_chains, multiplier_cores, use_pool, use_parallel_init,
mode_jump_threshold=mode_jump_threshold,
)
#
additional_opac = False
for name in ["k_opac", "k_cond", "lambda0_micron", "omega_scale_micron", "xi"]:
if name in list_bestpars or name in list_fixed:
additional_opac = True
#
self.atmosphere = Atmosphere()
self.atmosphere.update_params(
line_species, line_species_isotope, line_species_complete_name_hr, line_species_complete_name_lr, list_condensed_molecules,
rayleigh_species, range_min_pressures, range_max_pressures,
layers, lbl_high_res=lbl_hr, lbl_low_res=lbl_lr,
resolution=resolution, instruments_HR=instrument_HR, instruments_LR=instrument_LR,
continum_opacity=continum_opacity,
mass_vector=mass_vector,
target=target,
initial_params=self.initial_param_array,
path_default=self.path_default,
rad_mode=rad_mode,
retrieval_model=retrieval_model, chemistry=chemistry, rv_sampling=rv_sampling, start_molecs=self.start_molec,
additional_opac=additional_opac
)
self.retrieval_data = Retrieval(
scale,
ecc_opi=ecc_opi, tell_rm_method=tell_rm_method,
model_reprocessing=model_reprocessing,
standardize_PCA=standardize_PCA,
format_temperature=format_temperature,
table_output_file=table_output_file,
maxsteps=maxsteps, order_sel=order_sel, id_process=id_process, path_results=path_results,
retrieval_model=retrieval_model, total_nights=total_nights, total_hr_instruments=total_hr_instruments,
mask_phase_enabled=mask_phase_enabled, mask_night_index=mask_night_index, min_phase=min_phase, max_phase=max_phase,
mask_inside_range=mask_inside_range, transit_limit_mode=transit_limit_mode,
nthin=nthin,
init_scatter_factor=init_scatter_factor,
epsilon_scale_divisor=epsilon_scale_divisor,
use_map_optimizer=use_map_optimizer,
map_maxiter_global=map_maxiter_global,
map_popsize=map_popsize,
map_tol_global=map_tol_global,
map_maxiter_local=map_maxiter_local,
init_mode=init_mode,
mcmc_sampler=mcmc_sampler,
range_min_vector=range_min_vector,
range_max_vector=range_max_vector,
gelman_rubin_threshold=gelman_rubin_threshold,
maxnpass=maxnpass,
save_chain_debug=save_chain_debug,
)
self.night_data_extraction()
self.atmosphere.set_wl_range(
minwlen_lr=minwlen_lr, maxwlen_lr=maxwlen_lr,
minwlen_hr=minwlen_hr, maxwlen_hr=maxwlen_hr,
)
if "Retrieval" in self.model_type:
with open(table_output_file, "a") as f:
f.write(f"Initializing {retrieval_model} module\n")
if retrieval_model == 'petitRADTRANS':
self.atmosphere.read_opacities(
table_output_file=table_output_file , path_target=path_target, stellar_spectrum_type_lr=stellar_spectrum_type_lr,
stellar_spectrum_type_hr=stellar_spectrum_type_hr, min_hwhm_hr=min_hwhm_hr,
instrument_LR=instrument_LR, stellar_spectrum_type_lr_mode=stellar_spectrum_type_lr_mode,
use_hr_linelists_for_lr=use_hr_linelist_for_lr
)
elif retrieval_model == 'PyratBay':
output_folder = str(Path(table_output_file).parent)
self.atmosphere.init_pyratbay(output_folder, self.initial_param_array)
[docs]
def night_data_extraction(self):
"""
Extract and process night observation data for each instrument.
This method handles the extraction and organization of night observation data,
including path setup, file archiving (for retrieval mode), and creation of
Night objects for each observation date. The method handles different
processing modes (Model, Manual, Retrieval) with appropriate file handling.
"""
tar_folder = None
# Set up temporary folder for model plotting
if "Model" in self.model_type:
username = os.environ.get("USER")
tar_folder = str(
Path(self.path_default, "GUIBRUSHR", "Files", "Temp_Files_GUI",
username, "tar_unpack")
)
# Initialize wavelength list for overplotting
self.wlen_list_overplot = []
# Exit if no high-resolution instruments available
if self.atmosphere.resolution_obj.instruments_HR is None:
return
# Process each high-resolution instrument
for instrument in self.atmosphere.resolution_obj.instruments_HR:
for i in range(instrument.ndate):
# Handle Model and Manual modes (use extracted tar files)
if "Model" in self.model_type or "Manual" in self.model_type:
path_night_reference = str(
Path(tar_folder, instrument.name)
)
path_order_tell = str(
Path(path_night_reference, instrument.dates[i])
)
if "Model" in self.model_type:
# Clean and extract tar archive for model plotting
os.system("rm -rf " + path_order_tell)
os.system("mkdir -p " + path_order_tell)
os.system(
"tar -xzf "
+ self.retrieval_data.path_results
+ instrument.dates[i]
+ "_" + instrument.name
+ ".tar.gz --transform='s/.*\///' -C "
+ path_order_tell
)
path_order_test = path_order_tell
path_good_order = str(
Path(path_order_test, "good_order_" + self.retrieval_data.order_sel + ".pkl")
)
else:
# Handle Retrieval mode (create and archive files)
# Set up paths for order selection and processing
path_night_reference = instrument.path_instrument
path_order = str(Path(
instrument.path_instrument, instrument.dates[i], "orders_selection", self.retrieval_data.order_sel
))
# Define paths for different processing stages
path_good_order = str(
Path(path_order, "good_order_" + self.retrieval_data.order_sel + ".pkl")
)
path_order_tell = str(
Path(path_order, self.retrieval_data.tell_rm_method)
)
path_order_test = str(Path(path_order_tell, "test"))
# Define paths for observation files
path_aligned = str(
Path(instrument.path_instrument, instrument.dates[i], "aligned.fits")
)
path_wlen = str(
Path(instrument.path_instrument, instrument.dates[i], "wlen_refined.fits")
)
path_fase = str(
Path(instrument.path_instrument, instrument.dates[i], "phase.pkl")
)
if not os.path.exists(path_fase):
path_fase = str(
Path(instrument.path_instrument, instrument.dates[i], "fase.pkl")
)
# Create archive with all necessary files
final_folder = str(
Path(self.retrieval_data.path_results,
instrument.dates[i] + "_" + instrument.name + ".tar.gz")
)
# Archive telluric removal results and observation files
os.system(
f"tar -zcvf - {path_order_tell} {path_good_order} "
f"{path_aligned} {path_wlen} {path_fase} > {final_folder}"
)
# Copy resolution string file
os.system(
"cp " + str(Path(path_order_tell, "string_res.txt"))
+ " " + self.retrieval_data.path_results
)
# Create Night object for this observation
target_limit_phase = (
self.atmosphere.target.limit_phase_t14 if self.retrieval_data.transit_limit_mode == "T14"
else self.atmosphere.target.limit_phase_t23
)
night = Night(
target_limit_phase,
path_night_reference,
instrument.dates[i],
self.atmosphere.rad_mode,
path_good_order,
path_order_tell,
path_order_test,
table_output_file=self.table_output_file, CC=False, simulation=False, synthetic=False,
retrieval_standardize=self.retrieval_data.standardize_PCA
)
# Store wavelength list and add night to instrument
self.wlen_list_overplot.append(night.wlen_list_overplot)
instrument.append_night_obj(night)