#!/usr/bin/env python3
import sys
from os import path, getcwd, popen
from os import chdir, system
import numpy as np
from scipy.stats import sem
import matplotlib.pyplot as plt
import matplotlib as mpl

#folders = ['0.007_fs_TN','0.008_fs_TN','0.012_fs_TN','0.021_fs_TN','0.030_fs_TN','0.044_fs_TN','0.053_fs_TN','0.069_fs_TN','0.089_fs_TN','0.108_fs_TN','0.125_fs_TN','0.141_fs_TN','0.144_fs_TN']
folders = ['0.007_fs_TN','0.030_fs_TN','0.144_fs_TN']

here=getcwd()
def collectdisp(folder):
	xy = []
	chdir(folder)

	Ntraj = int( popen('ls | egrep "^[0-9].....$" | wc -l').read() )
	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').readlines()
		if 'REACTION ' in outcome[0]:
			post = open('PostAnalysis.dat', 'r').readlines()
			xy.append( float( post[8].split()[0] ) )

	chdir(here)

	return xy

mpl.rc("text", usetex=True)
fig = plt.figure(figsize=(3.69,3.))

for i in range(len(folders)):
	disp = collectdisp(folders[i])
	avg = np.average( disp )
	err = sem( disp )

	print('{0:s} {1:4.1f} {2:4.1f}'.format( folders[i], avg, err ))

	plt.hist( disp, bins=np.linspace(0.,100.,25), density=True, histtype='step', label='{:s}'.format(folders[i][:5]) )

# Terrace length of (433)
plt.ylim( plt.ylim() )
plt.plot( [14.3231367, 14.3231367], plt.ylim(), color='k', ls='--' )

plt.legend(loc='best', numpoints=1, frameon=True, fontsize=8)
plt.xlim(0., 100.)
#plt.ylim( 0., plt.ylim()[1] )
#plt.ylim(0.0, 0.99)
#plt.ylim(1e-7, 1.)
#plt.yscale('log')
plt.tick_params(length=6, width=1, direction='in', top=True, right=True)
plt.xlabel(r'Displacement in $X$ (\r{A})')
plt.ylabel(r'Distribution (a.u.)')

plt.tight_layout()
plt.savefig('Xdisp_histo.pdf')