read_models_obs.py

from nugridpy import nugridse as mp
import matplotlib.pyplot as plt
import matplotlib as mpl
from numpy import log10
import pandas as pd
from collections import OrderedDict

#run1=mp.se('/run/media/s1785107/Seagate Backup Plus Drive/eddie3/NuclearPhysics/NuPPN/frames/mppnp/RUN_set1upd_m3z3m2/H5_surf','surf.h5')
run1=mp.se('/media/ashley/Seagate Backup Plus Drive/eddie3/NuclearPhysics/NuPPN_Ne22_newer_update4/frames/mppnp/RUN_set1upd_m3z3m2_hCBM_TDUx4p3_ll22/H5_surf','surf.h5')
run2=mp.se('/media/ashley/Seagate Backup Plus Drive/eddie3/NuclearPhysics/NuPPN_Ne22_newer_update5/frames/mppnp/RUN_set1upd_m3z3m2_hCBM_TDUx4p3_ChetecNe22a_stdstd/H5_surf','surf.h5')
run3=mp.se('/media/ashley/Seagate Backup Plus Drive/eddie3/NuclearPhysics/NuPPN_Ne22_newer_update1/frames/mppnp/RUN_set1upd_m3z3m2_hCBM_TDUx4p3_newerNe22a_stdstd/H5_surf','surf.h5')
#run4=mp.se('/run/media/s1785107/Seagate Backup Plus Drive/eddie3/NuclearPhysics/NuPPN_Ne22_newer_update3/frames/mppnp/RUN_set1upd_m3z3m2_hCBM_TDUx4p3_newerNe22a_stdhi/H5_surf','surf.h5')
#run6=mp.se('/run/media/s1785107/Seagate Backup Plus Drive/eddie3/np/NuPPN_ll12/frames/mppnp/RUN_set1upd_m3z3m2_hCBM_TDUx4p3_95Zrll12/H5_surf','surf.h5')
#run7=mp.se('/run/media/s1785107/Seagate Backup Plus Drive/eddie3/np/NuPPN_ll23/frames/mppnp/RUN_set1upd_m3z3m2_hCBM_TDUx4p3_95Zrll23/H5_surf','surf.h5')

#run_mass=[run1,run2,run3,run4,run5,run6,run7,run8,run9,run10,run11,run12,run13]
run_mass=[run1,run2,run3]#,run4]
metallicity_label=[
#'m3z2m2','m3z2m2 RI18']
#'m3z3m2','m3z3m2 Reifarth18']
#'m3z1m2','m3z1m2\_hCBM\_TDUx7','m3z2m2','m3z2m2\_hCBM\_TDUx7']
#'m3z3m2',
#'m3z3m2-hCBM this work', 'm3z3m2-hCBM this work lostd','m3z3m2-hCBM this work stdhi']
'm3z3m2-hCBM longland','m3z3m2-hCBM This work no TAMU','m3z3m2-hCBM This work with TAMU']
#'m2z1m2','m2z2m2','m2z3m2','m3z1m2','m3z2m2','m3z3m2','m2z3p5m2','m3z3p5m2']
#'m3z2m2','m3z2m2-hCBM','RI18']
#'m2z1m2','m2z2m2','m2z3m2','m3z1m2','m3z2m2','m3z3m2','m3z2m2.rotmix.stx5','m3z2m2.rotmix.st']#,'M3.z2m2.rotmix.stx2','M3.z2m2.rotmix.std2'] ## legend labels

