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!")