Solar Calibration Script using PySEP

This Python script performs solar calibration for stellar models using the PySEP package. It reads composition, atmospheric, and opacity model files from input CSV files, generates pre-main sequence (preMS) models, and carries out solar calibration. The calibrated models are saved to a specified output directory. The script includes functionalities for parsing atmospheric model files, low and high temperature opacity files, and composition files. Error handling is implemented for situations where calibration fails for certain models.

from pysep.dsep import calibrate
from pysep.dsep import stellarModel as sm
from pysep.dm import filetypes as ft
from pysep.api.starting import get_basic_stellar_model
 
from pysep.io.nml.control.defaults import solar
from pysep.io.nml.physics.defaults import phys1
 
from pysep.opac.tops import open_and_parse
 
from pysep.newpoly.popFileToPolyModel import generate_prems_models
 
import os
 
import pandas as pd
 
import pickle as pkl
 
OPLIB = 
 
if __name__ == "__main__":
    import argparse
    parser = argparse.ArgumentParser(description="Run Solar Calibration for models")
    parser.add_argument("inputs", help="path to csv with model file paths listed", type=str)
    parser.add_argument("outputDir", help="path to save calibrated results too", type=str)
    parser.add_argument("-t", "--tau", help="atmospheric fitting point", type=str, default="010")
 
    args = parser.parse_args()
 
    if os.path.exists(args.outputDir):
        raise RuntimeError("Error! Will not run calibration if outputDir exists")
    os.mkdir(args.outputDir)
 
    modelsToBuild = pd.read_csv(os.path.join(args.inputs, "composition.dat"), comment='#')
 
    for idx, modelFiles in modelsToBuild.iterrows():
        atmDirPath = os.path.join(args.inputs, modelFiles['atm'])
        files = os.listdir(atmDirPath)
        atmFilePath = next(filter(lambda x: f"tau{args.tau}.sbc.nc" in x, files))
        atmFile = ft.bcphx95_file(atmFilePath)
 
        opacHFile = ft.opacf_file(os.path.join(args.inputs, modelFiles['opacH']))
 
 
        lowTempOpacFileNames = os.listdir(os.path.join(args.inputs, modelFiles['opacL']))
        lowTempOpacPaths = [os.path.join(args.inputs, modelFiles['opacL'],  x) for x in lowTempOpacFileNames]
 
        parsedAbunFile = open_and_parse(os.path.join(args.inputs, modelFiles['comp']))
 
        controlNML = solar.copy()
        physicsNML = phys1.copy()
 
        controlNML.update_low_temp_opacity_files(lowTempOpacPaths)
 
        premsModelPath = generate_prems_models(parsedAbunFile, [1.0], "premsTMP/", compName=modelFiles['name'], DDAGE=1000, fnameT="temp_m%e%m")
        premsFile = ft.premsf_file(os.path.join('premsTMP/', premsModelPath[1.0]))
 
        modelOutputPath = os.path.join(args.outputDir, modelFiles['name'])
        os.mkdir(modelOutputPath)
 
        model = sm(modelOutputPath, controlNML, physicsNML, opacHFile, premsFile)
        model.update_composition(X=parsedAbunFile['AbundanceRatio']['X'], Z=parsedAbunFile['AbundanceRatio']['Z'])
 
 
 
        try:
            calibResult = calibrate.LRZX_calib_solar(model, 0.7, 0.02, 1.9)
 
            with open(os.path.join(args.outputDir, f"{modelFiles['name']}_calib.pkl"), 'wb') as f:
                pkl.dump(calibResult,f)
        except:
            print(f"Error Calibrating {idx}, moving on!")