sparcity_sindex=2000   ## sparcity to adopt reading s-process index data
sparcity_isoratio=2000 ## sparcity to adopt reading isotopic ratio data
markers=['tab:red','xkcd:indigo','tab:green','tab:brown','k','tab:cyan', 'tab:gray','tab:olive', 'y', 'tab:purple', 'tab:pink', 'tab:orange','tab:cyan']## markers to use while plotting s-process indices
symbols=['o','v','s','o','h','^','>','D','<','p','d','^','s','*'] ## linestyles to use while plotting s-process indices
line=['solid','solid','solid','solid','solid','solid','solid','solid','solid','solid','dashed','dotted']
ls_element = ['Sr-84','Sr-86','Sr-87','Sr-88','Y-89','Zr-90','Zr-91','Zr-92','Zr-94','Zr-96']
hs_element = ['Ba-130','Ba-132','Ba-134','Ba-135','Ba-136','Ba-137','Ba-138','La-138','La-139','Nd-142','Nd-143','Nd-144','Nd-145','Nd-146','Nd-148','Nd-150','Sm-144','Sm-147','Sm-148','Sm-149','Sm-150','Sm-152','Sm-154']
# notice that Pb (3rd s-process index) is not included here.
s_element =  ['Sr-84','Sr-86','Sr-87','Sr-88','Y-89','Zr-90','Zr-91','Zr-92','Zr-94','Zr-96','Ba-130','Ba-132','Ba-134','Ba-135','Ba-136','Ba-137','Ba-138','La-138','La-139','Nd-142','Nd-143','Nd-144','Nd-145','Nd-146','Nd-148','Nd-150','Sm-144','Sm-147','Sm-148','Sm-149','Sm-150','Sm-152','Sm-154']	

t = mpl.markers.MarkerStyle(marker='*')
t._transform = t.get_transform().rotate_deg(180)

# '[Rb/M]', '[Sr/M]', '[Y/M]', '[Zr/M]', '[Ba/M]', '[La/M]',
# '[Nd/M]', '[Sm/M]', '[ls/M]', '[hs/M]', '[hs/ls]', '[s/M]', '[M/H]'
iso_ratio_16='[hs/M]'
iso_ratio_15='[hs/ls]'
iso_ratio_11='[ls/M]'
iso_ratio_14='[s/M]'

# '[M/H]', '[Rb/M]', '[Sr/M]', '[Y/M]', '[Zr/M]', '[Ba/M]',
# '[La/M]', '[Ce/M]', '[Nd/M]', '[Sm/M]'
iso_ratio_3='[Rb/M]'
iso_ratio_4='[Sr/M]'
iso_ratio_5='[Y/M]'
iso_ratio_6='[Zr/M]'
iso_ratio_7='[Ba/M]'
iso_ratio_8='[La/M]'
iso_ratio_9='[Nd/M]'
iso_ratio_10='[Sm/M]'
iso_ratio_12='[Ce/M]'
iso_ratio_13='[M/H]'

#read data from zamora
dfl = pd.read_csv('Zamora-09.txt',skiprows=[1,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,25,26],sep='\t',engine='python') 
dfe = pd.read_csv('Zamora-09.txt',skiprows=[1,2,3,4,5,6,7,8,26],sep='\t',engine='python') 
df5 = pd.read_csv('Zamora-09.txt',skiprows=[1,2,3,4,5,6,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,25],sep='\t',engine='python')
dfl = dfl[~(dfl[iso_ratio_13] <= -0.3)] 
dfe = dfe[~(dfe[iso_ratio_13] <= -0.3)] 

xhsl = dfl[iso_ratio_16].values
xhse = dfe[iso_ratio_16].values
exhsl = df5[iso_ratio_16].values[0]
exhse =  df5[iso_ratio_16].values[1]
xhslsl = dfl[iso_ratio_15].values
xhslse = dfe[iso_ratio_15].values
exhslsl = df5[iso_ratio_15].values[0]
exhslse = df5[iso_ratio_15].values[1]
xlsl = dfl[iso_ratio_11].values
xlse = dfe[iso_ratio_11].values
exlsl = df5[iso_ratio_11].values[0]
exlse =  df5[iso_ratio_11].values[1]
xrbl = dfl[iso_ratio_3].values
xrbe = dfe[iso_ratio_3].values
exrbl = df5[iso_ratio_3].values[0]
exrbe = df5[iso_ratio_3].values[1]
xsl = dfl[iso_ratio_14].values
xse = dfe[iso_ratio_14].values
exsl = df5[iso_ratio_14].values[0]
exse = df5[iso_ratio_14].values[1]
xml = dfl[iso_ratio_13].values
xme = dfe[iso_ratio_13].values
exml = df5[iso_ratio_13].values[0]
exme = df5[iso_ratio_13].values[1]

# read data from abia
df1 = pd.read_csv('Abia-02.txt',skiprows=[30],usecols=[1,2,3,4,5,6,7,8,9,10],sep='\s\s',engine='python')
df2 = df1.fillna(value=0.00)
df2 = df2[~(df2[iso_ratio_13] <= -0.3)] 
df3 = df2.notnull()
df4 = pd.read_csv('Abia-02.txt',skiprows=[1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,25,26,27,28,29,30],sep='\s\s',engine='python')

ls = (df2[iso_ratio_4].values+df2[iso_ratio_5].values+df2[iso_ratio_6].values)\
/(df3[iso_ratio_4].values*1+df3[iso_ratio_5].values*1+df3[iso_ratio_6].values*1) #issues with zero values
els = df4[iso_ratio_11].values[0]
hs = (df2[iso_ratio_7].values+df2[iso_ratio_8].values+df2[iso_ratio_9].values+df2[iso_ratio_10].values)\
/(df3[iso_ratio_7].values*1+df3[iso_ratio_8].values*1+df3[iso_ratio_9].values*1+df3[iso_ratio_10].values*1) #issues with zero values
ehs = df4[iso_ratio_16].values[0]
s = (ls+hs)/2
es = df4[iso_ratio_14].values[0]
hsls = hs-ls
ehsls = df4[iso_ratio_15].values[0]
m = df2[iso_ratio_13].values
em = df4[iso_ratio_13].values[0]
rb = df2[iso_ratio_3].values 
erb = df4[iso_ratio_3].values[0]

df6 = pd.read_csv('Karinkuzhi-18.txt',sep='\t',engine='python') 

ls2 = df6[iso_ratio_11].values[1:]
els2 = df6[iso_ratio_11].values[0]
hs2 = df6[iso_ratio_16].values[1:]
ehs2 = df6[iso_ratio_16].values[0]
s2 = (ls2+hs2)/2
es2 = df6[iso_ratio_15].values[0]
hsls2 = df6[iso_ratio_15].values[1:]
ehsls2 = df6[iso_ratio_15].values[0]
m2 = df6[iso_ratio_13].values[1:]
em2 = df6[iso_ratio_13].values[0]

df7 = pd.read_csv('Kong-18.txt',usecols=[0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18],index_col=[0],sep='\t',engine='python').T
df8 = df7.fillna(value=0.00)
df9 = df8.notnull()

ls3 = (df8[iso_ratio_4].values+df8[iso_ratio_5].values+df8[iso_ratio_6].values)\
/(df9[iso_ratio_4].values*1+df9[iso_ratio_5].values*1+df9[iso_ratio_6].values*1)
els3 = df8['Error[Y/M]'].values
hs3 = (df8[iso_ratio_7].values+df8[iso_ratio_8].values+df8[iso_ratio_9].values+df8[iso_ratio_12].values)\
/(df9[iso_ratio_7].values*1+df9[iso_ratio_8].values*1+df9[iso_ratio_9].values*1+df9[iso_ratio_12].values*1) 
ehs3 = df8['Error[Ba/M]'].values
s3 = (ls3+hs3)/2
es3 = df8['Error[Ba/M]'].values
hsls3 = hs3-ls3
ehsls3 = df8['Error[Ba/M]'].values
m3 = df8[iso_ratio_13].values
em3 = df8['Error[M/H]'].values

df10 = pd.read_csv('elementi_m3p0z2m2_000_20190111_135823.txt',skiprows=[38,56],sep='\s+',engine='python').T
df10.columns = df10.iloc[0]
df10 = df10[1:]

df11 = pd.read_csv('sproc_m3p0z2m2_000_20190114_163131.txt',sep='\s+',engine='python')

s_1 = (df11['ls/Fe'].values + df11['hs/Fe'].values)/2
rb_1 = df10['Rb'].values[2:]
rb_3 = rb_1-0.068
c = df10['C'].values[2:]
o = df10['O'].values[2:]
co = c - o 
co_ini = (2.26e-03+2.75e-05)/(6.07e-03+2.40e-06)
c_o = (10**co)*co_ini*(4./3.)
rb_2=[]
s_2=[]
for i in range(len(s_1)):
    if c_o[i] >= 1.:
        rb_2.append(rb_3[i])
        s_2.append(s_1[i])

