How to implement your own seed service¶
Overview¶
A log-likelihood reconstructed must be seeded with a hypothesis (e.g.
the result of a first-guess reconstruction) relatively close to the
global minimum in the likelihood space. For Gulliver-based
reconstructions, a hypothesis is represented as an
I3EventHypothesis
object. Obviously, the seed hypothesis
must be an object of the same type. The I3SeedServiceBase
class collects one or several first-guess reconstruction results
(I3Particle
objects) from an event, forms seed event
hypotheses, and hands a seed to I3Gulliver
via the
I3SeedServiceBase::GetSeed()
method. Usually, this happens
inside a Guilliver-based module like I3SimpleFitter
(gulliver-modules
).
For most trivial event hypotheses that do not use
I3EventHypothesis::nonstd
, you can use the
I3BasicSeedSerive
class (lilliput
) if you want to
implement your own Gulliver reconstruction. This seed service collects a
list of first-guess reconstruction results and provides different
options to correct the seed’s vertex position and time (for infinite
tracks). If your event hypothesis needs other event data that is not
provided by the I3Particle
class, you definitely need to
implement your own seed service. Projects that define their own seed
service are: toprec.
API¶
All seed services have to be derived from I3SeedServiceBase
and should implement the following pure virtual methods:
- virtual unsigned int SetEvent (const I3Frame &f)
Provide event data: get the first-guess tracks; fill blank data members; do tweaks like space and time corrections.
The generated seeds should all have an OK fit status and non-NAN position, direction, and time data members.
- virtual I3EventHypothesis GetSeed (unsigned int iseed) const
Return a seed as a deep copy of some internally stored seed because
I3Gulliver
will most likely change its values.- virtual I3EventHypothesis GetDummy () const
Return a dummy seed; useful in case all first guesses fail. In particular if the seeds normally have a non-NULL
nonstd
datamember, then this should be initialized in this dummy seed as well.- virtual I3EventHypothesis GetCopy (const I3EventHypothesis &eh) const
Return a (deep) copy of an existing event hypothesis.
- virtual void Tweak (const I3EventHypothesis &eh) const
When a seed service generates a seed, it may perform some tweaks or transformations on it, intended to facilitate the minimization. The generic example is the vertex correction for an infinite track.
- virtual void FillInTheBlanks (const I3EventHypothesis &eh) const
Many first guess algorithms do not provide the complete description for the kind of event we are trying to reconstruct. A simple example is the energy of a cascade or muon. The seed service might take care of filling in a sensible initial value for such a quantity.
Example¶
Example for an I3SeedServiceBase
implementation in Python
that seeds a reconstruction with different directions based on a
healpy grid:
import healpy
import numpy
import icecube
import icecube.icetray
import icecube.dataclasses
import icecube.gulliver
class HEALPixSeedService(icecube.gulliver.I3SeedService):
"""Seed service implementation based on HEALPix
Seed reconstruction with directions based on a full-sky HEALPix
grid.
"""
def __init__(self, vertex_fit="", nside=32):
super(HEALPixSeedService, self).__init__()
self.key = vertex_fit
self.nside = nside
self.npix = healpy.nside2npix(self.nside)
self.pixels = numpy.arange(self.npix)
def SetEvent(self, frame):
self.vertex_fit = frame[self.key]
return self.npix
def GetSeed(self, i):
eh = icecube.gulliver.I3EventHypothesis(self.vertex_fit)
eh.particle.shape = icecube.dataclasses.I3Particle.InfiniteTrack
eh.particle.speed = icecube.dataclasses.I3Constants.c
eh.particle.length = numpy.nan
eh.particle.dir = icecube.dataclasses.I3Direction(
*healpy.pix2ang(self.nside, self.pixels[i]))
return eh
def GetCopy(self, eh):
return icecube.gulliver.I3EventHypothesis(eh.particle)
def GetDummy(self):
p = icecube.dataclasses.I3Particle()
return icecube.gulliver.I3EventHypothesis(p)
def Tweak(self, eh):
return
def FillInTheBlanks(self, eh):
return