In [1]:
%pylab inline


Populating the interactive namespace from numpy and matplotlib

In [2]:
import urllib
import os
import time

In [3]:
datadir = 'sage_output/' # Change me to whatever you like
# Retrieve output of SAGE from Mini Millennium with fiducial parameters if directory is empty
# Run this for an empty directory first to get the necessary plotting scripts and test the notebook!
if not os.path.isfile(datadir+'model_z0.000_0'):
	urllib.urlretrieve('https://bitbucket.org/arhstevens/sage-mm-output/get/5a222f5b1053.zip', 'sage_mm.zip')
	os.system('unzip -j sage_mm.zip -d '+datadir)
	os.system('rm sage_mm.zip')
if os.path.isfile(datadir+'plot_read_routines.py'):
    os.system('cp '+datadir+'plot_read_routines.py .')

In [4]:
# Import routines for reading and plotting
import plot_read_routines as pr

In [5]:
# Read SAGE data
M = pr.read_sagesnap(datadir+'model_z0.000')

In [6]:
# Plot settings
fsize = 28
matplotlib.rcParams.update({'font.size': fsize, 'xtick.major.size': 10, 'ytick.major.size': 10, 'xtick.major.width': 1, 'ytick.major.width': 1, 'ytick.minor.size': 5, 'xtick.minor.size': 5, 'axes.linewidth': 1, 'text.usetex': True, 'font.family': 'serif', 'font.serif': 'Times New Roman', 'figure.figsize': (9.0, 9.0*2/3), 'legend.numpoints': 1})

In [7]:
h = 0.73 # Choose value (default is for the Millennium simulation)
Lbox = 62.5/h # Box length (default is for the Millennium simulation)

In [8]:
# Plot stellar mass function
plt.figure()
pr.massfunction_obsdata(h)
pr.massfunction(M.StellarMass*1e10/h, Lbox, Nbins=30)
plt.yscale('log', nonposy='clip')
plt.axis([8.5, 11.9, 1e-6, 0.2])
plt.xlabel(r'$\log_{10}(m_*\ [\mathrm{M}_{\odot}])$')
plt.ylabel(r'$\Phi\ [\mathrm{Mpc}^{-3}\ \mathrm{dex}^{-1}]$')
plt.legend(fontsize=fsize-6, loc='best', frameon=False)


Out[8]:
<matplotlib.legend.Legend at 0x10d35a410>

In [9]:
# Plot baryonic Tully-Fisher relation
M_filt = (M.Type==0) * (M.StellarMass+M.ColdGas > 0) * (M.BulgeMass > 0.1*M.StellarMass) * (M.BulgeMass < 0.5*M.StellarMass)
V_rot, M_bary = M.Vmax[M_filt], (M.StellarMass[M_filt]+M.ColdGas[M_filt])*1e10/h
plt.figure()
pr.btf_obsdata(h)
pr.contour(np.log10(V_rot), np.log10(M_bary), Nbins=50, Nlevels=5)
plt.axis([1.4, 3.1, 7.1, 11.4])
plt.xlabel(r'$\log_{10} (V_{\rm max}\ [\mathrm{km\ s}^{-1}])$')
plt.ylabel(r'$\log_{10} (m_*+m_{\rm cold}\ [\mathrm{M}_{\bigodot}])$')
plt.plot([-1,-1], [-1,-2], 'k-', lw=2, label=r'\textsc{sage}') # Just for legend
plt.legend(fontsize=fsize-6, loc='lower right', frameon=False)


Out[9]:
<matplotlib.legend.Legend at 0x10c493bd0>

In [10]:
# Plot black hole -- bulge mass relation
M_filt = (M.BlackHoleMass/h > 1e-5) * (M.BulgeMass/h > 10**-1.5)
M_BH, M_bulge = M.BlackHoleMass[M_filt]*1e10/h, M.BulgeMass[M_filt]*1e10/h
plt.figure()
pr.bhbulge_obsdata(h)
pr.contour(np.log10(M_bulge), np.log10(M_BH), Nbins=50, Nlevels=5)
plt.axis([8.5, 12.1, 5.1, 9.4])
plt.xlabel(r'$\log_{10} (m_{\rm *,bulge}\ [\mathrm{M}_{\bigodot}])$')
plt.ylabel(r'$\log_{10} (m_{\rm BH}\ [\mathrm{M}_{\bigodot}])$')
plt.plot([-1,-1], [-1,-2], 'k-', lw=2, label=r'\textsc{sage}') # Just for legend
plt.legend(fontsize=fsize-6, loc='best', frameon=False, numpoints=1)


Out[10]:
<matplotlib.legend.Legend at 0x10d04f190>

In [11]:
# Plot the stellar mass -- gas metallicity relationship
M_filt = (M.Type==0) * (M.ColdGas > 0) * (M.MetalsColdGas > 0) * (M.StellarMass > 0.01)
M_star, OonH = M.StellarMass[M_filt]*1e10/h, M.MetalsColdGas[M_filt]/M.ColdGas[M_filt]/0.02
plt.figure()
pr.massmet_obsdata(h)
#pr.contour(np.log10(M_star), np.log10(OonH)+9, range=[[8,11],[7.5,10.5]], Nbins=50, Nlevels=5)
m_plot, met_high, met_med, met_low = pr.percentiles(np.log10(M_star), np.log10(OonH)+9, xrange=[8,12])
plt.errorbar(m_plot, met_med, yerr=[met_high-met_med, met_med-met_low], color='k', lw=2, ms=0, ls='none', capsize=0)
plt.xlabel(r'$\log_{10}(m_*}\ [\mathrm{M}_{\bigodot}])$')
plt.ylabel(r'$12 + \log_{10}(\mathrm{O} / \mathrm{H})$')
plt.plot([-1,-1], [-1,-2], 'k|', ms=20, mew=2, label=r'\textsc{sage}') # Just for legend
plt.axis([8.1, 11.9, 8, 9.5])
plt.legend(fontsize=fsize-6, loc='lower right', frameon=False)


Out[11]:
<matplotlib.legend.Legend at 0x10da54210>

In [12]:
# Plot the quiescent fraction of galaxies
plt.figure()
pr.quiescent_obs(h)
pr.quiescent(M.StellarMass*1e10/h, (M.SfrBulge+M.SfrDisk), h=h, Nbins=15)
plt.axis([9.05, 11.4, 0, 1])
plt.xlabel(r'$\log_{10}(m_*\ [\mathrm{M}_{\bigodot}])$')
plt.ylabel(r'Quiescent fraction')
plt.legend(fontsize=fsize-6, loc='lower right', frameon=False)


Out[12]:
<matplotlib.legend.Legend at 0x10db19fd0>

In [ ]: