#!/usr/bin/env python
import matplotlib.pyplot as plt
import matplotlib as mpl
import numpy as np

def baule(E, mass):
	mu = (12. + 1.008 + 3 * 2.014) / mass

	ET = E * 2.4 * mu / (1 + mu)**2
	stderr = ( ET / 2.4 * 4. - ET )

	return ET

#Ei, Et, standard deviation
pd_LO = np.array(
[[89.2, 25.8],
[102.4, 29.7],
[111.8, 31.6],
[120.0, 33.3]]
)

pd_v1 = np.array(
[[89.2, 24.2],
[97.4, 29.2],
[102.4, 30.6],
[111.8, 31.6],
[120.0, 33.4]]
)

pt_LO = np.array(
[[81.7, 14.574],
[89.3, 16.103],
[97.4, 16.683],
[102.5, 18.053],
[111.9, 19.651],
[120.1, 21.767]]
)

pt_v1 = np.array(
[[60.7, 10.968],
[71.4, 13.128],
[81.9, 14.172],
[92.2, 15.846],
[104.6, 16.161]]
)

#Ei, Et
ni_LO = np.array(
[[101.1, 27.3249792731],
[112.3, 29.290332734827608],
[121.2, 32.59538852616818],
[130.7, 34.5546815954],
[136.4, 37.59329017442232],
[146.6, 40.94244770219179],
[160.4, 43.8412167303]]
)

ni_v1 = np.array(
[[101.1, 25.1936192292],
[112.3, 29.8106900874],
[121.2, 32.1593831303],
[130.7, 35.05502347494413],
[136.4, 39.0269222599],
[146.6, 39.5114032316],
[160.4, 43.6122275741]]
)

# Careful, the following is not so exactly an average and is taken from different states
cu = np.array(
[[83.4, 33.17198969],
[97.0, 38.6089964288],
[111.6, 44.4981275495],
[127.3, 51.1326663715],
[143.9, 58.0248047304],
[160.3, 64.7952813537],
[170.8, 69.1141081705],
[181.3, 73.5616236633]]
)
mpl.rc("text", usetex=True)
fig = plt.figure(figsize=(3.69,3))

plt.plot( ni_LO[:,0], ni_LO[:,1], marker='s', color='b', linestyle='--', label='Ni(111)' )
plt.plot( ni_v1[:,0], ni_v1[:,1], marker='s', color='b', linestyle=':', fillstyle='none' )
plt.plot( [10., 200.], [baule(10., 58.690), baule(200., 58.690)], color='b', linestyle='-' )
plt.plot( pd_LO[:,0], pd_LO[:,1], marker='o', label='Pd(111)', color='k', linestyle='--' )
plt.plot( pd_v1[:,0], pd_v1[:,1], marker='o', color='k', fillstyle='none', linestyle=':' )
plt.plot( [10., 200.], [baule(10., 106.42), baule(200., 106.42)], color='black', linestyle='-' )
plt.plot( pt_LO[:,0], pt_LO[:,1], marker='^', label='Pt(111)', color='r', linestyle='--' )
plt.plot( pt_v1[:,0], pt_v1[:,1], marker='^', color='red', linestyle=':', fillstyle='none' )
plt.plot( [10., 200.], [baule(10., 195.084), baule(200., 195.084)], color='r', linestyle='-' )
plt.plot( cu[:,0], cu[:,1], marker='d', color='g', linestyle=':', label='Cu(111)', fillstyle='none' )
plt.plot( [10., 200.], [baule(10., 63.546), baule(200., 63.546)], color='g', linestyle='-' )

plt.legend(loc='upper left', numpoints=1, scatterpoints=1)
plt.xlim(50,170)
plt.ylim(0,140)
plt.tick_params(length=6, width=1)
plt.ylabel('Energy transfer (kJ/mol)')
plt.xlabel('Incidence energy (kJ/mol)')

plt.tight_layout()
plt.savefig('energytransfer_comparison.pdf')
#plt.show()