Getting Started
The source code for simweights is available on github.
Prerequisites
Installation
To install from pypi run:
pip install simweights
Alternatively, if you need to install unreleased code from main you can run:
pip install git+https://github.com/icecube/simweights.git
On certain installs of python on cvmfs you might get the following error:
ModuleNotFoundError: No module named 'glob'
. If this happens you can add the
following option --no-build-isolation
to the above command.
If you want to develop simweights you can install directly with flit.
The -s
option will symlink the module into site-packages rather than copying it,
so that you can test changes without reinstalling the module:
pip install flit
git clone git@github.com:icecube/simweights.git
cd simweights
flit install [--user] -s
Basic Usage
For triggered CORSIKA or CORSIKA produced by corsika-reader
with S-Frames files use
CorsikaWeighter()
without any additional arguments:
>>> import simweights, pandas
>>> simfile = pandas.HDFStore("Level2_IC86.2016_corsika.021889.000000.hdf5", "r")
>>> flux_model = simweights.GaisserH4a()
>>> weight_obj = simweights.CorsikaWeighter(simfile)
>>> weights = weight_obj.get_weights(flux_model)
>>> print(f"Rate = {weights.sum():5.2f} Hz")
Rate = 122.84 Hz
The value returned by get_weights()
is the rate of events in Hz
For traditional CORSIKA files made with corsika-reader
you will also use
simweights.CorsikaWeighter()
, but you need to know the number of .i3
files that contributed to create this hdf5 file and pass it as the nfiles
parameter.
For neutrino-generator
you can use NuGenWeighter()
which also
requires you to know the number of files.
Flux models from nuflux can be used:
>>> import nuflux
>>> simfile = pandas.HDFStore("Level2_IC86.2016_NuMu.020878.000000.hdf5")
>>> flux_model = nuflux.makeFlux("CORSIKA_GaisserH3a_QGSJET-II")
>>> weight_obj = simweights.NuGenWeighter(simfile, nfiles=1)
>>> weights = weight_obj.get_weights(flux_model)
>>> print(f"Rate = {weights.sum():5.2e} Hz")
Rate = 1.41e-02 Hz
To weight a spectrum with a function you can also pass a callable to get_weights()
>>> weights = weight_obj.get_weights(lambda energy: 7.2e-8 * energy**-2.2)
>>> print(f"Rate = {weights.sum():5.2e} Hz")
Rate = 2.34e-05 Hz
You can also pass flux values as a numpy array with the same length as the sample
>>> fluxes = 7.2e-8 * simfile["I3MCWeightDict"]["PrimaryNeutrinoEnergy"] ** -2.2
>>> weights = weight_obj.get_weights(fluxes)
>>> print(f"Rate = {weights.sum():5.2e} Hz")
Rate = 2.34e-05 Hz
You can also pass a scalar to weight all events with the same flux. Passing
a value of 1.0
will result in the well known quantity OneWeight divided
by the number of events.
>>> OneWeight = weight_obj.get_weights(1.0)
>>> OldOneWeight = simfile["I3MCWeightDict"]["OneWeight"] / (simfile["I3MCWeightDict"]["NEvents"] / 2)
>>> (OneWeight - OldOneWeight).median()
0.0
Simulation created with genie-reader
can be weighted with GenieWeighter()
:
>>> simfile = pandas.HDFStore("genie_reader_NuE.hdf5")
>>> flux_model = nuflux.makeFlux("IPhonda2014_spl_solmax")
>>> weight_obj = simweights.GenieWeighter(simfile)
>>> weights = weight_obj.get_weights(flux_model)
>>> print(f"Rate = {weights.sum():5.2e} Hz")
Rate = 3.78e+00 Hz
Also note that these examples use pandas
. SimWeights will work equally well with
pandas
, h5py
, or pytables
.
Documentation
Full documentation is available on the IceCube Documentation Server.
Getting Help
Please direct any questions to @kjm
on the slack channel
#software.