from matplotlib.pylab import *
import numpy
import sys
import os

fig_size = [800/72.27 ,700/72.27]
params = {'backend': 'ps', 'axes.labelsize': 25, 'font.size': 25, 'legend.fontsize': 14,
          'xtick.labelsize': 25, 'ytick.labelsize': 25, 
          'xtick.major.size': 15,'ytick.major.size': 15,
          'xtick.minor.size': 8,'ytick.minor.size': 8,
          'text.usetex': False, 'figure.figsize': fig_size, 'xtick.major.pad': 11,
          'figure.subplot.top': 0.94,'figure.subplot.left': 0.105,'figure.subplot.bottom': 0.095,'figure.subplot.right': 0.978}

rcParams.update(params)
a=subplot(111)

title_ = "First-Principles Equation of State Database";
if (len(sys.argv)>1):
   title_ = str(sys.argv[1])
title(title_)

scriptName = sys.argv[0][:-3]; # [:3] removes ".py" for file name 
try:
   os.remove(scriptName+'.pdf')
   os.remove(scriptName+'.png')
except:
   pass

d  = numpy.loadtxt('FPEOS_isochores_ideal_Debye.txt',usecols=(1,3,7,25,27,29,31)); # index,rho,T,1-P/PIdeal,1-E/EIdeal,1-PDebye/PIdeal,1-EDebye/EIdeal
dp = numpy.loadtxt('FPEOS_isochor_points.txt',     usecols=(1,3,7,25,27,29,31)); # index,rho,T,1-P/PIdeal,1-E/EIdeal,1-PDebye/PIdeal,1-EDebye/EIdeal

iMin = int(min(d[:,0]))
iMax = int(max(d[:,0]))
#print(iMin,iMax)
#iMin = min(dp[:,0])
#iMax = max(dp[:,0])
#print(iMin,iMax)

dashes     = ((),(2,1),(4,0.5,0.5,0.5),(8,1,1,1,1,1,1,1,1,1),(10,1,2,1))
markers    = ('s','o','D','^','v','>','<','d')
linewidths = (2,2)
mec        = ('k','r','b','darkgreen','k')
mfc        = ('w','pink','lightblue','lightgreen','lightgrey')

shift = 0
for i in range(iMin,iMax+1):
#for i in range(iMax,iMin-1,-1):
    shift = 0.2*int(iMax-i)
    ii = (d[:,0]==i)
    iii = np.where(ii==True)[0][0]

    color    = mec[i % len(mec)]
    linewidth= linewidths[i % len(linewidths)]
    dashes_  = dashes[i % len(dashes)]
    marker   = markers[i % len(markers)]
    mec_     = mec[i % len(mec)]
    mfc_     = mfc[i % len(mfc)]

    plot(d[ii,2],d[ii,3]+shift,'-',color=color,linewidth=linewidth, dashes=dashes_, zorder=10);

    ii = (dp[:,0]==i)
    iii = np.where(ii==True)[0][0]
    plot(dp[ii,2],dp[ii,3]+shift,marker=marker, lw=0, markersize=12, mec=mec_, mfc=mfc_, zorder=20);

    rho_str  =  "{:6.3f}".format(  dp[iii,1] )   # 3 decimals for rho
    shift_str =  "{:2.1f}".format( shift )         # 1 decimal for shift
    plot(-dp[ii,2],dp[ii,3],marker=marker, markersize=12, mec=mec_, mfc=mfc_, linewidth=linewidth, color=mec_, dashes=dashes_, label=r'$\rho$[g/cc]='+rho_str+' shift='+shift_str, zorder=20);
    
    #shift += 0.5

#############################################################################

a.set_ylabel(r'(P-P$_0$)/P$_0$')
#xlabel(r'Density (g$\,$cm$^{-3}$)')
#xlabel(r'Compression ratio $\rho / \rho_0$')
a.set_xlabel(r'Temperature (K)')
#ylabel(r'Pressure (GPa)')

a.set_xscale('log')
#a.yaxis.set_minor_locator(MultipleLocator(0.01))

a.xaxis.set_ticks_position('both')
a.get_xaxis().set_tick_params(which='both', direction='in')
a.yaxis.set_ticks_position('both')
a.get_yaxis().set_tick_params(which='both', direction='in')

a.legend(loc='lower right',frameon=False,handlelength=4.5,ncol=1,handletextpad=0.4,numpoints=1)
a.legend(frameon=False,handlelength=4.5,ncol=1,handletextpad=0.4,numpoints=1,bbox_to_anchor=(1.03, 1), loc='upper left', borderaxespad=0.)

#a.set_xlim(1e6,1e9)
#plot([1e6,1e9],[0,0],'--', color='grey',linewidth=1);#print DFT
#a.set_ylim(-0.25,0.01)

savefig(scriptName+'.pdf', bbox_inches='tight')
savefig(scriptName+'.png', bbox_inches='tight', dpi=200)