Usage¶
IceTray Service¶
The intended usage of rpdf is to provide the likelihood service for a Gulliver reconstruction. As shown in the following example:
tray.AddService("I3RecoLLHFactory", "spe_likelihood",
InputReadout="InIcePulses",
JitterTime=15*I3Units.ns,
NoiseProbability=10.*I3Units.hertz
)
tray.Add("I3SimpleFitter",
OutputName="SPESingleFit",
SeedService="seed",
Parametrization="simpletrack",
LogLikelihood="spe_likelhood",
Minimizer="minuit",
)
In the example an instance of I3RecoLLH
with the name
"spe_likelihood"
is placed in the tray’s
context by I3RecoLLHFactory
. This example reads an
I3RecoPulsesSeriesMap
from the frame and uses a jitter time of
15 ns and a noise rate of 10 Hz. A Gulliver reconstruction
I3SimpleFitter
is then added to the tray and passed the name of the
rpdf service to the LogLikelihood
parameter.
Python¶
The low level functions are also available via a python interface. An example calculation is below:
from icecube import rpdf
icemodel = rpdf.H2
peprob = rpdf.FastConvolutedPandel()
jitter = 15*I3Units.ns
noise = 10*I3Units.hz
llhtot = 0
for omkey, pulseseries in self.pulses:
#find the position of the OM
ompos = geometry[omkey].position
#calculate the geometrical parameters tgeo and deff
geo_params = rdf.muon_geometry(ompos,part,icemodel)
#residual time is the pulse time minus tgeo
t_res= pulseseries[0].time-geo_params.first
#calculate the convoluted Pandel function of t_res
llh = peprob.pdf(t_res,geo_params.second,jitter,icemodel)
#add noise and sum the log of the likelihoods
llhtot += math.log(self.noise + llh)
print llhtot
This is equivalent to calculating the SPE1st likelihood in I3RecoLLH
.
A easier method if you want to calculate the likelihood of a single event
hypothesis is to instantiate a I3RecoLLH
object in python land:
recollh = rpdf.I3RecoLLH("SRTHVInIcePulses","SPE1st","GaussConvoluted",
15*I3Units.ns,10.*I3Units.hertz,rpdf.H2)
particle = phys["SPEFitSingle_TWHV"]
hypothesis = gulliver.I3EventHypothesis(particle)
recollh.set_geometry(i3geometry)
recollh.set_event(physics_frame)