# grab all current files, voltage files
# from voltage file determine where cutoff steps are
# for each current file, find maximum activation, then for 2nd step fraction of that activation
# print voltage, delay time, fraction

from sys import argv
from glob import glob
import numpy as np
from scipy.optimize import curve_fit



prefix = 'RECOVERY'
OUTFILE=None
ACTTOG = 1
DIRECTION=1 # 1 = outward, 0 = inward
SKIP=int(0.005/0.0001) # for excluding capacitive current
#print SKIP
EREV=-90
if len(argv) > 1:
	for args in argv[1:]:
		temp = args.split('=')
		if len(temp) == 2:
			try: 
				exec(args)
			except:
				exec('%s =\'%s\'' % (temp[0],temp[1]))
		
files =  glob('*%s*.dat' %(prefix))
ACTTOG=int(ACTTOG)

def triexponential(t,tau1,tau2,tau3,amp1,amp2):
	out = amp1*(1.0-np.exp(-t/tau1))*np.exp(-t/tau2)+amp2*(1.0-np.exp(-t/tau1))*np.exp(-t/tau3)
	return out

def boltz(v,vhalf,vslope,A):
	out = A/(1.0+np.exp(-(v-vhalf)/vslope))
	return out

volts = []
times = []
peakcurrs = []
endcurrs = []
taus = []
for things in files:
	time, current, voltage = np.loadtxt(things,unpack=True,usecols=[0,1,2])
	dt = time[1]-time[0]
	size = len(voltage)
	edges = [i for i in range(size-1) if abs(voltage[i]-voltage[i+1]) > 0.0025] #checks when voltage steps at least 2.5 mV
	if len(edges) < 2:
		continue
	edges.append(size-1)
	print edges, current[edges[0]]
	quit()
	vval = voltage[int(edges[2+ACTTOG]/2+edges[3+ACTTOG]/2)]
	tval = 1e9
	for j in range(len(edges)-1):
		temp = edges[j+1] - edges[j]
		if  temp < tval:
			tval = temp
	volts.append(vval*1000.0)
	times.append(tval*dt*1000)
	delay = time[edges[2]]
	
	if DIRECTION:
		peak1 = np.max(current[edges[1]+SKIP:edges[2]])
		peak2 = np.max(current[edges[3]+SKIP:edges[4]])
	else:
		peak1 = np.min(current[edges[1]+SKIP:edges[2]])
		peak2 = np.min(current[edges[3]+SKIP:edges[4]])

	
	end1 = np.mean(current[edges[2]-11:edges[2]-1])
	end2 = np.mean(current[edges[4]-11:edges[4]-1])
	
	
	peakcurrs.append([peak1,peak2])
	endcurrs.append([end1,end2])
	
	print tval*dt*1000, peak1, peak2, end1, end2
	#print diffcur,difftime
		
		
	#print vval*1000, peak-end
	#for p in popt:
	#	print p,
	#print


nsteps = len(volts)

diffcur = []

for i in range(nsteps):
	
	maxcurr = peakcurrs[i][0]-endcurrs[i][1]
	postcurr = peakcurrs[i][1]-endcurrs[i][1]
	precurr = peakcurrs[i][1]-endcurrs[i][0]
	
	preact = precurr/maxcurr
	postact = postcurr/maxcurr
	
	maxrecov = 1.0-preact
	recovered = 1.0-postact
	
	frac_recovered = recovered/maxrecov
	diffcur.append(frac_recovered)
	
	
data = np.array([times, diffcur,volts])

data = np.transpose(data)

data = np.sort(data,axis=0)


print data

if OUTFILE != None:
	try:
		fp = open(OUTFILE,'w')
	except:
		print 'OUTFILE must be a string, ideally <something>.dat'
		quit()
	
	for j in range(nsteps):
		fp.write('%e  %e  %e\n' % (data[j][0],data[j][1],data[j][2]))
	fp.close()
	
else:
	for j in range(nsteps):
		print '%e  %e  %e' % (data[j][0],data[j][1],data[j][2])