Extract average cross-sections for both polarization states in parser (b93429ed) · Commits · Giacomo Mulas / NP_TMcode

src/scripts/parse_output.py

+87 −31

Original line number	Diff line number	Diff line
		@@ -133,25 +133,31 @@ def parse_legacy_oclu(config):
		errors = 0
		oclu_name = config['result_file']
		root_name = config['output_name']
		out20 = None # Cluster average cross-sections
		out21 = None # Cluster average cross-sections in state -1
		out22 = None # Cluster average cross-sections in state +1
		out31 = None # Cluster differential cross-sections in state -1
		out32 = None # Cluster differential cross-sections in state +1
		out40 = None # Cluster integrated asymmetry parameter and radiation pressure
		out41 = None # Cluster integrated asymmetry parameter and radiation pressure in state -1
		out42 = None # Cluster integrated asymmetry parameter and radiation pressure in state +1
		out51 = None # Cluster differential asymmetry parameter and radiation pressure forces in state -1
		out52 = None # Cluster differential asymmetry parameter and radiation pressure forces in state +1
		oclu_file = open(oclu_name, "r") # open the OCLU file for reading
		file_line = "first line" # a string to parse the OCLU file lines
		if ('ALL' in config['selection'] or 'ICS' in config['selection']):
		out20 = open(root_name + "_ics.csv", "w") # open a file for integrated cross sections
		out20.write("Wavelength_m,ScaSec_m2,AbsSec_m2,ExtSec_m2\n")
		out21 = open(root_name + "_ics1.csv", "w") # open a file for integrated cross sections in state -1
		out22 = open(root_name + "_ics2.csv", "w") # open a file for integrated cross sections in state +1
		out21.write("Wavelength_m,ScaSec_m2,AbsSec_m2,ExtSec_m2\n")
		out22.write("Wavelength_m,ScaSec_m2,AbsSec_m2,ExtSec_m2\n")
		if ('ALL' in config['selection'] or 'DCS' in config['selection']):
		out31 = open(root_name + "_dcs1.csv", "w") # open a file for differential cross-sections in state -1
		out31.write("Wavelength_m,THi_deg,THs_deg,PHi_deg,PHs_deg,ScaSec_m2,AbsSec_m2,ExtSec_m2\n")
		out32 = open(root_name + "_dcs2.csv", "w") # open a file for differential cross-sections in state +1
		out32.write("Wavelength_m,THi_deg,THs_deg,PHi_deg,PHs_deg,ScaSec_m2,AbsSec_m2,ExtSec_m2\n")
		if ('ALL' in config['selection'] or 'IRP' in config['selection']):
		out40 = open(root_name + "_irp.csv", "w") # open a file for integrated radiation pressure forces
		out40.write("Wavelength_m,CosAv,RaPr_m2\n")
		out41 = open(root_name + "_irp1.csv", "w") # open a file for integrated radiation pressure forces in state -1
		out42 = open(root_name + "_irp2.csv", "w") # open a file for integrated radiation pressure forces in state +1
		out41.write("Wavelength_m,CosAv,RaPr_m2\n")
		out42.write("Wavelength_m,CosAv,RaPr_m2\n")
		if ('ALL' in config['selection'] or 'DRP' in config['selection']):
		out51 = open(root_name + "_drp1.csv", "w") # open a file for differential radiation pressure forces in state -1
		out51.write("Wavelength_m,THi_deg,THs_deg,PHi_deg,PHs_deg,CosAv,RaPr_m2,Fl_m2,Fr_m2,Fk_m2,Fx_m2,Fy_m2,Fz_m2,TQEl_m2,TQEr_m2,TQEk_m2,TQEx_m2,TQEy_m2,TQEz_m2,TQSl_m2,TQSr_m2,TQSk_m2,TQSx_m2,TQSy_m2,TQSz_m2\n")
		@@ -213,13 +219,13 @@ def parse_legacy_oclu(config):
		alam = 2.0 * math.pi / vk
		if ("CLUSTER (ENSEMBLE AVERAGE, " in file_line):
		# we are in average section. We know the average cross-sections
		# are after 3 more lines.
		# are after 3 more lines for state -1.
		for i in range(3):
		file_line = oclu_file.readline()
		# the following check is needed to parse C++ output
		if ("INSERTION" in file_line):
		file_line = oclu_file.readline()
		# we know we are in LIN -1 because it is the first one
		# we know we are in state -1 because it is the first one
		# we also know that the next line contains the values
		# we are looking for, so we read it
		scasm = float(file_line[1:15].replace("D", "E"))
		@@ -228,7 +234,7 @@ def parse_legacy_oclu(config):
		# we can write the average values, similarly to fort.22
		# note that \n puts a new line at the end of the string
		output_line = "{0:.7E},{1:.7E},{2:.7E},{3:.7E}\n".format(alam, scasm, abssm, extsm)
		if (out20 is not None): out20.write(output_line)
		if (out21 is not None): out21.write(output_line)
		# we know that the asymmetry parameter and the radiation
		# pressure forces will be after 8 more lines
		for i in range(8):
		@@ -236,7 +242,28 @@ def parse_legacy_oclu(config):
		cosav = float(file_line[8:23].replace("D", "E"))
		rapr = float(file_line[31:46].replace("D", "E"))
		output_line = "{0:.7E},{1:.7E},{2:.7E}\n".format(alam, cosav, rapr)
		if (out40 is not None): out40.write(output_line)
		if (out41 is not None): out41.write(output_line)
		# now we look for state +1 cross-sections, after 4 lines
		for i in range(4):
		file_line = oclu_file.readline()
		# the following check is needed to parse C++ output
		if ("INSERTION" in file_line):
		file_line = oclu_file.readline()
		scasm = float(file_line[1:15].replace("D", "E"))
		abssm = float(file_line[17:30].replace("D", "E"))
		extsm = float(file_line[32:45].replace("D", "E"))
		# we can write the average values, similarly to fort.22
		# note that \n puts a new line at the end of the string
		output_line = "{0:.7E},{1:.7E},{2:.7E},{3:.7E}\n".format(alam, scasm, abssm, extsm)
		if (out22 is not None): out22.write(output_line)
		# we know that the asymmetry parameter and the radiation
		# pressure forces will be after 8 more lines
		for i in range(8):
		file_line = oclu_file.readline()
		cosav = float(file_line[8:23].replace("D", "E"))
		rapr = float(file_line[31:46].replace("D", "E"))
		output_line = "{0:.7E},{1:.7E},{2:.7E}\n".format(alam, cosav, rapr)
		if (out42 is not None): out42.write(output_line)
		# the averages were written. We look for differential section
		for di in range(ndirs):
		while ("JTH =" not in file_line):
		@@ -369,10 +396,12 @@ def parse_legacy_oclu(config):
		# closes di for loop
		# The parsing loop ends here

		if (out20 is not None): out20.close()
		if (out21 is not None): out21.close()
		if (out22 is not None): out22.close()
		if (out31 is not None): out31.close()
		if (out32 is not None): out32.close()
		if (out40 is not None): out40.close()
		if (out41 is not None): out41.close()
		if (out42 is not None): out42.close()
		if (out51 is not None): out51.close()
		if (out52 is not None): out52.close()
		oclu_file.close() # close the OCLU file
		@@ -386,31 +415,37 @@ def parse_legacy_oinclu(config):
		errors = 0
		oclu_name = config['result_file']
		root_name = config['output_name']
		out20 = None # Inclusion average cross-sections
		out21 = None # Inclusion average cross-sections in state -1
		out22 = None # Inclusion average cross-sections in state +1
		out31 = None # Inclusion differential cross-sections in state -1
		out32 = None # Inclusion differential cross-sections in state +1
		out40 = None # Inclusion integrated asymmetry parameter and radiation pressure
		out41 = None # Inclusion integrated asymmetry parameter and radiation pressure in state -1
		out42 = None # Inclusion integrated asymmetry parameter and radiation pressure in state +1
		out51 = None # Inclusion differential asymmetry parameter and radiation pressure forces in state -1
		out52 = None # Inclusion differential asymmetry parameter and radiation pressure forces in state +1
		try:
		oclu_file = open(oclu_name, "r") # open the OINCLU file for reading
		file_line = "first line" # a string to parse the OINCLU file lines
		if ('ALL' in config['selection'] or 'ICS' in config['selection']):
		out20 = open(root_name + "_ics.csv", "w") # open a file for integrated cross sections
		out20.write("Wavelength,ScaSec,AbsSec,ExtSec\n")
		out21 = open(root_name + "_ics1.csv", "w") # open a file for integrated cross sections in state -1
		out21.write("Wavelength_m,ScaSec_m2,AbsSec_m2,ExtSec_m2\n")
		out22 = open(root_name + "_ics2.csv", "w") # open a file for integrated cross sections in state +1
		out22.write("Wavelength_m,ScaSec_m2,AbsSec_m2,ExtSec_m2\n")
		if ('ALL' in config['selection'] or 'DCS' in config['selection']):
		out31 = open(root_name + "_dcs1.csv", "w") # open a file for differential cross-sections in state -1
		out31.write("Wavelength,THi,THs,PHi,PHs,ScaSec,AbsSec,ExtSec\n")
		out31.write("Wavelength_m,THi_deg,THs_deg,PHi_deg,PHs_deg,ScaSec_m2,AbsSec_m2,ExtSec_m2\n")
		out32 = open(root_name + "_dcs2.csv", "w") # open a file for differential cross-sections in state +1
		out32.write("Wavelength,THi,THs,PHi,PHs,ScaSec,AbsSec,ExtSec\n")
		out32.write("Wavelength_m,THi_deg,THs_deg,PHi_deg,PHs_deg,ScaSec_m2,AbsSec_m2,ExtSec_m2\n")
		if ('ALL' in config['selection'] or 'IRP' in config['selection']):
		out40 = open(root_name + "_irp.csv", "w") # open a file for integrated radiation pressure forces
		out40.write("Wavelength,CosAv,RaPr\n")
		out41 = open(root_name + "_irp1.csv", "w") # open a file for integrated radiation pressure forces in state -1
		out41.write("Wavelength_m,CosAv,RaPr_m2\n")
		out42 = open(root_name + "_irp1.csv", "w") # open a file for integrated radiation pressure forces in state +1
		out42.write("Wavelength_m,CosAv,RaPr_m2\n")
		if ('ALL' in config['selection'] or 'DRP' in config['selection']):
		out51 = open(root_name + "_drp1.csv", "w") # open a file for differential radiation pressure forces in state -1
		out51.write("Wavelength,THi,THs,PHi,PHs,CosAv,RaPr,Fl,Fr,Fk,Fx,Fy,Fz,TQEl,TQEr,TQEk,TQEx,TQEy,TQEz,TQSl,TQSr,TQSk,TQSx,TQSy,TQSz\n")
		out51.write("Wavelength_m,THi_deg,THs_deg,PHi_deg,PHs_deg,CosAv,RaPr_m2,Fl_m2,Fr_m2,Fk_m2,Fx_m2,Fy_m2,Fz_m2,TQEl_m2,TQEr_m2,TQEk_m2,TQEx_m2,TQEy_m2,TQEz_m2,TQSl_m2,TQSr_m2,TQSk_m2,TQSx_m2,TQSy_m2,TQSz_m2\n")
		out52 = open(root_name + "_drp2.csv", "w") # open a file for differential radiation pressure forces in state +1
		out52.write("Wavelength,THi,THs,PHi,PHs,CosAv,RaPr,Fl,Fr,Fk,Fx,Fy,Fz,TQEl,TQEr,TQEk,TQEx,TQEy,TQEz,TQSl,TQSr,TQSk,TQSx,TQSy,TQSz\n")
		out52.write("Wavelength_m,THi_deg,THs_deg,PHi_deg,PHs_deg,CosAv,RaPr_m2,Fl_m2,Fr_m2,Fk_m2,Fx_m2,Fy_m2,Fz_m2,TQEl_m2,TQEr_m2,TQEk_m2,TQEx_m2,TQEy_m2,TQEz_m2,TQSl_m2,TQSr_m2,TQSk_m2,TQSx_m2,TQSy_m2,TQSz_m2\n")

		# Define the quantities that you need to extract
		alam = 0.0
		@@ -463,14 +498,14 @@ def parse_legacy_oinclu(config):
		# replacing FORTRAN's D with E
		vk = float(file_line[5:20].replace("D", "E"))
		alam = 2.0 * math.pi / vk
		if ("ENSEMBLE AVERAGE, MODE 0" in file_line):
		if ("ENSEMBLE AVERAGE, MODE" in file_line):
		# we are in average section. We start a nested loop to
		# extract the average values
		found_averages = False
		while (not found_averages):
		file_line = oclu_file.readline()
		if ("----- SCS ----- ABS ----- EXS ----- ALBEDS --" in file_line):
		# we know we are in LIN -1 because it is the first one
		# we know we are in state -1 because it is the first one
		# we also know that the next line contains the values
		# we are looking for, so we read it
		file_line = oclu_file.readline()
		@@ -481,7 +516,7 @@ def parse_legacy_oinclu(config):
		# we can write the average values, similarly to fort.22
		# note that \n puts a new line at the end of the string
		output_line = "{0:.7E},{1:.7E},{2:.7E},{3:.7E}\n".format(alam, scasm, abssm, extsm)
		if (out20 is not None): out20.write(output_line)
		if (out21 is not None): out21.write(output_line)
		# we know that the asymmetry parameter and the radiation
		# pressure forces will be after 5 more lines
		for i in range(5):
		@@ -489,7 +524,26 @@ def parse_legacy_oinclu(config):
		cosav = float(file_line[8:23].replace("D", "E"))
		rapr = float(file_line[31:46].replace("D", "E"))
		output_line = "{0:.7E},{1:.7E},{2:.7E}\n".format(alam, cosav, rapr)
		if (out40 is not None): out40.write(output_line)
		if (out41 is not None): out41.write(output_line)
		# we know that state +1 starts after 4 more lines
		for i in range(4):
		file_line = oclu_file.readline()
		# we now extract the values from string sections
		scasm = float(file_line[1:15].replace("D", "E"))
		abssm = float(file_line[17:30].replace("D", "E"))
		extsm = float(file_line[32:45].replace("D", "E"))
		# we can write the average values, similarly to fort.22
		# note that \n puts a new line at the end of the string
		output_line = "{0:.7E},{1:.7E},{2:.7E},{3:.7E}\n".format(alam, scasm, abssm, extsm)
		if (out22 is not None): out22.write(output_line)
		# we know that the asymmetry parameter and the radiation
		# pressure forces will be after 5 more lines
		for i in range(5):
		file_line = oclu_file.readline()
		cosav = float(file_line[8:23].replace("D", "E"))
		rapr = float(file_line[31:46].replace("D", "E"))
		output_line = "{0:.7E},{1:.7E},{2:.7E}\n".format(alam, cosav, rapr)
		if (out42 is not None): out42.write(output_line)
		found_averages = True # terminate the inner loop
		# the averages were written. We look for SINGLE SCATTERER section
		# using another inner loop
		@@ -609,10 +663,12 @@ def parse_legacy_oinclu(config):
		found_differentials = True # terminate the inner loop
		# The parsing loop ends here

		if (out20 is not None): out20.close()
		if (out21 is not None): out21.close()
		if (out22 is not None): out22.close()
		if (out31 is not None): out31.close()
		if (out32 is not None): out32.close()
		if (out40 is not None): out40.close()
		if (out41 is not None): out41.close()
		if (out42 is not None): out42.close()
		if (out51 is not None): out51.close()
		if (out52 is not None): out52.close()
		oclu_file.close() # close the OINCLU file
		@@ -637,13 +693,13 @@ def parse_legacy_osph(config):
		file_line = "first line" # a string to parse the OSPH file lines
		if ('ALL' in config['selection'] or 'ICS' in config['selection']):
		out20 = open(root_name + "_ics.csv", "w") # open a file for integrated cross sections
		out20.write("Wavelength,ScaSec,AbsSec,ExtSec\n")
		out20.write("Wavelength_m,ScaSec_m2,AbsSec_m2,ExtSec_m2\n")
		if ('ALL' in config['selection'] or 'IRP' in config['selection']):
		out30 = open(root_name + "_irp.csv", "w") # open a file for integrated radiation pressure forces
		out30.write("Wavelength,CosAv,RaPr\n")
		out30.write("Wavelength_m,CosAv,RaPr_m2\n")
		if ('ALL' in config['selection'] or 'DRP' in config['selection']):
		out40 = open(root_name + "_drp.csv", "w") # open a file for differential radiation pressure forces in state -1
		out40.write("Wavelength,THi,THs,PHi,PHs,Fx,Fy,Fz\n")
		out40.write("Wavelength_m,THi_deg,THs_deg,PHi_deg,PHs_deg,Fx_m2,Fy_m2,Fz_m2\n")

		# Define the quantities that you need to extract
		alam = 0.0