s_ini = []
ls_ini = []
hs_ini = []
fe_ini = []
rb_ini = []
for i in run_mass:
    dum_s_ini = 0.
    dum_ls_ini = 0.
    dum_hs_ini = 0.
    dum_fe_ini = 0.
    dum_rb_ini = 0.
    dum_fe_ini = float(i.se.get(min(i.se.cycles),'iso_massf','Fe-54'))+float(i.se.get(min(i.se.cycles),'iso_massf','Fe-56'))+float(i.se.get(min(i.se.cycles),'iso_massf','Fe-57'))+float(i.se.get(min(i.se.cycles),'iso_massf','Fe-58'))	
    dum_rb_ini = float(i.se.get(min(i.se.cycles),'iso_massf','Rb-85'))+float(i.se.get(min(i.se.cycles),'iso_massf','Rb-87'))
    for j in ls_element:
        dum_ls_ini = dum_ls_ini + float(i.se.get(min(i.se.cycles),'iso_massf',j))
    for j in hs_element:
        dum_hs_ini = dum_hs_ini + float(i.se.get(min(i.se.cycles),'iso_massf',j))
    for j in s_element:
        dum_s_ini = dum_s_ini + float(i.se.get(min(i.se.cycles),'iso_massf',j))
    fe_ini.append(dum_fe_ini)
    ls_ini.append(dum_ls_ini)    #  /float(len(ls_element)))
    hs_ini.append(dum_hs_ini)     #/float(len(hs_element)))
    s_ini.append(dum_s_ini)  #/float(len(s_element)))				
    rb_ini.append(dum_rb_ini)
# sparcity for cycles I am looking at.
sparsity = sparcity_sindex
s_fe = []
ls_fe = []
hs_fe = []
hs_ls = []
fe_h = []
rb_fe = []
s_fe_tps = []
ls_fe_tps = []
hs_fe_tps = []
hs_ls_tps = []
fe_h_tps = []
rb_fe_tps = []
s_fe_tps_co  = []
ls_fe_tps_co = []
hs_fe_tps_co = []
hs_ls_tps_co = []
fe_h_tps_co = []
rb_fe_tps_co = []
k = 0

