#For plotting the HRDs

import nugrid_set as set1

path='/apod2/NuGrid/data/set1ext'

symbs=[]
end_model=20*[-1]
symbs=['o','s','D','+','<','*','p','x']
linestyle=['--','-','-.',':','-','-.','--',':']


#t0_model=4849 for M3



from imp import *
ms=load_source('mesa','/rpod2/critter/PPN/forum.astro.keele.ac.uk_update/utils/pylib/mesa.py')

#######M3
historydata=ms.history_data(path+'/set1.5a/see_wind/M3.000Z0.0001/LOGS')

historydata.kippenhahn(xax='logtimerev',t0_model=0,num_frame=1,plot_star_mass=False)
plt.xlim(8.435083281747918,0)
plt.xlim(8.435083281747918,2.696)
plt.ylabel('$M/M_{\odot}$')
plt.ylabel('$M [M_{\odot]$')

plt.xlabel('$log(t_{final}-t)\ [yr]$',fontsize=20,labelpad=1)
plt.legend().set_visible(False)
plt.ylabel('$M$ $[M_{\odot}]$',fontsize=20,)
plt.legend(loc=5,fontsize=12)
plt.legend(loc=6,fontsize=12)
plt.legend(loc=6,fontsize=14)
plt.savefig('M3_set1_5a_kippcont_preAGB.png',dpi=300) #changed for now
'''



################### version 2 of figure
historydata.kippenhahn(xax='time',t0_model=4849,num_frame=1,plot_star_mass=False)
plt.xlim(0,90000)
plt.ylabel('$M-0.836\ [10^{-3}M_{\odot}]$',fontsize=20,labelpad=1)
ax=plt.gca()
plt.tick_params(axis='both', which='major', labelsize=16)
ax.set_xticks([2e4,4e4,6e4,8e4])
plt.ylim(0.836,0.855)
ax.set_xticklabels([2,4,6,8])
ax.set_yticks([0.836,0.84,0.844,0.848,0.852])
ax.set_yticklabels([0,4,8,12,16])
plt.xlabel('$t-t_0\ [10^4yr]$',fontsize=20,labelpad=0)
plt.ylabel('$M-0.836\ [10^{-3}M_{\odot}]$',fontsize=20,labelpad=0)
plt.legend().set_visible(False)
plt.savefig('M3_set1_5a_kippcont_AGB.png')


historydata=ms.history_data(path+'/set1.5a/see_wind/M5.000Z0.0001/LOGS')

mesaset=set1.mesa_set(path+'/set1.5a/see_wind',multi_dir=['M5.000Z0.0001'])
beginTP=mesaset.set_find_first_TP()[0]

historydata.kippenhahn(xax='time',t0_model=beginTP,num_frame=1,plot_star_mass=False)
plt.xlim(0,64000)
plt.title('')
plt.minorticks_on()
ax=plt.gca()
ax.set_yticks([0,4,8,12])
plt.yticks([0.978,0.982,0.986,0.99])
plt.ylim(0.978,0.992)
ax.set_yticklabels([0,4,8,12])
plt.xlabel('$t-t_0\ [10^4yr]$',fontsize=20,labelpad=1)
ax.set_xticklabels([0,1,2,3,4,5,6])
plt.tick_params(axis='both', which='major', labelsize=16)
plt.legend().set_visible(False)
plt.ylabel('$M-0.978\ [10^{-3}M_{\odot}]$',fontsize=20,labelpad=1)
plt.savefig('kipp_TPs_M5_0_Z0_0001.png')

mesaset=set1.mesa_set(path+'/set1.5a/see_wind',multi_dir=['M7.000Z0.0001'])
beginTP=mesaset.set_find_first_TP()[0]
historydata=ms.history_data(path+'/set1.5a/see_wind/M7.000Z0.0001/LOGS')

historydata.kippenhahn(xax='time',t0_model=beginTP,num_frame=1,plot_star_mass=False)
plt.xlim(0,3500)
plt.title('')
ax=plt.gca()
plt.ylabel('$M-1.2714\ [10^{-4}M_{\odot}]$',fontsize=20,labelpad=1)
plt.ylim(1.2714,1.27225)
ax.set_yticks([1.2714,1.2716,1.2718,1.272,1.2722])
ax.set_yticklabels([0,2,4,6,8])
plt.xlabel('$t-t_0\ [10^2yr]$',fontsize=20,labelpad=1)
ax.set_xticklabels([5,10,15,20,25,30,35])
plt.tick_params(axis='both', which='major', labelsize=16)
plt.xlabel('$t-t_0\ [10^2yr]$',fontsize=20,labelpad=1)
plt.legend().set_visible(False)
plt.minorticks_on()
plt.savefig('kipp_TPs_M7_0_Z0_0001.png')

'''