# John F. Raffensperger, 2019. This code changes the input to Hector,
# with a pulse increase in each chemical, one chemical at a time,
# and then records the change in temperature.
from math import log10, floor
import os
import subprocess


def roundToSignificantFigures(valueToRound, sigFigures):
	return valueToRound
	if valueToRound == 0.0:
		return 0.0
	else: return round(valueToRound, sigFigures - int(floor(log10(abs(valueToRound)))) - 1)
	# return float('{:g}'.format(float('{:.{p}g}'.format(valueToRound, p=sigFigures))))
# print (roundToSignificantFigures(1.23456e-6, 5))
# -----------------------------------------------
def increaseEmission(inputfilename, outputfilename, chemical):
	with open (inputfilename, 'r') as inputfile, open (outputfilename, 'w') as outputfile:
		counter = 0
		previousline = []
		for line in inputfile:
			linelist = line.split(',')
			if counter == 3:
				if chemical not in linelist:
					print("Fail with chemical "+ chemical)
					return 'fail', 0.0
				else:
					chemindex = linelist.index(chemical)
					units = previousline[chemindex]
			if linelist[0] == '2020':
				oldvalue = linelist[chemindex]
				# pulse = roundToSignificantFigures(1.5*float(linelist[chemindex]),8)
				pulse = roundToSignificantFigures(15.0*float(linelist[chemindex]),9)
				linelist[chemindex] = str(pulse)
				print ("Changed " + chemical + " from " + oldvalue + " to " + linelist[chemindex])
				outputfile.write(','.join(linelist))
			else: outputfile.write(line)
			counter += 1
			previousline = linelist
	return units, pulse
# -----------------------------------------------
def readTemperatureOutput(outputfilename):
	yearlist = []
	templist = []
	with open (outputfilename, 'r') as youroutputfile:
		for line in youroutputfile:
			fields = line.split(',')
			if len(fields) >= 5: # temperature global average
				if fields[4] == 'Tgav':
					# print (fields[0],fields[5])
					yearlist.append(int(fields[0]))
					templist.append(float(fields[5]))
	return  yearlist, templist
# -----------------------------------------------
def savePulseTemperature(chemical, pulse, units, pulselist, outputfilename):
	with open (outputfilename, 'a+') as outputfile:
		outputfile.write(chemical + ',' + str(pulse) + ',' + units + ',' + ','.join([str(p) for p in pulselist]) + '\n')
# -----------------------------------------------
# pathname = "D:/Work documents/2 Research/Global warming/Numerical Simulation/hector-2.0.1-Windows/"
pathname = "C:/Users/johnr/Documents/Work documents/2 Research/Global warming/Numerical Simulation/hector-2.0.1-Windows/"
os.chdir(pathname)

# 1. Read original temperatures.
originalyearlist, originaltemplist = readTemperatureOutput(pathname + "output/sample_outputstream_rcp26.csv")

# 2. Create pulse.
chemicals = ['ffi_emissions', 'CH4_emissions', 'N2O_emissions', 'CF4_emissions', 'C2F6_emissions', 'HFC125_emissions', 'HFC134a_emissions', 'HFC143a_emissions', 'SF6_emissions']
for chemical in chemicals:
	units, pulse = increaseEmission(pathname + 'input/emissions/RCP26_emissions - backup.csv', pathname + 'input/emissions/RCP26_emissions.csv', chemical)
	if pulse > 0.0:
		# run Hector.
		returnvalue = subprocess.call(pathname + 'run_hector', shell=True)

		# 3. Read new temperature.
		newyearlist, newtemplist = readTemperatureOutput(pathname + "output/outputstream_rcp26.csv")
		pulselist = []
		for i in range(len(originaltemplist)):
			if newyearlist[i] >= 2020:
				# pulselist.append(roundToSignificantFigures (newtemplist[i] - originaltemplist[i], 8))
				pulselist.append(roundToSignificantFigures (newtemplist[i] - originaltemplist[i], 8))
		print ("Total pulse effect = %f" % sum(pulselist))
		if abs(sum(pulselist)) > 0.0005:
			savePulseTemperature(chemical, pulse, units, pulselist, 'pulseoutput_bigchange.txt')
			print ("Saved")
		else: print("Chemical " + chemical + " was ignored.")
	# else: print ("Died with chemical " + chemical)
print ("Done")