This notebook demonstrates the method to produce the most basic simulation under consideration using the package varsim. The inputs here are:
varsim.BasePopulation. This class must implement the abstract methods of varsim.BasePopulation, and provide an index for each source and the model parametersvarsim.BaseModel. Again this subclass must implement all of the abstract methods and properties of varsim.BaseModel. The essential functionality of this class is to represent an astrophysical source as a model with model parameters, given which, this class has methods of predicting the model flux as a function of time at the top of the earth's atmosphere (ie. no sky noise included).varsim.BasicSimulation which is a subclass of the abstract base class varsim.BaseSimulation and implements concrete methods and properties necessary for the simulation. These methods and properties only use the abstract properties and methods of BaseModel and BasePopulation, and are therefore guaranteed to with any subclass.
In [1]:
import os
In [2]:
from opsimsummary import HealpixTiles, OpSimOutput
In [3]:
import numpy as np
import pandas as pd
In [5]:
from lsst.sims.photUtils import BandpassDict
In [6]:
from varsim import BasePopulation, BasicSimulation, BaseModel
In [7]:
#hptiles = HealpixTiles(nside=256,
# preComputedMap='/Users/rbiswas/data/LSST/OpSimData/healpixelized_MINION_1016_256_64_indexed.db')
In [8]:
from lsst.sims.catUtils.supernovae import SNObject
In [9]:
import analyzeSN as ans
In [10]:
import varsim
In [11]:
# The set of pointings randomly kept for convenience.
example_data = varsim.example_data
pointings_File = os.path.join(example_data, 'example_pointings.csv')
pointings = pd.read_csv(pointings_File, index_col='obsHistID')
In [12]:
# To show what this looks like
print(len(pointings))
pointings.head()
Out[12]:
In [13]:
# We did not need to have all these columns. The essential columns are ra, dec, filter,
# fivesigmadepth. It is also good to have fieldID, etc.
Here I use supernovae with the SALT model. I keep two supernovae in a completely random way without attempting to make sense. The important parts are that no matter how the Population model is implemented, it has:
modelparams : the method which takes the unique index of a supernova in the population and provides its model parameters as a dictionary. An important requirement is that the dictionary as the keys ra, dec, while the other parameters are completely user dependentidxvalues : a property which is a sequence of indices. Here the sequence is implemented as a tuple, which is perhaps how it should be for large simulations.
In [14]:
class SALTPopulation(BasePopulation):
def __init__(self):
self.x0 = [5.0e-2, 3.e-5]
self.x1 = [0, .1]
self.c = [-0.2, 0.5]
self.t0 = [59581., 59580. ]
self.z = [0.5, 0.6]
self.ra = [30., 30.]
self.dec = [-45., -45.]
def modelparams(self, idx):
return dict(x0=self.x0[idx], x1=self.x1[idx], c=self.c[idx], t0=self.t0[idx],
z=self.z[idx], ra=self.ra[idx], dec=self.dec[idx])
@property
def idxvalues(self):
return (x for x in (0, 1))
@property
def hasPositions(self):
return True
@property
def positions(self):
return (0, 1)
In [15]:
sp = SALTPopulation()
In [16]:
# The indices
idxs = tuple(sp.idxvalues)
print(idxs)
In [17]:
# The model parameters
sp.modelparams(1)
Out[17]:
and the other one
In [18]:
sp.modelparams(0)
Out[18]:
varsim.BaseModel with its methods implemented. Here I use lsst.sims.catUtils.supernovae.SNObject functionality to provide functionsminMjd and maxMjd to restrict the set of pointings over which we will build the light curve.modelFlux: concrete implementation of abstract method BaseModel.modelFlux providing the flux, given the model parameters and mjd, bandpass.
In [19]:
# We will need the LSST bandpasses. Let us load them using the catsim method
bandpassdict = BandpassDict.loadTotalBandpassesFromFiles()
In [20]:
class Model(SNObject, BaseModel):
def setModelParameters(self, **params):
paramDict = params.copy()
self.setCoords(ra=params['ra'], dec=params['dec'])
for key in ('ra', 'dec'):
params.pop(key)
self.set(**params)
self.mwEBVfromMaps()
@property
def minMjd(self):
return self.mintime()
@property
def maxMjd(self):
return self.maxtime()
def modelFlux(self, mjd, bandpassobj):
return self.catsimBandFlux(bandpassobject=bandpassobj,
time=mjd)
In [21]:
model = Model()
In [22]:
# Let us set parameters (using the method `setModelParameters` in the abstract class
# to the parameters of the first object in SALTParameters)
model.setModelParameters(**sp.modelparams(0))
In [23]:
# Now, this should predict fluxes (in maggies) given a time, and a bandpassobject
f = (model.modelFlux(mjd=59580, bandpassobj=bandpassdict['y']), model.modelFlux(mjd=59580, bandpassobj=bandpassdict['g']))
print(f, -2.5 * np.log10(f))
In [24]:
bsim = BasicSimulation(sp, model, pointings=pointings, rng=np.random.RandomState(0),
maxObsHistID=1000000, pointingColumnDict=None,
pruneWithRadius=False)
In [38]:
!rm *.hdf
In [39]:
bsim.write_simulation(phot_output='sim_phot_all.hdf', pop_output='sim_pop.hdf', method='hdf', key='sim1')
In [40]:
phot_sim_df = pd.read_hdf('sim_phot_all.hdf', key='sim1')
In [41]:
phot_sim_df
Out[41]:
In [42]:
pop_sim_df = pd.read_hdf('sim_pop.hdf', key='sim1')
In [43]:
pop_sim_df
Out[43]:
In [38]:
bsim.lc(1).index.values.size
Out[38]:
In [24]:
np.unique(bsim.lc(0).index.values).size
Out[24]:
In [25]:
bsim.lc(0).columns
Out[25]:
In [26]:
bsim.write_lc(0, 'test_0.hdf', 'hdf')
In [27]:
lc_0 = pd.read_hdf('test_0.hdf')
In [28]:
lc_0
Out[28]:
In [22]:
bsim.write_simulation('test_sim.hdf', 'hdf', key='test', clobber=False)
In [23]:
simdf = pd.read_hdf('test_sim.hdf')
In [24]:
simdf
Out[24]:
In [25]:
bsim.write_population()
In [23]:
df = pd.read_hdf('test.hdf', key='test')
In [36]:
bsim.lc(1).to_hdf('test.hdf', key='test', mode='w', append=True, format='t')
In [37]:
bsim.lc(0).to_hdf('test.hdf', key='test', mode='a', append=True, format='t')
In [38]:
df = pd.read_hdf('test.hdf', key='test')
In [39]:
df
Out[39]:
In [23]:
bsim.write_lc(0, 'test_0.csv', 'csv', key=None)
In [25]:
dfcsv = pd.read_csv('test_0.csv')
In [25]:
bsim.write_simulation('test.hdf', 'hdf', key=None)
In [40]:
df = pd.read_hdf('test.hdf', key='0')
In [41]:
df.to_sql('')
Out[41]:
In [24]:
bsim.write_population('pop.hdf', method='hdf')
In [ ]:
sp.idxvalues
In [28]:
dfcsv.equals(df, )
Out[28]:
In [24]:
bsim = BasicSimulation(population=sp, model=Model, pointings=)
In [21]:
bsim.write_lc(0, 'lc_0.hdf', 'hdf')
In [ ]:
pd.DataFrame.to_hdf()
In [28]:
bsim.write_simulation('test.hdf', 'hdf', key='test')
In [29]:
simdf = pd.read_hdf('test.hdf', key='test')
In [30]:
simdf
Out[30]:
In [ ]:
sp.positions
In [ ]:
model = Model()
In [ ]:
sp.modelparams(0)
In [ ]:
sp.modelparams(0)
In [ ]:
model.setModelParameters(**sp.modelparams(0))
In [ ]:
model.sn.SNstate
In [ ]:
sn = SNObject()
In [ ]:
params
In [ ]:
sn.setCoords(params['ra'], params['dec'])
In [ ]:
sn.set(**params)
In [ ]:
sn.mwEBVfromMaps()
In [ ]:
sn.SNstate
In [ ]:
model = Model()
In [ ]:
params = sp.modelparams(0).copy()
In [ ]:
for key in ('ra', 'dec'):
params.pop(key)
In [ ]:
params
In [ ]:
sn.set(**params)
In [ ]:
sn.SNstate
In [ ]:
model.setModelParameters(**sp.modelparams(0))
In [ ]:
model.sn.SNstate
In [ ]:
bandpassobj = bandpassdict['y']
In [ ]:
model.setModelParameters(**sp.modelparams(0))
In [ ]:
model.modelFlux(59580., bandpassobj=bandpassobj)
In [ ]:
bsim = BasicSimulation(sp, model, pointings, rng=np.random.RandomState(0),
maxObsHistID=1000000, pointingColumnDict=None,
timeRange=
pruneWithRadius=False)
In [ ]:
bsim.model.maxMjd
In [ ]:
bsim.model.r
In [ ]:
scrtch__
In [ ]:
opsout = _OpSimOutput.fromOpSimDB('/Users/rbiswas/data/LSST/OpSimData/minion_1016_sqlite.db',
)
In [ ]:
pointings = opsout.summary.query('expMJD < 59580.04')
In [ ]:
pointings.expMJD.min()
In [ ]:
len(pointings)
In [ ]:
pointings.to_csv('check_poinx``tings.csv')
In [ ]:
!head check_pointings.csv
In [ ]:
from var
In [ ]:
from varsim import BaseExcep
In [ ]:
import pandas as pd
In [ ]:
import numpy as np
In [ ]:
df = pd.DataFrame()
In [ ]:
df['objid'] = np.arange(10)
In [ ]:
df
In [24]:
280 *2 /4 + 20
Out[24]:
In [116]:
from collections import namedtuple, OrderedDict as odict
In [129]:
keys = 'abcde'
vals = np.arange(5)
idx = 'pq'
s1 = (dict( (k, v) for (i, k,v) in zip(idx, keys, vals)))
s2 = (dict( (k, v) for (i, k,v) in zip(idx, keys, vals)))
In [117]:
s = namedtuple(idx, s1)
In [123]:
s = list()
_ = list(s.append(x) for x in (s1, s2))
In [127]:
s1
Out[127]:
In [102]:
xx = list()
In [107]:
xx = tuple(s)
In [103]:
xx.append(s)
In [133]:
idx = (l for l in idx)
In [135]:
df = pd.DataFrame(s, index=idx)
In [136]:
df
Out[136]:
In [57]:
df = pd.DataFrame()
In [58]:
type(s)
Out[58]:
In [56]:
df.append(s)
In [ ]: