In [1]:
%pylab inline
from astropy.io import fits
from astropy.table import Table
from simqso.sqrun import buildWaveGrid,load_sim_output
from simqso import hiforest
from simqso.sqmodels import WP11_model
import bossqsos
In [2]:
M1450 = linspace(-30,-22,20)
zz = arange(0.7,3.5,0.5)
ple = bossqsos.BOSS_DR9_PLE()
lede = bossqsos.BOSS_DR9_LEDE()
for z in zz:
if z<2.2:
qlf = ple if z<2.2 else lede
plot(M1450,qlf(M1450,z),label='z=%.1f'%z)
legend(loc='lower left')
xlim(-21.8,-30.2)
xlabel("$M_{1450}$")
ylabel("log Phi")
Out[2]:
In [3]:
_ = bossqsos.qsoSimulation(bossqsos.simParams,saveSpectra=True)
In [4]:
wave,qsos = load_sim_output('boss_dr9qlf_sim','.')
the table of simulated quasars, including redshift, luminosity, synthetic flux/mags in nine bands, and "observed" photometry with errors included.
also includes details of the model inputs for each quasar: slopes is the set of broken power law slopes defining the continuum, emLines is the set of Gaussian parameters for each emission line (wave, EW, sigma) measured in the rest frame.
In [5]:
qsos[::40]
Out[5]:
the distribution in g-band magnitude:
In [6]:
_ = hist(qsos['obsMag'][:,1],linspace(17,22,20),log=True)
color-color diagram from observed magnitudes, including errors:
In [7]:
scatter(qsos['obsMag'][:,0]-qsos['obsMag'][:,1],qsos['obsMag'][:,1]-qsos['obsMag'][:,2],
c=qsos['z'],cmap=cm.autumn_r,alpha=0.7)
colorbar()
xlabel('u-g')
ylabel('g-r')
xlim(-0.75,3)
ylim(-0.5,1.5)
Out[7]:
the list of emission lines in the model:
In [8]:
qsodatahdr = fits.getheader('boss_dr9qlf_sim.fits',1)
for i,n in enumerate(qsodatahdr['LINENAME'].split(',')):
print('%d:%s, '% (i,n,),end=" ")
print()
broad CIV equivalent width, displaying the Baldwin Effect:
In [9]:
scatter(qsos['absMag'],qsos['emLines'][:,13,1],c=qsos['z'],cmap=cm.autumn_r)
colorbar()
xlabel("$M_{1450}$")
ylabel("CIV equivalent width $\AA$")
Out[9]:
In [10]:
figure(figsize=(14,4))
plot(wave/1e4,qsos['spec'][0])
yscale('log')
xlabel('wave [micron]')
Out[10]:
zoom in on the lyman alpha - CIV region:
In [11]:
figure(figsize=(14,4))
plot(wave,qsos['spec'][20])
xlim(3500,7500)
title('$z=%.3f$'%qsos['z'][20])
Out[11]:
an example of the forest transmission spectra at R=30,000 (the native resolution for the monte carlo forest spectra):
In [2]:
# XXX WARNING -- an ugly hack is needed here. Internally, a table of Voigt profiles is generated
# at startup in order to speed the forest spectra generation. This table is defined in terms of
# the wave dispersion the first time a simulation is run. Here we are changing the wavelength
# model, and thus before executing the next cells you must restart the kernel and execute only
# the first cell.
np.random.seed(12345)
wave = buildWaveGrid(dict(waveRange=(3500,4800),SpecDispersion=30000))
forest = hiforest.IGMTransmissionGrid(wave,WP11_model,1)
T = forest.next_spec(0,2.9)
In [3]:
figure(figsize=(14,4))
plot(wave,T)
Out[3]:
In [4]:
figure(figsize=(14,4))
plot(wave,T)
xlim(4300,4800)
Out[4]:
In [ ]: