Skip to content

Instantly share code, notes, and snippets.

@astro313
Created July 25, 2019 21:23
Show Gist options
  • Save astro313/72fe18920aec8c64d2e4ce555c741aa8 to your computer and use it in GitHub Desktop.
Save astro313/72fe18920aec8c64d2e4ce555c741aa8 to your computer and use it in GitHub Desktop.
check nHden from cloudy output hden and saved grid param file
gridfile = '/mnt/home/daisyleung/Downloads/SIGAME_dev/sigame/temp/z6_data_files/cloud_models/GMC/grids/GMCgrid_abun_z6_em.models'
import glob
import pandas as pd
import numpy as np
gridsDF = pd.read_pickle(gridfile)
nHmw = gridsDF['nHmw']
nHmin = gridsDF['nHmin']
#
outPath = '/mnt/home/daisyleung/Downloads/SIGAME_dev/sigame/temp/z6_data_files/cloud_models/GMC/output/'
def line_num_for_phrase_in_file(phrase='xxx', filename='file.txt'):
totNumLines = np.loadtxt(filename).shape[0]
with open(filename, 'r') as f:
i = 0
magicnumber_list = []
while i < totNumLines:
for (i, line) in enumerate(f):
if phrase in line:
magicnumber_list.append(i - 2)
return magicnumber_list
return -1
def read_ovrfile(ovrfile):
# read in overview file
try:
ovrdata = np.loadtxt(ovrfile, skiprows=1)
except:
return ovrfile, None
f = open(ovrfile, 'r')
ovrhdr = f.readline().split('\t')
f.close()
ovrhdr[0] = ovrhdr[0][1:]
col = dict.fromkeys(ovrhdr)
for ii, jj in enumerate(ovrhdr): # clean off any \n \r\n
ovrhdr[ii] = jj.rstrip()
for i in range(ovrdata.shape[-1]):
col[ovrhdr[i]] = ovrdata[:, i]
magicnumber_list = line_num_for_phrase_in_file("####", ovrfile)
return col, magicnumber_list
import _pickle as pick
param = pick.load(open('/mnt/home/daisyleung/Downloads/SIGAME_dev/sigame/temp/z6_data_files/cloud_models/GMC/grids/GMCgrid_abun_z6.param', "rb"))
pc2cm = 3.086e+18
for iii, prefix in enumerate(glob.glob(outPath + '*.in')):
prefix = prefix[:prefix.find('.in')]
prefix = prefix.split('/')[-1]
ovrfile = outPath + prefix + '.ovr'
col, magicnumber_list = read_ovrfile(ovrfile)
if type(col) is str and magicnumber_list is None:
pass
elif len(magicnumber_list) == 0:
# haven't converged
pass
else:
Rgmc = param['Rgmcs'][iii] # pc
col=pd.DataFrame(col)
import pdb; pdb.set_trace()
col = col[magicnumber_list[-1]-10:]
print(col['hden'])
col[col['depth'] > 0.01 * Rgmc * pc2cm].reset_index(drop=True)
nH = col['hden']
print(nH.min(), nHmin[int(prefix[prefix.find('_')+1:])])
import pdb; pdb.set_trace()
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment