In [1]:
from LSSTmetrics import PerSNMetric
from LSSTmetrics.efficiencyTable import EfficiencyTable
In [2]:
import LSSTmetrics as metrics
print metrics.__file__
In [3]:
from opsimsummary import summarize_opsim as oss
In [4]:
import pandas as pd
import sncosmo
In [5]:
%matplotlib inline
import matplotlib.pyplot as plt
import os
import numpy as np
import copy
In [6]:
from lsst.sims.photUtils import BandpassDict
In [7]:
# Catsim bandpasses
lsst_bp = BandpassDict.loadTotalBandpassesFromFiles()
In [8]:
# 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 [9]:
# You will have to download this file and change the value of the variable opsimHDF to the absolute path of this file on your disk: http://lsst.astro.washington.edu/simdata/SN_data/storage.h5
opsimHDF = os.path.join(os.getenv('HOME'), 'data', 'LSST', 'OpSimData', 'storage.h5')
summarydf = pd.read_hdf(opsimHDF, 'table')
# df = df.query('propID == [364, 366]')
In [10]:
# Look at the output
summarydf.head()
Out[10]:
In [11]:
# summarydf = summarydf.query('propID == [364, 366]').query('night < 365')
In [12]:
# Create the summary instance
so = oss.SummaryOpsim(summarydf)
In [13]:
_ = so.cadence_plot(fieldID=290, mjd_center=49540)
In [14]:
_ = so.cadence_plot(fieldID=309, mjd_center=49540)
In [16]:
# 1st Season for a WFD field
_ = so.cadence_plot(fieldID=309, sql_query='night < 366')
In [17]:
# 1st Season for a deep field
_ = so.cadence_plot(fieldID=290, sql_query='night < 366')
In [18]:
qm = PerSNMetric(fieldID=309, t0=49540, summarydf=summarydf, lsst_bp=lsst_bp)
In [19]:
qm.lcplot()
Out[19]:
In [20]:
# Coadd the light curve
qm.lcplot(nightlyCoadd=True)
Out[20]:
In [34]:
#Let us see how this looks compared to fiveSigmaDepth
xx = qm.lightcurve[['time', 'band','flux', 'fluxerr', 'fiveSigmaDepth']].copy()
In [26]:
xx['mag'] = -2.5 * np.log10(xx['flux'])
In [29]:
xx.groupby('band').get_group('i')
Out[29]:
In [29]:
# qm.writeLightCurve('coadded_lc.dat',nightlyCoadd=True)
In [35]:
cm = PerSNMetric(fieldID=290, t0=49540, summarydf=summarydf, lsst_bp=lsst_bp)
In [37]:
_ = cm.lcplot()
In [39]:
_ = cm.lcplot(nightlyCoadd=True)
In [43]:
yy = cm.lightcurve[['time', 'band','flux', 'fluxerr', 'fiveSigmaDepth']].copy()
In [45]:
yy['mag'] = -2.5 * np.log10(yy.flux)
In [46]:
yy
Out[46]:
In [47]:
cm.SN.SNstate
Out[47]:
In [48]:
# Write Light Curve to Disk
In [49]:
cm.writeLightCurve('DDF_290_lc.txt')
In [50]:
!head DDF_290_lc.txt
In [ ]: