Automating Solar Model Calibration with PySEP
This Python script is dedicated to automating the calibration of solar models using data from various sources, including atmospheric data, opacity files, and composition data. The script utilizes the PySEP library to calibrate stellar models based on input parameters and output the calibrated models into a specified directory. It includes functionality for handling command-line arguments, reading input files, adjusting model parameters, and performing calibration with error handling.
from pysep.dsep import calibrate from pysep.dsep import stellarModel as sm from pysep.dm import filetypes as ft 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 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!")