import numpy as np

tunit = 6.88027e+14
tyr = tunit/3.15e7
currenttime = 9.9284345716569
redshift = 17.91
omega_b = 0.0449
h = 0.71
baryonmass = 1.9e-29*h**2*omega_b*((8*3.08e24)**3*0.2745822715872526)/2e33
mass_solar = 1.998e33
L_mass = 1.12e46
L_ed = 1.3e38 * 6.24e11 # ev
epsilon = 0.2
T_ed = epsilon*mass_solar/1.44e17
e_ph0 =  10.2   # 4.202e3  # 1.2 ev Temperature, median wavelength
lifetime_max = 2.0e7*3.157e7 # 20 Myr
e_ph =  e_ph0 #  e_ph0*(10.0/(40.0))**0.25   #e_ph0*(10.0/(mass_max+mass_min)*2)**0.25
mass_min = 5

#f=open('PopIII_list_0002_type1.txt','r')
#f1=open('PopIII_list_0002_type5.txt','r')

fo = open('PopIII_LW_%5.2f.txt'%redshift,'w')

source_number = 0

haloposition=np.zeros(3)

for  file in ['PopIII_list_0011_type1.txt','PopIII_list_0011_type5.txt']:
 f =  open(file,'r')
 line = f.readline()
 line = f.readline()
 while line !="":
   starmass = float(line.split()[1])
   if starmass > 1.0e20: starmass=starmass/1.0e20 
   if starmass > mass_min:
        ct = float(line.split()[2])
        logm = np.log10(starmass)
        lifetime = 10**(9.785 - 3.759*logm + 1.413*logm*logm - \
                              0.186*logm*logm*logm)*3.1557e7
        if (currenttime - ct) < lifetime/tunit:
          haloposition[0] = float(line.split()[3])
          haloposition[1] = float(line.split()[4])
          haloposition[2] = float(line.split()[5])
          energyphoton=e_ph      
          halolifetime = lifetime/tunit
          luminosity = 10**(44.03 + 4.59*logm  - 0.77*logm*logm)
          halophotonenergy = e_ph0     
          halocreationtime = ct
          fo.write(str("PhotonTestSourceType[%d]       = 1"%source_number+"\n"))
          fo.write(str("PhotonTestSourcePosition[%d]   = %f %f %f"%(source_number, \
             haloposition[0],haloposition[1],haloposition[2])+"\n"))
          fo.write(str("PhotonTestSourceLuminosity[%d] = %e"%(source_number,luminosity)+"\n"))
          fo.write(str("PhotonTestSourceLifeTime[%d]   = %f"%(source_number,halolifetime)+"\n"))
          fo.write(str("PhotonTestSourceCreationTime[%d] = %f"%(source_number,halocreationtime)+"\n"))
          fo.write(str("PhotonTestSourceEnergyBins[%d] = 1"%source_number+"\n")) 
          fo.write(str("PhotonTestSourceEnergy[%d] = %f"%(source_number,halophotonenergy)+"\n\n"))     
          source_number = source_number + 1
   line = f.readline()
 f.close()

#f.close()
#f1.close()
fo.close()