for i in run_mass:
    jjjj=0
    dum_s_fe = []
    dum_ls_fe = []
    dum_hs_fe = []
    dum_hs_ls = []
    dum_fe_h = []
    dum_rb_fe = []
    dum_s_fe_tps = []
    dum_ls_fe_tps = []
    dum_hs_fe_tps = []
    dum_hs_ls_tps = []
    dum_fe_h_tps = []
    dum_rb_fe_tps = []
    dum_s_fe_tps_co = []
    dum_ls_fe_tps_co = []
    dum_hs_fe_tps_co = []
    dum_hs_ls_tps_co = []
    dum_fe_h_tps_co = []
    dum_rb_fe_tps_co = []
    dum_co=[]

    for j in i.se.cycles[0::sparsity]:
        dum_s = 0.
        dum_ls = 0.
        dum_hs = 0.
        dum_fe = 0.
        dum_rb = 0.
        dum_c = 0.
        dum_o = 0.
        print(j)
        dum_fe = float(i.se.get(j,'iso_massf','Fe-54'))+float(i.se.get(j,'iso_massf','Fe-56'))+float(i.se.get(j,'iso_massf','Fe-57'))+float(i.se.get(j,'iso_massf','Fe-58')) 
        dum_rb = float(i.se.get(j,'iso_massf','Rb-85'))+float(i.se.get(j,'iso_massf','Rb-87'))
        dum_c = float(i.se.get(j,'iso_massf','C-12'))
        dum_o = float(i.se.get(j,'iso_massf','O-16'))
        dum_c=(float((i.se.get((int(j)),'iso_massf','C-12')))+float((i.se.get(j,'iso_massf','C-13'))))
        dum_o=(float((i.se.get((int(j)),'iso_massf','O-16')))+float((i.se.get(j,'iso_massf','O-17'))))
        dum_co.append(((dum_c/dum_o)*(16./12.)))
        for jj in ls_element:
            dum_ls = dum_ls + float(i.se.get(j,'iso_massf',jj)) #/float(len(ls_element)))
        for jj in hs_element:
            dum_hs = dum_hs + float(i.se.get(j,'iso_massf',jj))  #/float(len(hs_element)))
        for jj in s_element:
            dum_s = dum_s + float(i.se.get(j,'iso_massf',jj)) #/float(len(s_element)))
        dum_s_fe.append(log10((dum_s/dum_fe)/(s_ini[k]/fe_ini[k])))
        dum_ls_fe.append(log10((dum_ls/dum_fe)/(ls_ini[k]/fe_ini[k])))
        dum_hs_fe.append(log10((dum_hs/dum_fe)/(hs_ini[k]/fe_ini[k])))
        dum_hs_ls.append(log10((dum_hs/dum_ls)/(hs_ini[k]/ls_ini[k])))
        dum_fe_h.append(log10((float(i.se.get('zini')))/(0.018)))
        #dum_rb_fe.append(log10((dum_rb/dum_fe)/(rb_ini[k]/fe_ini[k]))-0.2) ## correction for Lambert 1995 data
        dum_rb_fe.append(log10((dum_rb/dum_fe)/(rb_ini[k]/fe_ini[k]))-0.068) ## correction for Zamora 2009 data

        if (len(dum_co)>1):
            if (dum_co[len(dum_co)-1]>(dum_co[len(dum_co)-2]+0.02)):
                dum_s_fe_tps.append(log10((dum_s/dum_fe)/(s_ini[k]/fe_ini[k])))
                dum_ls_fe_tps.append(log10((dum_ls/dum_fe)/(ls_ini[k]/fe_ini[k])))
                dum_hs_fe_tps.append(log10((dum_hs/dum_fe)/(hs_ini[k]/fe_ini[k])))
                dum_hs_ls_tps.append(log10((dum_hs/dum_ls)/(hs_ini[k]/ls_ini[k])))
                dum_fe_h_tps.append(log10((float(i.se.get('zini')))/(0.018)))
                #dum_rb_fe_tps.append(log10((dum_rb/dum_fe)/(rb_ini[k]/fe_ini[k]))-0.2) ## correction for Lambert 1995 data
                dum_rb_fe_tps.append(log10((dum_rb/dum_fe)/(rb_ini[k]/fe_ini[k]))-0.068) ## correction for Zamora 2009 data

        if (dum_co[len(dum_co)-1]>1.):
            dum_s_fe_tps_co.append(log10((dum_s/dum_fe)/(s_ini[k]/fe_ini[k])))
            dum_ls_fe_tps_co.append(log10((dum_ls/dum_fe)/(ls_ini[k]/fe_ini[k])))
            dum_hs_fe_tps_co.append(log10((dum_hs/dum_fe)/(hs_ini[k]/fe_ini[k])))
            dum_hs_ls_tps_co.append(log10((dum_hs/dum_ls)/(hs_ini[k]/ls_ini[k])))
            dum_fe_h_tps_co.append(log10((float(i.se.get('zini')))/(0.018)))
            dum_rb_fe_tps_co.append(log10((dum_rb/dum_fe)/(rb_ini[k]/fe_ini[k]))-0.068) ## correction for Zamora 2009 data
        jjjj=jjjj+1
	
    s_fe.append(dum_s_fe)
    ls_fe.append(dum_ls_fe)
    hs_fe.append(dum_hs_fe)
    hs_ls.append(dum_hs_ls)
    fe_h.append(dum_fe_h)
    rb_fe.append(dum_rb_fe)
    s_fe_tps.append(dum_s_fe_tps)
    ls_fe_tps.append(dum_ls_fe_tps)
    hs_fe_tps.append(dum_hs_fe_tps)
    hs_ls_tps.append(dum_hs_ls_tps)
    fe_h_tps.append(dum_fe_h_tps)
    rb_fe_tps.append(dum_rb_fe_tps)
    s_fe_tps_co.append(dum_s_fe_tps_co)
    ls_fe_tps_co.append(dum_ls_fe_tps_co)
    hs_fe_tps_co.append(dum_hs_fe_tps_co)
    hs_ls_tps_co.append(dum_hs_ls_tps_co)
    fe_h_tps_co.append(dum_fe_h_tps_co)
    rb_fe_tps_co.append(dum_rb_fe_tps_co)
    k = k+1

