From 6cc54cc88ad08700581248725ff402703a7ee254 Mon Sep 17 00:00:00 2001 From: kpolsen Date: Mon, 11 Sep 2017 15:06:30 -0700 Subject: [PATCH] cleaned up analysis.py and plot.py in terms of comments. --- Tables/Observations/SFR_Mstar/Jiang16.txt | 28 ++++++ Tables/observations_z6 | Bin 0 -> 8922 bytes sigame/analysis.py | 57 ++++++----- sigame/plot.py | 114 +++++++++------------- 4 files changed, 102 insertions(+), 97 deletions(-) create mode 100644 Tables/Observations/SFR_Mstar/Jiang16.txt create mode 100644 Tables/observations_z6 diff --git a/Tables/Observations/SFR_Mstar/Jiang16.txt b/Tables/Observations/SFR_Mstar/Jiang16.txt new file mode 100644 index 0000000..0380b42 --- /dev/null +++ b/Tables/Observations/SFR_Mstar/Jiang16.txt @@ -0,0 +1,28 @@ +# ID age(M) mass(G) SFR_Lya SFR_UV E(B-V) + 03 4 0.717 9.5 20.9 0.04 + 04 4 0.336 8.0 9.9 0.00 + 15 8 2.113 19.0 33.9 0.08 + 20 492 3.282 15.1 10.5 0.02 + 23 444 5.794 4.6 16.0 0.00 + 24 4 0.377 12.4 15.2 0.00 + 25 956 6.653 13.3 15.9 0.02 + 27 324 1.735 6.3 6.9 0.00 + 28 4 0.268 16.6 10.4 0.00 + 30 912 2.684 4.2 7.8 0.00 + 31 4 0.238 17.1 8.4 0.04 + 34 992 25.000 3.7 35.0 0.12 + 35 4 8.106 16.5 19.5 0.34 + 36 988 4.639 2.5 16.1 0.00 + 43 972 3.076 4.7 8.6 0.02 + 44 872 17.291 18.1 32.6 0.08 + 47 184 39.111 4.9 48.4 0.18 + 49 212 3.556 5.4 9.7 0.08 + 50 4 0.897 8.7 7.7 0.22 + 54 460 5.154 11.8 14.1 0.00 + 58 960 6.162 4.9 17.5 0.04 + 61 32 6.634 11.4 11.0 0.20 + 62 472 11.501 13.7 23.9 0.04 + 63 8 0.945 16.5 22.4 0.02 + 64 200 4.132 11.8 20.2 0.00 + 66 4 0.135 7.6 7.1 0.00 + 67 4 3.201 20.8 43.3 0.14 diff --git a/Tables/observations_z6 b/Tables/observations_z6 new file mode 100644 index 0000000000000000000000000000000000000000..757a9c95c688d1f3cdb93a1381ca6f01abb5c262 GIT binary patch literal 8922 zcmb7K30xD$_YbEK6j1>=1o4hwC?OC6R@evRa?8o1%io&gjo`+D2P??Xw`cB zEGl0AdSBIQv1+R*+SdE1^(fx8`m?o)M_aZ2@IO1*#LWO&m(PbcGy8pS-n{qb&9SL) z&M#(+1x!g~zKPXF7P3r{R-DY3nG{q4PJJf<=Q7|Tp*(Yov}UFNz5=c^aHCvv&EtxJ zJLSnWAwJNCa?3>k4+o)oqgl%u8AFL!Z7}7JPG^kFNG&S{o=~<}LIW=e1~1T#)s{r& zF(q1Yu2K78ZkiFch(X)dVLYN;8Ss`+u1JUW`9@1o@i^SZ6 zLV2Y#dXqtKG;zX1DECY=!^sGxJTfi$I@pOEbfwy5X*nP0^YdYY0)$asIr%z9%W9$2 zSkR4XlT)NObMm@V;-Pwj!DNQDSkOaj&ZAt7FeE_F3`Nk(5(5<&u+<4c?$!A&CG{h)RKs@`QDsGn#UPm68G=r`!$J{$eOMPS+9;OS!XJ zro?0f3d-GLeHlkNmBQiZtOHSMX%O!)gaqVmO@}Hpi2f3Xmv3?}oho(E`-L{|4?y03 z)6x4P=zYpV1RrIMfrQ-&t(3@18Nfp-rX=DJRWnEigwTM{h?vxgk2%32k-gv`ofNS4 z{iEHUsD)2pEIL5c0*F;qHv)qs2_SOxbjBgc0*F;gBkYI2 z+O?-aD%Cc7*dUoaCN4o14dZekL$w*h7+^H(L7E-sAajsFX2^R5{janS#`V*15OJjr znmJ25C_QDv{Y_*+uY;?*t^A5SRhg#kcx^Yy?s5LrWz11Bc+`TjZxW1TLPMYEk!^n> zSAX}1>-rIOq|bqw+v8VmChH#$xlWyCU8>c%q?Ue;4(h2 zsg4}5{c!TqGwR3g#-MIj;Thy!Hxy}otZ!oMb`)m(o1bnG&6ZLgZ0 zMe15|+E!QP4exL8uF?fvb=tr0gTe93$dQYe44VPgLoH1uf0wT(oBDk*ZSt~OdsnWF z3!CP@wVSLNzj;mXlN-ppY}f6wBlAfKdn}`&^Oxj-F`udHj&4U0mLY*z&HCqLFNw4c z#?Rz*@I;m!xc|)HnAAyslj-KSyq6B?qwbn_=lHRL&UR{x46kiX>yJ&5azV|0-tM4i7WEg+ZUyXJ3wxFx=zME68l|EKoE9DB;Rb4B8fBWmuUsd9B zpBE1h-%@dJeK()u-;0U&XYc=0^8P&)dVTJpst`SQRW>JMU1jrK!Th1O<&T_J^G8+y zFD8(u9eFC}sS8hKJoV$L4VT2%*YOn7<{e&}95_MTd59ZNQSO3e1z)$3etf}(v(Zq# zJ&dPVSCpVDo7T=e)K<4?%{PnM@&y~tMt%AAV4iY0XUM;Z+MxnJpW1S5-<)7E(Qx+% zAI(JtaocNmX<0@%(dXcjDO)@Bg~q+uaEA%H+Hv~4RJ^{1WWKob=CAm;s{`2z0m1KO)=K*yqe?(EPJ$4LD(sw14%fRg^vvhGt1D?(b1ziQdajB>KK7ZofD&^UOoBRp1c!m2D^8#FlcuW0>HjVxRv)=S0 zh)zuXj}d;vw0!zvii;nS_}>c~pVs&i$sZeT{nOW<_(c&^H+;E2z8;L^x^hKXD-8PW z>iA0yoKbz8$m<-fzO%ZKqZNf4>TxY~j;wG|--!IN)ng}Ok$kon)#o-X#aZ>Q#Y%%w z_DI&+S3Y~e!AqkBlAo=+e>y!cOVt|1;?p2tg+msJbR!5vi>F#?a+Q^~&Bm=*gtHc> z@~lIxRF`F?Ba}F+@wHfKV0vK+hk+Lh-?;MoWRjkj^1o@8jpUK@CwizoH;{Fw zXKcQB+XB5|Hy#aug!0Wv&dCWI6s3%bj(~d&Wmvwcs2H%6Z&J#@x5EZS%c7#BvN&0+ zwY`M$$;!!u&2pJEI?CE+=G){^N@=vr+Gf$}TG^+J5$I|dTe_-WQUGzyYu>t#-VWdw zk04@L{PrKlGzH-1`jBIFM)B$+g3pu2T`0&_jy-_~*;pP5GdY5CaNC2qe%wo;G$`fX zJH0CQ#$&pPrqvUj3cXg$t*^$1F+l)bR{U%k^#(=g=k?gM?_jgoDic7hE9oc#(ee** zspUMB8*z7jsO@>mEeI1O8fT5O=$Ax?IGH4&+T|pt zguyLobc`%EG7fwsXoylq%E4rVp^L=vjw8tI!8?Pzt9Em z{Zj#ii-d|Pvv6z9;UO=v=G+z`G?@1q)|`)NA`%7S%KoRtZX5z>L5mO?eAbFJt10(^ zN$Ft<5-em|3(iFX%?6-NZL^?U%t2gGP@?_k*mb!m?NnrTM08p>z|C5Nv20nLilF>_ zL<9@kCI5>ut*^yOgD-8HB`3#4ZyP5Ey0VH}Wx}c-#>Y6{_a%Nv zs{JF`SGnKv99m-EJQ8lH~qPix3)o+hTHe2+chL1#AI>xoochVqRST!DY1! zpH_7O2zOWf{CgD+*(ZRYaH;DKAvD;padx~`t4BbM6Y{J_-PXr{xts}yP=N;zP;P0N z8mSB%WOS|MhVKO3S#L3#U8H?4fN&Ze-8&t({SQ3kCAL4*B7_EqUv2v%m;xeEv-`@U zIE2W=4SNiSh#ukd^>~XA8l1qdG*2R+&GvnIzI8~q?382fJ_VO!eUW~&37l@#_&*Bz zLsw2rv*q{<2jQ1vs;%42*`9OAxs#G7KlWKgPT2D-kgH>H_>+J|^K+X=t1D)ZJI&d% zuhy&}{THpOT?OkH%h`X$LW6UBmU&gWJui^==OB2XU*Jq^(q9HMGF*uA5aq~Zy$kN1 z3j&ro3Hx#eolYQT4?D6bVtapLY1!~z9Zp8$EI$ia(7GL%NT9pwZHd-%L!km^X%MiW zqs2|zB~CoQfWHA*XmHWGaAPNkO9JlD>GNaX$#|kVej~4ZL;Pdafp43tr776W_GJm> zm6Vy2o|z+)E8^n870NT?op+K`6!Mr@aFz0)l5=E==y(~p)@stcF6dKlYOFjozYB5e zA9&tB7DSB64D}y1vsuji!m$YU7n5Aaqs6+L@+;i&H+U$VxUV`n{`yiV4Q_JpUwLx8 zh3N+K3Be^5cp-z&#{UYSGcQ#AcJEd<;@RcNQES$axLa=vph#FnKLXJ#3WxsI4C+sy zRIHLB&;Q=S$?*=y15Jujt)Z>m<1rb`P;Qy`Zh<^BP`_>r=>C^_1i-KD%r zmD1FVWNB=4yaL?Qww8t)1zK6Zw*>P`W3#QIX4wS zIf`X~hcAWF;1QbnubTIN3Z&q7QCNq6380&9MTDGeOdxj5Iu!Kg^=21Ge+!^!=NOn| zoviOB;of~LfLiS+=Ticg3TW`81@Ar;@YMhHV6r|Vj4e=jS6I;{SQtx*InR<8`Tv;DMg#x= literal 0 HcmV?d00001 diff --git a/sigame/analysis.py b/sigame/analysis.py index ca37d82..516d6d7 100644 --- a/sigame/analysis.py +++ b/sigame/analysis.py @@ -27,8 +27,19 @@ """ Global line luminosities """ #------------------------------------------------------------------------------- -def ISM_line_contributions(split=True,line='CII'): - print('Split L_[CII] into contributions from different ISM phases!') +def ISM_line_contributions(line='CII'): + ''' + Purpose + --------- + Splits line luminosity into contributions from different ISM phases + + Arguments + --------- + line: which line to plot - str + default = 'CII' + ''' + + print('') plt.close('all') # close all windows # Load model results @@ -69,32 +80,18 @@ def ISM_line_contributions(split=True,line='CII'): plt.savefig('plots/line_emission/ISM_contributions/['+line+']_lum_fractions_'+z1+ext_DENSE+ext_DIFFUSE+'.pdf', format='pdf') # .eps for paper! - print('Now plotting L_CII fractions vs mass fractions') - - # Starting figure - fontsize = 13 - mpl.rcParams['xtick.labelsize'] = fontsize - mpl.rcParams['ytick.labelsize'] = fontsize - - xr = [0.1,100] - yr = [0.1,100] - # Set up figure - fig = plt.figure(1,figsize = (8,7)) - ax1 = fig.add_subplot(1,1,1) - # xlog='y',ylog='y',\ - # colored by ISM phase version - plot.simple_plot(add='y',xr=xr, yr=yr,\ - xlab='Mass fraction [%]',ylab='Fraction of total L$_{[\mathrm{CII}]}$ [%]',\ - x1=f_M[0],y1=f_L[0],ma1='x',fill1='y',ms1=6,mew1=2,col1='b',\ - x2=f_M[1],y2=f_L[1],ma2='x',fill2='y',ms2=6,mew2=2,col2='orange',\ - x3=f_M[2],y3=f_L[2],ma3='x',fill3='y',ms3=6,mew3=2,col3='r') - +def ISM_line_efficiency(line='CII'): + ''' + Purpose + --------- + Plots the line efficiency (line luminosity / gas mass) in different ISM phases - plt.show(block=False) - plt.savefig('plots/line_emission/ISM_contributions/['+line+']_lum_mass_fractions_'+z1+ext_DENSE+ext_DIFFUSE+'.pdf', format='pdf') # .eps for paper! + Arguments + --------- + line: which line to plot - str + default = 'CII' + ''' -def ISM_line_efficiency(line='CII'): - print('Split L_[CII] into contributions from different ISM phases!') plt.close('all') # close all windows # Load model results @@ -164,7 +161,13 @@ def ISM_line_efficiency(line='CII'): #------------------------------------------------------------------------------- def ISM_mass_contributions(): - print('Investigate mass from different ISM phases!') + ''' + Purpose + --------- + Plots the mass contribution from different ISM phases + + ''' + plt.close('all') # close all windows # Load model results diff --git a/sigame/plot.py b/sigame/plot.py index 667c67a..f5ce7d6 100644 --- a/sigame/plot.py +++ b/sigame/plot.py @@ -859,7 +859,6 @@ def SFR_Mstar(): Plot of SFR vs stellar mass ''' - print('** Plot SFR vs stellar mass for sample **') plt.close('all') models = aux.find_model_results() @@ -883,25 +882,18 @@ def SFR_Mstar(): age = 0.984 # age of universe, Omega_m = 0.3, Omega_lambda = 0.7, h = 0.65 MS = 10.**((0.84-0.026*age)*np.log10(xr)-(6.51-0.11*age)) # z ~ 6 LBGs and LAEs - # columns = ['No','M_1500','UV_slope','L_Lyalpha','ew_Lyalpha','SFR_UV','SFR_Lyalpha'] - # J16 = pd.read_table('/Users/Karen/Desktop/temp/Observations/Jiang16_paper',names=columns,skiprows=5,sep=r'\t*',engine='python') - # columns = ['ID','age(M)','mass(G)','SFR_Lya','SFR_UV','E(B-V)'] - # J16 = pd.read_table('sigame/temp/Observations/Jiang16.txt',names=columns,skiprows=1,sep=r'\s*',engine='python') - # SFR_J16 = J16['SFR_UV'].values - # No_J16 = J16['ID'].values - # ebv = J16['E(B-V)'].values - # age_J16 = J16['age(M)'].values - # M_star_J13 = np.array([7.2,3.4,21.1,32.8,57.9,3.8,66.5,17.4,2.7,26.8,2.4,250,81.1,46.4,30.8,172.9,391.1,35.6,9.0,51.5,61.6,66.3,115,9.4,41.3,1.4,32])*1e8 - # No_J13 = np.array([3,4,15,20,23,24,25,27,28,30,31,34,35,36,43,44,47,49,50,54,58,61,62,63,64,66,67]) - # SFR_J16 = [float(SFR_J16[No_J16 == No][0]) for No in No_J13] - # # From Linhua Jiang - # klambda = 8.5 #at 2200A - # SFR_J16 = SFR_J16 * 10.**(0.4*klambda*ebv) - - if z1 == 'z2': - # eq. 28 - age = 3.223 # age of universe [Gyr], Omega_m = 0.3, Omega_lambda = 0.7, h = 0.70 - SFR_MS = 10.**((0.84-0.026*age)*np.log10(xr)-(6.51-0.11*age)) + columns = ['ID','age(M)','mass(G)','SFR_Lya','SFR_UV','E(B-V)'] + J16 = pd.read_table(d_t+'Observations/SFR_Mstar/Jiang16.txt',names=columns,skiprows=1,sep=r'\s*',engine='python') + SFR_J16 = J16['SFR_UV'].values + No_J16 = J16['ID'].values + ebv = J16['E(B-V)'].values + age_J16 = J16['age(M)'].values + M_star_J13 = np.array([7.2,3.4,21.1,32.8,57.9,3.8,66.5,17.4,2.7,26.8,2.4,250,81.1,46.4,30.8,172.9,391.1,35.6,9.0,51.5,61.6,66.3,115,9.4,41.3,1.4,32])*1e8 + No_J13 = np.array([3,4,15,20,23,24,25,27,28,30,31,34,35,36,43,44,47,49,50,54,58,61,62,63,64,66,67]) + SFR_J16 = [float(SFR_J16[No_J16 == No][0]) for No in No_J13] + # From Linhua Jiang + klambda = 8.5 #at 2200A + SFR_J16 = SFR_J16 * 10.**(0.4*klambda*ebv) simple_plot(fignum=0,xlog='y',ylog='y',xr=xr,yr=yr,fontsize=16,\ xlab='M$_{\mathrm{*}}$ [M$_{\odot}$]',ylab='SFR [M$_{\odot}$ yr$^{-1}$]',legloc=[0.01,0.75],frameon=False,\ @@ -913,14 +905,14 @@ def SFR_Mstar(): x6=xr,y6=10.**(np.log10(MS)+3.*0.2),col6='k',lw6=1,ls6=':',legend=False) ax1 = plt.gca() - # if z1 == 'z6': - # simple_plot(add='y',x1=M_star_J13[age_J16 > 30],y1=SFR_J16[age_J16 > 30],fill1='y',ma1='x',col1='r',ms1=10,mew1=2,lab1='Old $z\sim6$ LBGs/LAEs [Jiang+16]',\ - # x2=M_star_J13[age_J16 < 30],y2=SFR_J16[age_J16 < 30],fill2='y',ma2='+',col2='b',ms2=11,mew2=2,lab2='Young $z\sim6$ LBGs/LAEs [Jiang+16]') + if z1 == 'z6': + simple_plot(add='y',x1=M_star_J13[age_J16 > 30],y1=SFR_J16[age_J16 > 30],fill1='y',ma1='x',col1='r',ms1=10,mew1=2,lab1='Old $z\sim6$ LBGs/LAEs [Jiang+16]',\ + x2=M_star_J13[age_J16 < 30],y2=SFR_J16[age_J16 < 30],fill2='y',ma2='+',col2='b',ms2=11,mew2=2,lab2='Young $z\sim6$ LBGs/LAEs [Jiang+16]',legloc='upper left') # Plot models SC = ax1.scatter(M_star,SFR,marker='o',lw=2,s=64,c=ages,edgecolor='black',cmap='viridis',alpha=1,label='',zorder=10) cbar = plt.colorbar(SC,pad=0) - cbar.set_label(label='Mass-weighted stellar age [Myr]')#,size=12) # colorbar in it's own axis + cbar.set_label(label='Mass-weighted stellar age [Myr]') # colorbar in it's own axis plt.tight_layout() plt.show(block=False) @@ -942,35 +934,31 @@ def CII_SFR_z6(plot_models=True,plot_fits=True,plot_obs=True,twopanel=True,obs_g models = aux.find_model_results() for key in models.keys(): exec(key + '=models[key].values') + L_CII_tot = L_CII_GMC+L_CII_DNG+L_CII_DIG z1 = 'z6' + # Plotting parameters mpl.rcParams['xtick.minor.size'] = 4 mpl.rcParams['xtick.major.size'] = 6 mpl.rcParams['xtick.minor.width'] = 1.5 mpl.rcParams['xtick.major.width'] = 1.5 fs_labels = 20 if twopanel: fs_labels = 17 - - # Plot parameters, quickly change y or x range here + # Change y or x range here xr = 10.**np.array([0,3.05]) # SFR range yr = 10.**np.array([6.5,9.8]) # [CII] range - if not twopanel: fig = plt.figure(figsize = (12,9)) - if twopanel: fig = plt.figure(figsize = (18,8)) - - # Load [CII] luminosities - sims = cl.sim(sim_paths) - L_CII_tot = L_CII_GMC+L_CII_DNG+L_CII_DIG - - if not twopanel: - ax1 = fig.add_subplot(1,1,1) if twopanel: + fig = plt.figure(figsize = (18,8)) gs = gridspec.GridSpec(1, 2, width_ratios=[6, 5]) ax1 = plt.subplot(gs[1]) + if not twopanel: + fig = plt.figure(figsize = (12,9)) + ax1 = fig.add_subplot(1,1,1) + ax1.set_xlabel('SFR [M$_{\odot}$ yr$^{-1}$]',fontsize=fs_labels) ax1.set_ylabel('L$_{[\mathrm{CII}]}$ [L$_{\odot}$]',fontsize=fs_labels) ax1.set_xlim(xr) ax1.set_ylim(yr) - ax1.set_xlabel('SFR [M$_{\odot}$ yr$^{-1}$]',fontsize=fs_labels) xlab = ax1.get_xticks() ylab = ax1.get_yticks() ax1.set_xticklabels(xlab,fontsize=fs_labels) @@ -984,46 +972,37 @@ def CII_SFR_z6(plot_models=True,plot_fits=True,plot_obs=True,twopanel=True,obs_g # Models! # ====================================================================== - # Powerlaw fit to high-z observations [de Looze et al. 2014] - # logL_CII = np.array([4,11]) - # logSFR = -8.52+1.18*logL_CII - # fit0 = ax1.plot(10.**logSFR,10.**logL_CII,c='purple',ls='--',lw=lw) - # Powerlaw fit to z=0 metal-poor dwarf galaxies + # Powerlaw fit to z=0 metal-poor dwarf galaxies [de Looze et al. 2014] logL_CII = np.array([4,11]) logSFR = -5.73+0.80*logL_CII ax1.fill_between(10.**logSFR, 10.**(logL_CII-0.37), 10.**(logL_CII+0.37), facecolor='lightgrey', linewidth=0, alpha=0.5) fit0 = ax1.plot(10.**logSFR,10.**logL_CII,c='grey',ls='--',lw=lw) - # Powerlaw fit to z=0 starburst galaxies + # Powerlaw fit to z=0 starburst galaxies [de Looze et al. 2014] logL_CII = np.array([4,11]) logSFR = -7.06+1.00*logL_CII ax1.fill_between(10.**logSFR, 10.**(logL_CII-0.27), 10.**(logL_CII+0.27), facecolor='lightgrey', linewidth=0, alpha=0.5) fit1 = ax1.plot(10.**logSFR,10.**logL_CII,c='grey',ls='--',dashes=(1,4),lw=lw) - print('\nAverage factor below de Looze+14 metal-poor dwarfs') - - # Powerlaw fit to z=0 observations [Malhotra et al. 2001] - # L_FIR = xr*3.4e9 - # L_CII = 10.**(-2.51+np.log10(L_FIR)) - # fit1 = ax1.plot(xr,L_CII,'--',dashes=(2,5),c='teal',lw=lw) # Powerlaw fit to model galaxies at z=6 slope_models,intercept_models,slope_dev,inter_dev = aux.lin_reg_boot(np.log10(SFR),np.log10(L_CII_tot)) + # Making different line styles for the SFR-range covered by the simulations and + # the SFR-range where an extrapolation is necessary L_CII_model = 10.**(slope_models*np.log10(xr)+intercept_models) xr2 = [min(SFR),max(SFR)] L_CII_model2 = 10.**(slope_models*np.log10(xr2)+intercept_models) # dex spread around relation: dex = np.std(np.abs(np.log10(L_CII_tot)-(slope_models*np.log10(SFR)+intercept_models))) print('Dex: %s' % dex) - if not twopanel: if plot_models: fit = ax1.plot(xr,L_CII_model,ls='--',color='purple',lw=1.8) fit = ax1.plot(xr2,L_CII_model2,ls='-',color='purple',lw=1.8) fit = ax1.plot(xr,L_CII_model,ls='--',color='purple',lw=1.8) fit = ax1.plot(xr2,L_CII_model2,ls='-',color='purple',lw=1.8) - if plot_fits: # Powerlaw fit to Vallini+15, z~7 eq. 8 + sims = cl.sim(sim_paths) Zfit = np.mean(sims.Zmw()) L_CII_V15 = 10.**(7.0 + 1.2*np.log10(xr) + 0.021*np.log10(Zfit) + \ 0.012*np.log10(xr)*np.log10(Zfit) - 0.74*np.log10(Zfit)**2.) @@ -1033,7 +1012,7 @@ def CII_SFR_z6(plot_models=True,plot_fits=True,plot_obs=True,twopanel=True,obs_g # PCA fit to L_[CII](SFR,Z) Zfit2 = np.mean(Zsfr) - L_CII_PCA = 10.**(7.26+0.35*np.log10(xr)+0.38*np.log10(Zfit2)) + L_CII_PCA = 10.**(7.17+0.55*np.log10(xr)+0.23*np.log10(Zfit2)) ax1.plot(xr,L_CII_PCA,c='purple',ls='--',dashes=(6,4,1,4,1,4),lw=lw) if len(galnames) > 1: @@ -1050,19 +1029,20 @@ def CII_SFR_z6(plot_models=True,plot_fits=True,plot_obs=True,twopanel=True,obs_g rms_dex = np.sqrt(np.sum(dev**2./len(dev))) print('rms error of: '+str(rms_dex)+' dex') - print('On average this far below observed relation:') - L_CII_obs = 10.**((np.log10(SFR)+7.06)/1.00) - print(str(np.mean((L_CII_obs-L_CII_tot)/L_CII_obs)*100.)+' %') - print(str(np.mean(np.log10(L_CII_obs)-np.log10(L_CII_tot)))+' dex') - print('Power law fit this far above Vallini+15:') - L_CII_V15 = 10.**(7.0 + 1.2*np.log10(min(SFR)) + 0.021*np.log10(Zfit) + \ - 0.012*np.log10(min(SFR))*np.log10(Zfit) - 0.74*np.log10(Zfit)**2.) - L_CII = 10.**(slope_models*np.log10(min(SFR))+intercept_models) - print('at min SFR : '+str(np.log10(L_CII)-np.log10(L_CII_V15))+' dex') - L_CII_V15 = 10.**(7.0 + 1.2*np.log10(max(SFR)) + 0.021*np.log10(Zfit) + \ - 0.012*np.log10(max(SFR))*np.log10(Zfit) - 0.74*np.log10(Zfit)**2.) - L_CII = 10.**(slope_models*np.log10(max(SFR))+intercept_models) - print('at max SFR : '+str(np.log10(L_CII)-np.log10(L_CII_V15))+' dex') + if len(galnames) > 1: + print('On average this far below observed relation:') + L_CII_obs = 10.**((np.log10(SFR)+7.06)/1.00) + print(str(np.mean((L_CII_obs-L_CII_tot)/L_CII_obs)*100.)+' %') + print(str(np.mean(np.log10(L_CII_obs)-np.log10(L_CII_tot)))+' dex') + print('Power law fit this far above Vallini+15:') + L_CII_V15 = 10.**(7.0 + 1.2*np.log10(min(SFR)) + 0.021*np.log10(Zfit) + \ + 0.012*np.log10(min(SFR))*np.log10(Zfit) - 0.74*np.log10(Zfit)**2.) + L_CII = 10.**(slope_models*np.log10(min(SFR))+intercept_models) + print('at min SFR : '+str(np.log10(L_CII)-np.log10(L_CII_V15))+' dex') + L_CII_V15 = 10.**(7.0 + 1.2*np.log10(max(SFR)) + 0.021*np.log10(Zfit) + \ + 0.012*np.log10(max(SFR))*np.log10(Zfit) - 0.74*np.log10(Zfit)**2.) + L_CII = 10.**(slope_models*np.log10(max(SFR))+intercept_models) + print('at max SFR : '+str(np.log10(L_CII)-np.log10(L_CII_V15))+' dex') # Powerlaw fit to model galaxies at z=2 [Olsen+15] L_CII_model_z2 = 0.78e7*xr**1.27 @@ -1113,11 +1093,9 @@ def CII_SFR_z6(plot_models=True,plot_fits=True,plot_obs=True,twopanel=True,obs_g ax1.fill_between(10.**logSFR, 10.**(logL_CII-0.27), 10.**(logL_CII+0.27), facecolor='lightgrey', linewidth=0, alpha=0.5,zorder=0) fit1 = ax1.plot(10.**logSFR,10.**logL_CII,c='grey',ls='--',dashes=(1,4),lw=lw,zorder=0) - # pdb.set_trace() if plot_obs: add_observations_to_plot(slope_models,intercept_models,mark_reasons,mark_uplim,mark_det,z1=z1) # Legends for observations! - # get handles if plot_models: ax1.scatter([1e20,1e30],[1e20,1e30],marker='o',s=64,color='lightseagreen',label=r'S$\mathrm{\'I}$GAME (this work)$_{}$',lw=2,edgecolor='black') handles2, labels2 = ax1.get_legend_handles_labels() # remove the errorbars @@ -1130,14 +1108,12 @@ def CII_SFR_z6(plot_models=True,plot_fits=True,plot_obs=True,twopanel=True,obs_g ax1.legend(handles2, labels2, loc='upper left',fontsize=11,frameon=False,numpoints=1,borderpad=1) # Model galaxies - # L_CII_tot = L_CII_tot*0.+10.**6.6 ax1.xaxis.set_major_formatter(ScalarFormatter()) if plot_models: plot = ax1.scatter(SFR,L_CII_tot,marker='o',edgecolor='black',c=np.log10(SFRsd),s=64,lw=2,cmap='viridis',zorder=200) cbar = fig.colorbar(plot,ax=ax1,pad=0)#,cax = cbar_ax) cbar.set_label(label=getlabel('lSFRsd'),size=13)#,weight='bold') - for axis in ['top','bottom','left','right']: ax1.spines[axis].set_linewidth(1.) @@ -1599,10 +1575,8 @@ def grid_parameters(histo_color='teal',FUV=5,ISM_phase='GMC',figsize=(10,7)): if histo_color == 'colsel': ax1.legend(loc='upper left') panel += 1 - # plt.suptitle('Grid points (dashed lines) on histograms (solid lines) in diffuse gas') plt.show(block=False) plt.savefig('plots/diffuse_gas/grid_parameters/difgrid_points_on_histos'+ext_DIFFUSE+'_'+z1+'.pdf',format='pdf') # .eps for paper! - # fig.savefig('plots/dif/grids/grid_points_on_histos_paper.png', dpi=fig.dpi) #=============================================================================== """ Gas distribution plots """