#!/usr/bin/env python
import sys
from os import path, getcwd, popen
from os import chdir, system
import numpy as np

#folders = [ '0.94_ms_v0j8', '1.18_ms_v0j8', '1.29_ms_v0j8', '1.55_ms_v0j8', '1.80_ms_v0j8', '2.12_ms_v0j8', '2.56_ms_v0j8' ]
folders = [ '2.56_ms_v1j0', '2.56_ms_v1j2', '2.56_ms_v1j4', '2.56_ms_v1j6', '2.56_ms_v1j8' ] 

here=getcwd()
def collectangles(folder):
	chdir(folder)

	Ntraj = int( popen('ls | egrep "^[0-9].....$" | wc -l').read() )
	angles = []
	angles_scat = []
	angles_trap = []
	here_tmp=getcwd()				#workdir
	for i in range(1,(Ntraj+1)):
		traj= '%.6d' % i		# define working directory
		newdir=path.join(here_tmp, traj)	# updating mypath

		#print( newdir )
		if not path.exists(newdir): 
			print( "folder", newdir, "not found. exiting..." )
			quit()

		chdir(newdir)

		if not path.exists('outcome') or not path.exists('PostAnalysis.dat'):
			continue

		outcome = open('outcome', 'r').readline()
		post_dat = open('PostAnalysis.dat', 'r').readlines()
		if 'REACTION ' in outcome:
			angles.append( [float(post_dat[25].split()[0]), float(post_dat[27].split()[0]),	float(post_dat[25].split()[2]), float(post_dat[27].split()[2])] )
		elif 'SCATTERING' in outcome:
			angles_scat.append( [float(post_dat[24].split()[0]), float(post_dat[26].split()[0])] )
		elif 'TRAPPING' in outcome:
			angles_trap.append( [float(post_dat[24].split()[0]), float(post_dat[26].split()[0])] )
	
	bins = np.arange(2.5,190.,15.)
	if len(angles) != 0:
		angles = np.array(angles)
		nreac0, tmp = np.histogram( angles[:,0], bins=bins, density=False ) # Reacted t=0
		nreac1, tmp = np.histogram( angles[:,1], bins=bins, density=False ) # Reacted r=2.2

	if len(angles_scat) != 0:
		angles_scat = np.array(angles_scat)
		nscat, tmp = np.histogram( angles_scat[:,0], bins=bins, density=False ) # Scattered t=0

	if len(angles_trap) != 0:
		angles_trap = np.array(angles_trap)
		ntrap, tmp = np.histogram( angles_trap[:,0], bins=bins, density=False )
	else:
		ntrap = np.zeros( len(nreac0) )
	
	S0 = nreac0 / ( nreac0 + nscat + ntrap )
	error = np.sqrt( S0 * (1-S0) / ( nreac0 + nscat + ntrap ) )
	#print( '{:s} {:5.4f} {:5.4f}'.format( folder, *S0, *error ) )
	sys.stdout.write('{:s},'.format(folder))
	for i in range( len( S0 ) ):
		sys.stdout.write(' {:5.4f},'.format(S0[i]))
	for i in range( len( error ) ):
		sys.stdout.write(' {:5.4f},'.format(error[i]))
	sys.stdout.write('\n')
	#print(tmp)

	chdir(here)

for i in range(len(folders)):
	collectangles(folders[i])