mass_label =[]
for i in run_mass:
    mass_label.append(float(i.se.get('mini')))
params = {'text.usetex': True,
        'xtick.direction': 'in',
        'ytick.direction': 'in',
        #'axes.linewidth' : 5,
        'xtick.major.size': 20,
        'ytick.major.size': 20,
        'xtick.labelsize': 40,
        'ytick.labelsize': 40,
        'ytick.major.pad': 5,
        'ytick.major.width': 2,#3,
        'xtick.major.pad': 5,
        'xtick.major.width': 2}#3}
plt.rcParams.update(params)
plt.tick_params(axis='both', pad=5)
# Axes object: one row, one column, first plot (one plot!)
fig = plt.figure(1)            # Figure object
ax = fig.add_subplot(1,1,1)
array_to_plot_x = hs_ls
array_to_plot_x_tps = hs_ls_tps
array_to_plot_x_co = hs_ls_tps_co
array_to_plot_y = hs_fe
array_to_plot_y_tps = hs_fe_tps
array_to_plot_y_co = hs_fe_tps_co

plt.axhline(y=0,linewidth=2, color='k')
plt.axvline(x=0,linewidth=2, color='k')
for k in range(0,len(array_to_plot_x)):
    plt.plot(array_to_plot_x_tps[k],array_to_plot_y_tps[k],c=markers[k],ls=line[k],markersize=12.,linewidth=2.)
    plt.plot(array_to_plot_x_co[k],array_to_plot_y_co[k],marker=symbols[k],c=markers[k],ls=line[k],markersize=20.,linewidth=2.,label=metallicity_label[k])
plt.errorbar(hsls,hs,xerr=ehsls,yerr=ehs,marker=symbols[-1],c=markers[-1],ls='None',markersize=8.,label='Abia et al. 2002')#'$Abia$ $et$ $al.$ $2002$')
plt.errorbar(xhslse,xhse,xerr=exhslse,yerr=exhse,marker=symbols[-2],c=markers[-2],ls='None',markersize=8.,label='Zamora et al. 2009')#'$Zamora$ $et$ $al.$ $2009,$ $Early$ $Stars$')
plt.errorbar(xhslsl,xhsl,xerr=exhslsl,yerr=exhsl,marker=symbols[-2],c=markers[-2],ls='None',markersize=8.)#,label='$Zamora$ $et$ $al.$ $2009,$ $Late$ $Stars$')
#plt.errorbar(hsls2,hs2,xerr=ehsls2,yerr=ehs2,marker=symbols[-4],c=markers[-4],ls='None',markersize=8.,label='$Karinkuzhi$ $et$ $al.$ $2018$')
#plt.errorbar(hsls3,hs3,xerr=ehsls3,yerr=ehs3,marker=symbols[-5],c=markers[-5],ls='None',markersize=8.,label='$Kong$ $et$ $al.$ $2018$')
handles, labels = plt.gca().get_legend_handles_labels()
by_label = OrderedDict(zip(labels, handles))
plt.legend(by_label.values(), by_label.keys(),prop={'size':25})
plt.xlabel('[hs/ls]', fontsize=40)
plt.ylabel('[hs/Fe]', fontsize=40)
plt.xlim(-1,0.7)
plt.ylim(-0.6,2)
plt.show()

fig = plt.figure(2)            # Figure object
ax = fig.add_subplot(1,1,1)
array_to_plot_y = hs_ls
array_to_plot_y_tps = hs_ls_tps
array_to_plot_y_co = hs_ls_tps_co
array_to_plot_x = fe_h
array_to_plot_x_tps = fe_h_tps
array_to_plot_x_co = fe_h_tps_co

