In [26]:
%load_ext autoreload
%matplotlib inline
import matplotlib.pyplot as plt
import numpy as np
import os
import pandas as pd
import sncosmo
import gedankenLSST
import seaborn as sns
sns.set()
In [27]:
from lsst.sims.photUtils import BandpassDict
In [28]:
lsst_bp = BandpassDict.loadTotalBandpassesFromFiles()
In [29]:
# sncosmo Bandpasses required for fitting
throughputsdir = os.getenv('THROUGHPUTS_DIR')
from astropy.units import Unit
bandPassList = ['u', 'g', 'r', 'i', 'z', 'y']
banddir = os.path.join(os.getenv('THROUGHPUTS_DIR'), 'baseline')
for band in bandPassList:
# setup sncosmo bandpasses
bandfname = banddir + "/total_" + band + '.dat'
# register the LSST bands to the SNCosmo registry
# Not needed for LSST, but useful to compare independent codes
# Usually the next two lines can be merged,
# but there is an astropy bug currently which affects only OSX.
numpyband = np.loadtxt(bandfname)
print band
sncosmoband = sncosmo.Bandpass(wave=numpyband[:, 0],
trans=numpyband[:, 1],
wave_unit=Unit('nm'),
name=band)
sncosmo.registry.register(sncosmoband, force=True)
In [ ]:
# Expected number of galaxies in Twinkles area
In [14]:
# Galaxies with z < 1.2
galaxyBase = pd.read_csv('/Users/rbiswas/doc/projects/supernovae/LSST/sgi/src/fatboy_galaxyBase.csv')
In [15]:
zmids = 0.5 * (galaxyBase['minz'] + galaxyBase['maxz'])
numGals = galaxyBase['numGals']* (10./60./ 4.0)**2
# Total number of galaxies
numGals.sum()
Out[15]:
In [16]:
fig, ax = plt.subplots()
ax.plot(zmids,numGals ,'o')
Out[16]:
In [17]:
lsstchar = gedankenLSST.LSSTReq
In [18]:
lsstchar['meanNumVisits'] = pd.Series(np.repeat(3650.,6), index=['u','g','r','i','z','y'])
In [19]:
lsstchar['meanNumVisits']
Out[19]:
In [20]:
sn = gedankenLSST.GSN_Obs(mjd_center=49530., lsstrequirements=lsstchar)
In [21]:
df = sn.summary
df[df['filter'] == 'u'].hist('night',bins=80)
Out[21]:
In [30]:
s = gedankenLSST.SNObs(summarydf=df, t0=49530, lsst_bp=lsst_bp)
In [55]:
a = []
for z in zmids.values:
s.snState = {'z': z}
lc = s.lightcurve
totalEpochs = len(lc)
highSNRlc = lc.query('SNR > 5')
highSNREpochs = len(highSNRlc)
highSNREpochs_u = len(highSNRlc.query("filter == 'u'"))
highSNREpochs_g = len(highSNRlc.query("filter == 'g'"))
highSNREpochs_r = len(highSNRlc.query("filter == 'r'"))
highSNREpochs_i = len(highSNRlc.query("filter == 'i'"))
highSNREpochs_z = len(highSNRlc.query("filter == 'z'"))
highSNREpochs_y = len(highSNRlc.query("filter == 'y'"))
a.append([z, highSNREpochs, highSNREpochs_u, highSNREpochs_g, highSNREpochs_r, highSNREpochs_i, highSNREpochs_z,
highSNREpochs_y, totalEpochs, -2.5 * np.log10(s.SN.get('x0'))])
In [70]:
clear(df)
In [76]:
df = pd.DataFrame(a, columns=['redshift', 'highSNREpochs', 'u', 'g', 'r', 'i', 'z', 'y', 'totalEpochs', 'mB'])
In [108]:
df['frac'] = df.highSNREpochs / df.totalEpochs
In [112]:
df['NumSNperzBin'] = 10 * 3650. / 80. / df.frac
In [114]:
df['numGals'] = numGals
In [116]:
df['SNperGals'] = df['NumSNperzBin'] / df['numGals']
In [117]:
df
Out[117]:
In [118]:
df.to_csv('Twinkles_roughPlan.csv')
In [107]:
fig, ax, ay = pu.settwopanel(setdifflimits=(-29.805, -29.795))
ay.plot(df.redshift, df.mB - cosmo.distmod(df.redshift.values).value)
ax.plot(df.redshift, cosmo.distmod(df.redshift.values).value)
ax.plot(df.redshift, df.mB +29.8, 'o')
Out[107]: