1+ import numpy as np
2+ import pandas as pd
3+ from astropy .io import fits
4+ from astropy .table import Table
5+ import os
6+ import hoki
7+ from hoki import load
8+ import glob
9+ import re
10+
11+ """Source: https://stackoverflow.com/questions/44369504/
12+ how-to-convert-entire-dataframe-values-to-float-in-pandas"""
13+
14+ colsToElim = [x for x in range (96 ) if x not in [hoki .dummy_dict [y ]
15+ for y in hoki .dummy_dict .keys ()]]
16+ def merger_related (f ):
17+ # Code that eventually returns whether
18+ # The given model (in form of dataframe)
19+ # has the conditions necessary
20+ # to be a merger model.
21+ return False
22+
23+ def find_mergers (f ):
24+ # returns when the merger occurs (first time)
25+ # assuming that f is a merger model.
26+ return False
27+
28+
29+ def reformatter (bevo_dir ):
30+ """I will expedite the process of writing and saving files by
31+ making the BPASS stellar input files binary fits tables This method
32+ creates a reformatted FITS version of
33+ all of the available BPASS data so that data can be loaded quickly
34+ and processed in the extractor method."""
35+
36+ # The list of possible BPASS metallicities.
37+ # HMG stands for Qasi-Chemically Homogeneous Evolution Models.
38+
39+ metallicities = ["zem5" , "zem4" , "z001" , "z002" , "z003" , "z004" , "z006" , "z008" ,
40+ "z010" , "z014" ,"z020" ,"z030" , "z040" , "zem5_hmg" ,
41+ "zem4_hmg" , "z001_hmg" , "z002_hmg" , "z003_hmg" , "z004_hmg" ]
42+ # I will use the Python HashSet since (assuming that strings
43+ # have a good HashCode, adding to the HashSet should be fast
44+
45+ setofAll = set ()
46+ for x in metallicity :
47+ a = set (glob .glob (hoki .MODELS_PATH + "/NEWSINMODS/" + x + "/*" ))
48+ b = set (glob .glob (hoki .MODELS_PATH + "/NEWBINMODS/NEWBINMODS/" + x + "/*" ))
49+ c = set (glob .glob (hoki .MODELS_PATH + "/NEWBINMODS/NEWSECMODS/" + x + "_2/*" ))
50+ setofAll = setofAll .union (a , b , c )
51+
52+ for x in setofAll :
53+
54+ f = pd .read_csv (hoki .MODELS_PATH + x , sep = "\s+" ,
55+ header = None )
56+ f .apply (pd .to_numeric , errors = 'coerce' )
57+ f .transpose ()
58+
59+ # Drop all of the columns that are not defined by the
60+ # HOKI package provided mapping of column name to
61+ # column number.
62+ f .drop (columns = [f .columns [x ] for x in colsToElim ],
63+ inplace = True )
64+ # Source: https://stackoverflow.com/questions/483666/
65+ # reverse-invert-a-dictionary-mapping
66+ # I create a mapping from column NUMBER to
67+ # column name
68+ invmap = {u : v for v , u in hoki .dummy_dict .items ()}
69+ f .rename (columns = invmap , inplace = True )
70+ f ["model_type" ] = [typ ] * len (f .index )
71+ astroTable = Table .from_pandas (f )
72+ # We will be recalculating the photometry related stuff
73+ # when we create the isochrone.
74+ # Will remove photometry related columns
75+ astroTable .remove_columns (['V-I' , 'U' ,
76+ 'B' , 'V' , 'R' , 'I' , 'J' , 'H' , 'K' ,'u' ,
77+ 'g' , 'r' , 'i' , 'z' , 'f300w' , 'f336w' ,
78+ 'f435w' , 'f450w' , 'f555w' , 'f606w' ,'f814w' ,
79+ 'U2' , 'B2' , 'V2' ,'R2' , 'I2' ,'J2' ,'H2' ,'K2' ,
80+ 'u2' ,'g2' ,'r2' ,'i2' ,'z2' ,'f300w2' ,'f336w2' ,
81+ 'f435w2' ,'f450w2' ,'f555w2' ,'f606w2' ,'f814w2' ,
82+ 'Halpha' ,'FUV' ,'NUV' ])
83+ # Merger Related indicates whether the star is
84+ # a model with a merger in it.
85+ astroTable ['Merger_Related' ] = [False ] * len (astroTable )
86+ # Creates container directory for stellar models grouped by metallicity
87+ # and whether they are single or binary related
88+ # the latter included mergers, QHE Models
89+ if not os .path .isdir (bevo_dir + "iso/" + metallicity + '/' +
90+ "sin" +
91+ "/" ):
92+ os .makedirs (bevo_dir +
93+ 'iso/' + metallicity +
94+ "/" + SinOrBin + "/" )
95+ if not os .path .isdir (bevo_dir + "iso/" + metallicity + '/' +
96+ "bin" +
97+ "/" ):
98+ os .makedirs (bevo_dir +
99+ 'iso/' + metallicity +
100+ "/" + "bin" + "/" )
101+ # Single star models that are not the hmg models go into the <Metallicity>/sin
102+ # directory.
103+ if (x [0 :10 ] == 'NEWSINMODS' and
104+ not re .matches (r"NEWSINMODS/z[0-9]+_hmg/.*" )):
105+ astroTable .write (bevo_dir + 'iso/' + met + '/' +
106+ "sin" + "/" + "FitsModels" +
107+ x [16 :] + ".fits" , format = 'fits' ,
108+ overwrite = True )
109+ print (x [16 :] + ".fits" )
110+ else :
111+ # Models that count as secondary star with compact remnants
112+ # go into the <metallicity>/bin subdirectory of iso.
113+ if (x [11 :21 ] == 'NEWSECMODS' ):
114+ astroTable .write (bevo_dir + 'iso/' + met + '/' +
115+ "bin" + "/" + "FitsModels" +
116+ x [29 :] + "sec.fits" , format = 'fits' ,
117+ overwrite = True )
118+ print (bevo_dir + 'iso/' + metallicity + '/' +
119+ SinOrBin + "/" + "FitsModels" +
120+ x [29 :] + "sec.fits" )
121+ else :
122+ # Check for merger related (models_type=0) models
123+ # and mark them.
124+ if (find_mergers (astroTable )):
125+ astroTable ['Merger_Related' ] = [True ] * len (astroTable )
126+ # Models that count as secondary star with compact remnants
127+ # go into the <metallicity>/bin subdirectory of iso.
128+ astroTable .write (bevo_dir + 'iso/' + metallicity + '/' +
129+ "bin" + "/" + "FitsModels" +
130+ x [27 :] + "bin.fits" , format = 'fits' ,
131+ overwrite = True )
132+ print (x [27 :] + ".fits" )
133+
134+ except FileNotFoundError :
135+ print ("File not Readable/Found" )
136+ notFound .append (x )
137+ #The following is commented out since I have yet to delete this version yet.
138+ """ still_to_do=set()
139+
140+
141+ if (x=="zem5" or x=="zem4" or x=="z001" or x=="z002" or x=="z003" or x=="z004"):
142+ c=set(c).union(set(glob.glob(hoki.MODELS_PATH+"/NEWSINMODS/"+x+"hmg"+"/*")))
143+ still_to_do=still_to_do.union(a, b, c)
144+ still_to_do=still_to_do.difference(included)
145+ lenModPath=len(HOKI.MODELS_PATH)
146+ for x in still_to_do:
147+ try:
148+ f = pd.read_csv(x, sep="\s+",
149+ header=None)
150+ f.apply(pd.to_numeric, errors='coerce')
151+ f.transpose()
152+ f.drop(columns=[f.columns[x] for x in colsToElim], inplace=True)
153+ included.add(hoki.MODELS_PATH+x)
154+ # Source: https://stackoverflow.com/questions/483666/
155+ # reverse-invert-a-dictionary-mapping
156+ invmap = {u: v for v, u in hoki.dummy_dict.items()}
157+ f.rename(columns=invmap, inplace=True)
158+ f["model_type"]=[typ]*len(f.index)
159+ astroTable = Table.from_pandas(f)
160+ astroTable.remove_columns(['?','modelimf','mixedimf', 'V-I', 'U',
161+ 'B', 'V', 'R', 'I', 'J', 'H', 'K','u',
162+ 'g', 'r', 'i', 'z', 'f300w', 'f336w',
163+ 'f435w', 'f450w', 'f555w', 'f606w','f814w',
164+ 'U2', 'B2', 'V2','R2', 'I2','J2','H2','K2',
165+ 'u2','g2','r2','i2','z2','f300w2','f336w2',
166+ 'f435w2','f450w2','f555w2','f606w2','f814w2',
167+ 'Halpha','FUV','NUV'])
168+ if not os.path.isdir(bevo_dir + "iso/" + metallicity + '/' +
169+ SinOrBin + "/"):
170+ os.makedirs(bevo_dir + 'iso/' + metallicity +
171+ "/" + SinOrBin + "/")
172+ if (x[lenModPath+0:lenModPath+10]=='NEWSINMODS'):
173+ nameMatcher=re.match(re".*/NEWSINMODS/[a-z0-9]+hmg/.*")
174+ if (bool(nameMatcher)):
175+ astroTable["model_type"]=[4]*len(astroTable)
176+ else:
177+ astroTable["model_type"]=[-12]*len(astroTable)
178+ astroTable.write(bevo_dir + 'iso/' + met + '/' +
179+ SinOrBin + "/" + "FitsModels" +
180+ x[lenModPath+16:] + ".fits", format='fits',
181+ overwrite=True)
182+ print(x[lenModPath+16:]+".fits")
183+ else:
184+ print(x[lenModPath+11:lenModPath+21])
185+
186+ if (x[lenModPath+11:lenModPath+21]=='NEWSECMODS'):
187+ astroTable.write(bevo_dir + 'iso/' + met + '/' +
188+ SinOrBin + "/" + "FitsModels" +
189+ x[lenModPath+29:] + "sec.fits", format='fits',
190+ overwrite=True)
191+ print(bevo_dir + 'iso/' + metallicity + '/' +
192+ SinOrBin + "/" + "FitsModels" +
193+ x[lenModPath+29:] +"sec.fits")
194+ # Added designation for secondary models that can be
195+ # either regular or single star
196+ astroTable["model_type"]=[-15]*len(astroTable)
197+
198+ else:
199+ astroTable.write(bevo_dir + 'iso/' + metallicity + '/' + SinOrBin +
200+ "/" + "FitsModels" +
201+ x[lenModPath+27:] +"bin.fits", format='fits',
202+ overwrite=True)
203+ # Added designation for secondary models that can be
204+ # either regular or single star
205+ astroTable["model_type"]=[-10]*len(f.index)
206+ print(x[lenModPath+27:] + ".fits")
207+ """
208+
209+ except FileNotFoundError :
210+ print ("File not Readable/Found" )
211+ notFound .append (x )
212+
213+ """ The follwing method takes in the reformatted stellar model files listed in
214+ the corresponding model-metallicity's cluster input files originally meant
215+ for color magnitude diagram makers. This is suppoesed to create the isochrone
216+ files out of which BPASS creates its isochrone objects. For each stellar model
217+ file, we try to find the row in the stellar model representing a model with the
218+ age closest to the inputted age and is within the specified margin error for
219+ log(age in years). Then that row is concattenated with other similar rows from
220+ different stellar model files to create the isochrone table. The age closest
221+ method may get rid of some data but is meant to improve accuracy of the age of the
222+ stellar model and improve computational efficiency to some extent."""
223+ def extractor (SinOrBin , model , age , metallicity , input_dir , bpass_evo_dir ,
224+ margin ):
225+ entries = glob .glob (input_dir + "/*" )
226+ print (len (entries ))
227+ stopLen = 4
228+ bigOne = None ;
229+ initlMass = np .nan
230+ initlMass2 = np .nan
231+ for x in entries :
232+ if x [0 ] != '.' and not os .path .isdir (x ):
233+ # To Obtain initial mass I will use Regex
234+ # to extract the initial mass of a model.
235+ try :
236+ org = Table .read (x , format = 'fits' )
237+ org = org .to_pandas ()
238+ starIndex = "Not Applicable"
239+ if (x [- 8 :- 5 ] == 'bin' and org ['model_type' ][0 ] == 0 ):
240+ starIndex = findLastRelMax (org )
241+ f = org [np .where (np .abs (np .log10 (org ['age' ])- age ) <= margin )[0 ]]
242+ if (len (f ) != 0 ):
243+ f = f [np .where (np .abs (f ['age' ] - 10 ** age ) ==
244+ np .min (np .abs (f ['age' ] - 10 ** age )))[0 ]]
245+ if len (f ) != 0 :
246+ index = f .index
247+ indexlen = len (index )
248+ if initial :
249+ initial = False
250+ bigOne = f
251+ # mass and mass2 stand for the initial masses
252+ # of the primary and secondary respectively.
253+ bigOne ['mass' ] = [initlMass ] * indexlen
254+ bigOne ['mass2' ] = [initlMass2 ] * indexlen
255+ bigOne ['single' ] = [True ] * indexlen
256+ bigOne ['mergered?' ] = [False ] * indexlen
257+ if (x [- 8 :- 5 ] == 'bin' or x [- 8 :- 5 ] == 'sec' ):
258+ bigOne ['single' ] = [False ] * indexlen
259+ # I only need to check one row since
260+ # whether a model is a merger model
261+ # is the same for all rows of the model
262+ if bigOne ['Merger_Related' ][0 ]:
263+ # Find the row in which the model merges.
264+ # Note that I pass in the original grid.
265+ merge_pt = find_mergers (org )
266+ bigOne ['mergered?' ] = [x > merge for x in bigOne .index ]
267+ print (bigOne .shape )
268+ else :
269+ f ['mass' ] = [initlMass ] * indexlen
270+ f ['mass2' ] = [initlMass2 ] * indexlen
271+ f ['single' ] = [True ] * indexlen
272+ f ['Merger Related' ]= [False ] * indexlen
273+ if (x [- 8 :- 5 ] == 'bin' ):
274+ f ['single' ] = [False ] * indexlen
275+ if bigOne ['Merger_Related' ][0 ]:
276+ merge_pt = find_mergers (org )
277+ bigOne ['mergered?' ] = [x > merge for x
278+ in bigOne .index ]
279+ bigOne = pd .concat ([f , bigOne ], axis = 0 )
280+ print (bigOne .shape )
281+ except FileNotFoundError :
282+ print ('File Not Found/Readable' )
283+
284+ try :
285+ if not isinstance (bigOne , type (None )) and not (bigOne ['age' ].empty ):
286+ bigOne = bigOne .apply (pd .to_numeric , errors = 'coerce' )
287+ reduced = Table .from_pandas (bigOne )
288+ if not os .path .isdir (bpass_evo_dir +
289+ 'iso/' + model + '/' + metallicity + '/' ):
290+ os .makedirs (bpass_evo_dir +
291+ 'iso/' + model + '/' + metallicity + '/' )
292+ print (reduced .columns )
293+ reduced .write (bpass_evo_dir + 'iso/' + model + "/" + metallicity +
294+ '/' + 'iso' + str (age ) + SinOrBin +
295+ '.fits' , format = 'fits' , overwrite = True )
296+ print (bigOne .shape )
297+ except IndexError :
298+ print ('It looks like there are no stars in out input file' +
299+ ' that have the specified age...' )
0 commit comments