def set_plot_tcrhoc(self,linestyle=['-'],burn_limit=0.997,marker=['o'],markevery=500,end_model=[-1],deg_line=True): 
		'''
			Plots HRDs
			end_model - array, control how far in models a run is plottet, if -1 till end
			symbs_1  - set symbols of runs
		'''
		if len(linestyle)==0:
			linestyle=200*['-']

                for i in range(len(self.runs_H5_surf)):
                        m1p65_last=se(self.runs_H5_out[i])
			t1_model=-1
			if end_model[i] != -1:
				t1_model=end_model[i]

			model=m1p65_last.se.cycles
			model_list=[]
			for k in range(0,len(model),5):
				model_list.append(model[k])
                        rho1=m1p65_last.get(model_list,'rho')   #[:(t1_model-t0_model)]
                        T1=m1p65_last.get(model_list,'temperature')#[:(t1_model-t0_model)]
			T_unit=m1p65_last.get('temperature_unit')
			mass=m1p65_last.get(model_list,'mass')
			mini=m1p65_last.get('mini')
			zini=m1p65_last.get('zini')
			#info=m1p65_last.burnstage_upgrade()
			#burn_info=[]
			#burn_cycle=[]
			#for k in range(len(info[0])):
			#    if 'start' in info[3][k]:
			#	burn_info.append( info[3][k])
			#	burn_cycle.append(info[0][k])
			#print burn_info
			#print burn_cycle
			'''
			#H
			#Maybe use get_elemental_abunds(), too slow if individually taken
			h=m1p65_last.get(model_list,'H-1')
			h_ini=mini*h[0][-1]
			h_limit=h_ini*burn_limit
			print h_limit
			#He
			he=m1p65_last.get(model_list,'He-4')
			he+=m1p65_last.get(model_list,'He-3')
			he_ini=mini*he[0][-1]
			he_limit=he_ini*burn_limit
			print he_limit
			if mini>5:
				c=m1p65_last.get(model_list,'C-12')
				c_ini=mini*c[0][-1]
				c_limit=c_ini*burn_limit
			if mini>10:
				ne=m1p65_last.get(model_list,'Ne-20')
				ne_ini=	mini*ne[0][-1]
				ne_limit=ne_ini*burn_limit
				o=m1p65_last.get(model_list,'O-16')
				o_ini=mini*o[0][-1]
				o_limit=o_ini*burn_limit
				si=m1p65_last.get(model_list,'Si-28')
				si_ini=mini*si[0][-1]
				si_limit=si_ini*burn_limit
			print '#############################3'
			print h_ini,he_ini
			rho=[]
			T=[]
			he_depl_idx=-1
			h_depl_idx=-1
			c_depl_idx=-1
			ne_depl_idx=-1
			o_depl_idx=-1
			si_depl_idx=-1
			'''
			rho=[]
			T=[]
			for k in range(len(model_list)):
				rho_center=rho1[k][0]
				T_center=T1[k][0]
				rho.append(rho_center)
				T.append(T_center)
				'''
				mass11=np.array([0]+list(mass[k]))
				delta_mass=mass11[1:]-mass11[:-1]
				#for the low-mass + massive AGB
				if (sum(np.array(he[k])*np.array(delta_mass))<he_limit) and  (he_depl_idx==-1):
					he_depl_idx=k			
                                if (sum(np.array(h[k])*np.array(delta_mass))<h_limit) and  (h_depl_idx==-1):
                                        h_depl_idx=k
				#for the SAGB + massive
				if mini>5:
					if (sum(np.array(c[k])*np.array(delta_mass))<c_limit) and  (c_depl_idx==-1):
                                        	c_depl_idx=k
				#for the massive stars
				if mini>10:
					if (sum(np.array(ne[k])*np.array(delta_mass))<ne_limit) and  (ne_depl_idx==-1):
                                        	ne_depl_idx=k
					if (sum(np.array(o[k])*np.array(delta_mass))<o_limit) and  (o_depl_idx==-1):
                                        	o_depl_idx=k
					if (sum(np.array(si[k])*np.array(delta_mass))<si_limit) and  (si_depl_idx==-1):
                                        	si_depl_idx=k
				'''

			T=np.log10(np.array(T)*T_unit)
			#T_degeneracy=np.log10(np.array(rho)**(2./3.))
			rho=np.log10(np.array(rho))
                        figure(1)
                        #plot(logTeff,logL,marker=symbs[i],label=self.run_label[i],linestyle=linestyle[i],markevery=markevery)
			#pl.plot(rho,T,marker=symbs[i],label=self.run_label[i],linestyle=linestyle[i],markevery=markevery)
			label=str(mini)+'$M_{\odot}$, Z='+str(zini)
			plt.plot(rho,T,label=label,linestyle=linestyle[i],marker=marker[i],markevery=markevery)
			'''
			#plt.plot(rho,T_degeneracy,color='b',linestyle='-',label='$P_e = P_{e,deg}$')
			print burn_info
			print burn_cycle
			burn_cycle=np.array(burn_cycle)/100.
			print 'changed cycles:', burn_cycle
			print len(rho)
			#For different burning stages:
			if 'H_start' in burn_info:
				plt.plot(rho[burn_cycle[0]],T[burn_cycle[0]],marker='o',color='b')
			if 'He_start' in burn_info:
				plt.plot(rho[burn_cycle[1]],T[burn_cycle[1]],marker='o',color='r')
			if 'C_start' in burn_info:
				plt.plot(rho[burn_cycle[2]],T[burn_cycle[2]],marker='o',color='g')
			if 'Ne_start' in burn_info:
				plt.plot(rho[burn_cycle[3]],T[burn_cycle[3]],marker='D',color='b')
			if 'O_start' in burn_info:
				plt.plot(rho[burn_cycle[4]],T[burn_cycle[4]],marker='D',color='r')
			if 'Si_start' in burn_info:
				plt.plot(rho[burn_cycle[5]],T[burn_cycle[5]],marker='D',color='g')
			'''
		
			#ax = plt.gca()
			#ax.invert_xaxis()
			plt.rcParams.update({'font.size': 16})
			plt.rc('xtick', labelsize=16)
			plt.rc('ytick', labelsize=16)
			legend(loc=4)
		        plt.xlabel('log $\\rho_{\\rm c}$',fontsize=18)
        		plt.ylabel('log $T_{\\rm c}$',fontsize=18)
			'''
			#plt.gca().invert_xaxis()	
			figure(0)
			#pl.plot(rho,T,marker=symbs[i],label=self.run_label[i],linestyle=linestyle[i],markevery=markevery)
                        plt.plot(rho,T,label=self.extra_label[i],linestyle=linestyle[i])
			plt.plot(rho,T_degeneracy,color='b',linestyle='-',label='$P_e = P_{e,deg}$')
			#For different burning stages:
			if not h_depl_idx ==-1:
				plt.plot(rho[h_depl_idx],T[h_depl_idx],marker='o',color='b')
			if not he_depl_idx ==-1:
				plt.plot(rho[he_depl_idx],T[he_depl_idx],marker='o',color='r')
			if not c_depl_idx ==-1:
				plt.plot(rho[c_depl_idx],T[c_depl_idx],marker='o',color='g')
			if not ne_depl_idx ==-1:
				plt.plot(rho[ne_depl_idx],T[ne_depl_idx],marker='D',color='b')
			if not o_depl_idx ==-1:
				plt.plot(rho[o_depl_idx],T[o_depl_idx],marker='D',color='r')
			if not si_depl_idx ==-1:
				plt.plot(rho[si_depl_idx],T[si_depl_idx],marker='D',color='g')
			

			#ax = plt.gca()
			#ax.invert_xaxis()
			plt.rcParams.update({'font.size': 16})
                	plt.rc('xtick', labelsize=16)
                	plt.rc('ytick', labelsize=16)
                	legend(loc=4)
		        plt.xlabel('log $\\rho_{\\rm c}$',fontsize=18)
        		plt.ylabel('log $T_{\\rm c}$',fontsize=18)
			#plt.gca().invert_xaxis()
			'''
			i+=1	
		if deg_line==True:
			rho=np.arange(0,9,0.01)
			T_degeneracy=2./3. *rho +np.log10(1.207e5 * 1.8/(2.**(5./3.)))
			#T_degeneracy=np.log10( (10**np.array(rho))**(2./3.))
			plt.plot(rho,T_degeneracy,color='b',linestyle='-',label='$P_e = P_{e,deg}$')
		plt.legend(loc=2)