plt.axhline(y=0,linewidth=2, color='k')
plt.axvline(x=0,linewidth=2, color='k')
for k in range(0,len(array_to_plot_x)):
    plt.plot(array_to_plot_x_tps[k],array_to_plot_y_tps[k],c=markers[k],ls=line[k],markersize=12.,linewidth=2.)
    plt.plot(array_to_plot_x_co[k],array_to_plot_y_co[k],marker=symbols[k],c=markers[k],ls=line[k],markersize=20.,linewidth=2.,label=metallicity_label[k])
plt.errorbar(m,hsls,xerr=em,yerr=ehsls,marker=symbols[-1],c=markers[-1],ls='None',markersize=8.,label='Abia et al. 2002')#'$Abia$ $et$ $al.$ $2002$')
plt.errorbar(xme,xhslse,xerr=exme,yerr=exhslse,marker=symbols[-2],c=markers[-2],ls='None',markersize=8.)#,label='$Zamora$ $et$ $al.$ $2009,$ $Early$ $Stars$')
plt.errorbar(xml,xhslsl,xerr=exml,yerr=exhslsl,marker=symbols[-2],c=markers[-2],ls='None',markersize=8.,label='Zamora et al. 2009')#'$Zamora$ $et$ $al.$ $2009,$ $Late$ $Stars$')
#plt.errorbar(m2,hsls2,xerr=em2,yerr=ehsls2,marker=symbols[-4],c=markers[-4],ls='None',markersize=8.,label='$Karinkuzhi$ $et$ $al.$ $2018$')
#plt.errorbar(m3,hsls3,xerr=em3,yerr=ehsls3,marker=symbols[-5],c=markers[-5],ls='None',markersize=8.,label='$Kong$ $et$ $al.$ $2018$')
handles, labels = plt.gca().get_legend_handles_labels()
by_label = OrderedDict(zip(labels, handles))
plt.legend(by_label.values(), by_label.keys(),prop={'size':25})
plt.ylabel('[hs/ls]', fontsize=40)
plt.xlabel('[Fe/H]', fontsize=40)
plt.xlim(-0.55,0.65)
plt.ylim(-1,0.75)
plt.show()
       
fig = plt.figure(3)            # Figure object
ax = fig.add_subplot(1,1,1)
array_to_plot_x = hs_ls
array_to_plot_x_tps = hs_ls_tps
array_to_plot_x_co = hs_ls_tps_co
array_to_plot_y = ls_fe
array_to_plot_y_tps = ls_fe_tps
array_to_plot_y_co = ls_fe_tps_co

plt.axhline(y=0,linewidth=2, color='k')
plt.axvline(x=0,linewidth=2, color='k')
for k in range(0,len(array_to_plot_x)):
    plt.plot(array_to_plot_x_tps[k],array_to_plot_y_tps[k],c=markers[k],ls=line[k],markersize=12.,linewidth=2.)
    plt.plot(array_to_plot_x_co[k],array_to_plot_y_co[k],marker=symbols[k],c=markers[k],ls=line[k],markersize=20.,linewidth=2.,label=metallicity_label[k])
plt.errorbar(hsls,ls,xerr=ehsls,yerr=els,marker=symbols[-1],c=markers[-1],ls='None',markersize=8.,label='Abia et al. 2002')#'$Abia$ $et$ $al.$ $2002$')
plt.errorbar(xhslse,xlse,xerr=exhslse,yerr=exlse,marker=symbols[-2],c=markers[-2],ls='None',markersize=8.,label='Zamora et al. 2009')##'$Zamora$ $et$ $al.$ $2009,$ $Early$ $Stars$')
plt.errorbar(xhslsl,xlsl,xerr=exhslsl,yerr=exlsl,marker=symbols[-2],c=markers[-2],ls='None',markersize=8.)#,label='$Zamora$ $et$ $al.$ $2009,$ $Late$ $Stars$')
#plt.errorbar(hsls2,ls2,xerr=ehsls2,yerr=els2,marker=symbols[-4],c=markers[-4],ls='None',markersize=8.,label='$Karinkuzhi$ $et$ $al.$ $2018$')
#plt.errorbar(hsls3,ls3,xerr=ehsls3,yerr=els3,marker=symbols[-5],c=markers[-5],ls='None',markersize=8.,label='$Kong$ $et$ $al.$ $2018$')
handles, labels = plt.gca().get_legend_handles_labels()
by_label = OrderedDict(zip(labels, handles))
plt.legend(by_label.values(), by_label.keys(),prop={'size':25})
plt.xlabel('[hs/ls]', fontsize=40)
plt.ylabel('[ls/Fe]', fontsize=40)
plt.xlim(-1,0.7)
plt.ylim(-0.6,1.6)
plt.show()

