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

# Q, Eb, Easymp, Zb
nn_cl = np.array(
[[ -0.2 ,  1.06471976512 ,  2.3342807907786294 ],
[ -0.15000000000000002 ,  1.05572422513 ,  2.3653861676098678 ],
[ -0.10000000000000003 ,  1.0483398713 ,  2.3977591411472354 ],
[ -0.050000000000000044 ,  1.04178219938 ,  2.4306237699106727 ],
[ -5.551115123125783e-17 ,  1.03533446147 ,  2.463924728780647 ],
[ 0.04999999999999993 ,  1.02811854612 ,  2.4972946663619195 ],
[ 0.09999999999999992 ,  1.01966941076 ,  2.5305398012055065 ],
[ 0.1499999999999999 ,  1.00971385307 ,  2.563432810470523 ],
[ 0.1999999999999999 ,  0.998111734417 ,  2.5962195717003964 ]]
)

nn_h = np.array(
[[ -0.2 ,  1.06674810183 ,  2.4074582239856874 ],
[ -0.15000000000000002 ,  1.06058254599 ,  2.4225150433392475 ],
[ -0.10000000000000003 ,  1.05373425647 ,  2.436984962691941 ],
[ -0.050000000000000044 ,  1.04552485342 ,  2.4505711485303414 ],
[ -5.551115123125783e-17 ,  1.03533418935 ,  2.463936196821342 ],
[ 0.04999999999999993 ,  1.02227571732 ,  2.476918650055769 ],
[ 0.09999999999999992 ,  1.00586072081 ,  2.4895157898987033 ],
[ 0.1499999999999999 ,  0.985717220198 ,  2.501987566438015 ],
[ 0.1999999999999999 ,  0.961469970458 ,  2.514050941963971 ]]
)

nn_hcl = np.array(
[[ -0.2 ,  1.11970173114 ,  2.2848390574691875 ],
[ -0.15000000000000002 ,  1.09338151801 ,  2.3278382908139803 ],
[ -0.10000000000000003 ,  1.07195554477 ,  2.372456732012188 ],
[ -0.050000000000000044 ,  1.05321533533 ,  2.4180944495413765 ],
[ -5.551115123125783e-17 ,  1.03533500569 ,  2.4638631261997546 ],
[ 0.04999999999999993 ,  1.01623397328 ,  2.5106966996861413 ],
[ 0.09999999999999992 ,  0.994830313376 ,  2.5575216589438776 ],
[ 0.1499999999999999 ,  0.9703958493 ,  2.6053063295476298 ],
[ 0.1999999999999999 ,  0.942333291125 ,  2.6540731273851095 ]]
)

nn_hcl_fcc = np.array(
[[ -0.2 ,  1.32965009152 ,  2.2366583846552683 ],
[ -0.15000000000000002 ,  1.32391202655 ,  2.284235310495274 ],
[ -0.10000000000000003 ,  1.3234162351 ,  2.330484053510158 ],
[ -0.050000000000000044 ,  1.32600648693 ,  2.3749412332899387 ],
[ -5.551115123125783e-17 ,  1.33001907792 ,  2.418113257913953 ],
[ 0.04999999999999993 ,  1.33358948387 ,  2.459704679923133 ],
[ 0.09999999999999992 ,  1.33587170282 ,  2.4995492355792104 ],
[ 0.15000000000000002 ,  1.33635280013 ,  2.5373951767007337 ],
[ 0.1999999999999999 ,  1.33473372266 ,  2.5747548541000365 ]]
)

nn_hcl_bridge = np.array(
[[ -0.2 ,  1.15155157002 ,  2.150458853581671 ],
[ -0.15000000000000002 ,  1.13905392465 ,  2.209700617451087 ],
[ -0.10000000000000003 ,  1.13252373623 ,  2.268755786479021 ],
[ -0.050000000000000044 ,  1.1286455695 ,  2.3255772770023775 ],
[ 0.0 ,  1.12481012464 ,  2.3803905616009806 ],
[ 0.04999999999999993 ,  1.11830497069 ,  2.4336517885503532 ],
[ 0.1 ,  1.10791239815 ,  2.4846915745453226 ],
[ 0.15000000000000002 ,  1.09306668228 ,  2.533231708431221 ],
[ 0.20000000000000004 ,  1.07357843184 ,  2.5791859080852984 ]]
)

nn_hcl_hcp = np.array(
[[ -0.2 ,  1.27724341107 ,  2.2552012245246726 ],
[ -0.15000000000000002 ,  1.28834701714 ,  2.2992735273785985 ],
[ -0.10000000000000003 ,  1.30108303425 ,  2.3415590436252938 ],
[ -0.050000000000000044 ,  1.31428137281 ,  2.3811800412611994 ],
[ -5.551115123125783e-17 ,  1.32700405634 ,  2.4194666952677464 ],
[ 0.04999999999999993 ,  1.33805296752 ,  2.45642987213921 ],
[ 0.09999999999999992 ,  1.34708224963 ,  2.492286698251086 ],
[ 0.1499999999999999 ,  1.35403911258 ,  2.5267105176920426 ],
[ 0.1999999999999999 ,  1.35895103987 ,  2.5596687588564437 ]]
)

nn_hcl_top = np.array(
[[ -0.2 ,  1.1992686402 ,  2.398993591967966 ],
[ -0.15000000000000002 ,  1.19261001404 ,  2.449773465855162 ],
[ -0.1 ,  1.19207367763 ,  2.4985596054669936 ],
[ -0.05 ,  1.19419779842 ,  2.5464029792317424 ],
[ 0.0 ,  1.19658641388 ,  2.59301275833329 ],
[ 0.05 ,  1.19695621662 ,  2.6388382863465454 ],
[ 0.10000000000000002 ,  1.19467998418 ,  2.683719001467474 ],
[ 0.15 ,  1.18989758309 ,  2.7281846478884377 ],
[ 0.2 ,  1.18318480627 ,  2.771490786969015 ]]
)

nn_cl[:,1] *= 96.4853075
nn_h[:,1] *= 96.4853075
nn_hcl[:,1] *= 96.4853075

nn_hcl_fcc[:,1] *= 96.4853075
nn_hcl_bridge[:,1] *= 96.4853075
nn_hcl_hcp[:,1] *= 96.4853075
nn_hcl_top[:,1] *= 96.4853075

mpl.rc("font", size=12)
fig, ax = plt.subplots(nrows=1, ncols=2, figsize=(3.69,3))

ax2 = plt.subplot(1,2,1)

plt.scatter( nn_cl[:,0], nn_cl[:,1], label='Cl', color='b', marker='o' )
slope, intercept, r_value, p_value, std_err = linregress( nn_cl[:,0], nn_cl[:,1] )
plt.plot( [-1., 1.], [-slope + intercept, slope + intercept], color='b' )
plt.annotate('{0:2.2f}'.format( slope ), xy=(0.05, 100), size=8, color='b')

plt.scatter( nn_h[:,0], nn_h[:,1], label='H', color='r', marker='o' )
slope, intercept, r_value, p_value, std_err = linregress( nn_h[:,0], nn_h[:,1] )
plt.plot( [-1., 1.], [-slope + intercept, slope + intercept], color='r' )
plt.annotate('{0:2.2f}'.format( slope ), xy=(-0.2, 99), size=8, color='r')

plt.scatter( nn_hcl[:,0], nn_hcl[:,1], label='HCl', color='orange', marker='o' )
slope, intercept, r_value, p_value, std_err = linregress( nn_hcl[:,0], nn_hcl[:,1] )
plt.plot( [-1., 1.], [-slope + intercept, slope + intercept], color='orange' )
plt.annotate('{0:2.2f}'.format( slope ), xy=(-0.18, 107), size=8, color='orange')

plt.annotate(r'$\beta$', xy=(-0.2, 92), size=12)

plt.legend(loc='upper right', numpoints=1, frameon=False, borderaxespad=0, handletextpad=0.)
plt.xlim(-0.25,0.25)
plt.ylim(90., 110.)
plt.tick_params(length=6, width=1, direction='in', top=True, right=True)
plt.tick_params(axis='x', labelsize=9)
#ax2.set_xticks([-0.3, -0.1, 0.1, 0.3])
plt.ylabel('$E_b$ (kJ/mol)')
plt.xlabel(r'Q ($\mathrm{\AA}$)')
#plt.annotate('(a)', xy=(0.15, 175.), size=12)

ax3 = plt.subplot(1,2,2)

plt.scatter( nn_cl[:,0], nn_cl[:,2], label='Cl', color='blue', marker='o' )
slope, intercept, r_value, p_value, std_err = linregress( nn_cl[:,0], nn_cl[:,2] )
plt.plot( [-1., 1.], [-slope + intercept, slope + intercept], color='b' )
plt.annotate('{0:2.2f}'.format( slope ), xy=(-0.1, 2.5), size=8, color='b')

plt.scatter( nn_h[:,0], nn_h[:,2], label='H', color='r', marker='o' )
slope, intercept, r_value, p_value, std_err = linregress( nn_h[:,0], nn_h[:,2] )
plt.plot( [-1., 1.], [-slope + intercept, slope + intercept], color='r' )
plt.annotate('{0:2.2f}'.format( slope ), xy=(0.1, 2.45), size=8, color='r')

plt.scatter( nn_hcl[:,0], nn_hcl[:,2], label='HCl', color='orange', marker='o' )
slope, intercept, r_value, p_value, std_err = linregress( nn_hcl[:,0], nn_hcl[:,2] )
plt.plot( [-1., 1.], [-slope + intercept, slope + intercept], color='orange' )
plt.annotate('{0:2.2f}'.format( slope ), xy=(0.05, 2.65), size=8, color='orange')

plt.annotate(r'$\alpha$', xy=(0.15, 2.25), size=12)

plt.xlim(-0.25,0.25)
plt.ylim(2.2,2.7)
plt.tick_params(length=6, width=1, direction='in', top=True, right=True)
plt.tick_params(axis='y', labelleft=False, labelright=True)
plt.tick_params(axis='x', labelsize=9)
plt.ylabel(r'$Z_{Cl}$ ($\mathrm{\AA}$)')
ax3.yaxis.set_label_position("right")
#ax3.set_xticks([-0.3, -0.1, 0.1, 0.3])
plt.xlabel(r'Q ($\mathrm{\AA}$)')
#plt.annotate('(b)', xy=(-0.25, 2.28), size=12)

plt.tight_layout()
plt.subplots_adjust(hspace = 0.0, wspace = 0.08)
plt.subplots_adjust(left=0.18, right=0.84)
plt.savefig('coupling.pdf')
#plt.show()