fig = plt.figure(4)            # Figure object
ax = fig.add_subplot(1,1,1)
array_to_plot_x = s_fe
array_to_plot_x_tps = s_fe_tps
array_to_plot_x_co = s_fe_tps_co
array_to_plot_y = rb_fe
array_to_plot_y_tps = rb_fe_tps
array_to_plot_y_co = rb_fe_tps_co

plt.axhline(y=0,linewidth=2, color='k')
plt.axvline(x=0,linewidth=2, color='k')
for k in range(0,len(array_to_plot_x)):
    plt.plot(array_to_plot_x_tps[k],array_to_plot_y_tps[k],c=markers[k],ls=line[k],markersize=8.,linewidth=2.,zorder=10)
    plt.plot(array_to_plot_x_co[k],array_to_plot_y_co[k],marker=symbols[k],c=markers[k],ls=line[k],markersize=12.,linewidth=2.,label=metallicity_label[k],zorder=10)
plt.errorbar(s,rb,xerr=es,yerr=erb,marker=t,c=markers[-1],ls='None',markersize=10.,label='Abia et al. 2002')#'$Abia$ $et$ $al.$ $2002$') 
plt.errorbar(xse,xrbe,xerr=exsl,yerr=exrbe,marker=symbols[-2],c=markers[-2],ls='None',markersize=10.,fmt=symbols[-1], mfc='white', zorder=1,label='Zamora et al. 2009')#'$Zamora$ $et$ $al.$ $2009,$ $Early$ $Stars$')
plt.errorbar(xsl,xrbl,xerr=exsl,yerr=exrbl,marker=symbols[-2],c=markers[-2],ls='None',markersize=10.,fmt=symbols[-1], mfc='white', zorder=1)#'$Zamora$ $et$ $al.$ $2009,$ $Late$ $Stars$')
##plt.plot(s_1,rb_3,c=markers[-7],ls='-')
##plt.plot(s_2,rb_2,marker=symbols[-5],c=markers[-7],ls='-',markersize=20.,label='FRUITY')
#for i in s2:
#    plt.axvline(x=i,linewidth=1, color=marker[-4])
#for j in s3:
#    plt.axvline(x=j,linewidth=1, color='y')
#plt.errorbar(s2,[0.25,0.25,0.25],xerr=es2, c=marker[-4],marker=symbols[-4],markersize=8.,label='$Karinkuzhi$ $et$ $al.$ $2018$') 
#plt.errorbar(s3,[0.25,0.25,0.25,0.25,0.25,0.25,0.25,0.25,0.25,0.25,0.25,0.25,0.25,0.25,0.25,0.25,0.25,0.25],xerr=es3, color='k',fmt='y>',markersize=8.,label='$Kong$ $et$ $al.$ $2018$')     
handles, labels = plt.gca().get_legend_handles_labels()
by_label = OrderedDict(zip(labels, handles))
plt.legend(by_label.values(), by_label.keys(),prop={'size':30})
plt.xlabel('[s/Fe]', fontsize=40)
plt.ylabel('[Rb/Fe]', fontsize=40)
plt.xlim(-0.6,1.7)
plt.ylim(-0.6,1)
plt.subplots_adjust(left= 0.115, bottom=0.125, right=0.99, top=0.98)
plt.show()