gulliver-modules C++ API Reference¶
-
class I3IterativeFitter : public I3ConditionalModule¶
- #include <I3IterativeFitter.h>
Gulliver-based module to perform iterative reconstructions.
The I3IterativeFitter is a reconstruction module which specializes in reconstructing the direction of a tracks or showers. The module works very similar to I3SimpleFitter: it uses the Gulliver toolkit to reconstruct e.g. a muon track, using the likelihood, parametrization of fit variables and the minimizer algorithm from external (Gulliver) services.
For each seed (as provided by a seed service, see I3SeedServiceBase in gulliver or I3BasicSeedService in lilliput) first regular fit is performed, just like I3SimpleFitter does. This fit is then redone several times, each time seeding the minimizer with the result of the previous iteration except that the direction is replaced with some (quasi-) random other direction. In most iterations the result will be worse than for the fit with the original seed, but for some events the original fit may have been found just a local minimum, while during the iterations the global minimum of -log(L) is found (or at least a deeper local minimum).
See the documentation of I3SimpleFitter module for more information about the configuration of the services, and about how the results of the reconstruction are stored in the frame.
The number of iterations is set with the “NIterations” option. The generated directions can be restricted to fall within a limited range of cos(zenith), specified with the “CosZenithRange” option. You can either use a random service (specified with the “RandomService” option) purely random directions for the iterations, or use a socalled Sobol/Niederreiter2 sequence (see below) to generate a pseudorandom grid of directions. You will get the Sobol sequence if you leave the “RandomService” option empty, or if you fill out “SOBOL”. You will get the Niederreiter2 sequence if you fill out “NIEDERREITER2”.
About the Sobol and Niederreiter2 Grids
A 2-dimensional Sobol grid is a “quasirandom” sequence of pairs of numbers between 0 and 1. These numbers are uniformly distributed over the unit square (the square with corners at (0,0), (0,1), (1,0) and (1,1)). The sequence is designed to maximally “cover” the phase space. It’s actually a pretty simple sequence, the first few points are:
(0.500, 0.500)
(0.750, 0.250)
(0.250, 0.750)
(0.375, 0.375)
(0.875, 0.875)
(0.625, 0.125)
(0.125, 0.625)
….Using the Sobol sequence (instead of a proper random number generator) to generate the directions of the iterations of the iterative likelihood has two advantages: (1) the already mentioned uniform coverage of the phase space (with a true random number generator, it can in principle happen that all generated directions for that event are in the lower half of the phase space) (2) the generated directions are much more easily reproducible (a problem with truely random iterative fit in AMANDA was that redoing the fit sometimes led to vastly different results).
We compute the sequence Sobol grid points only once and use the same sequence for every seed in every event. However, the actual directions are chosen to be (kind of) relative to the original seed direction. This is done by rescaling cos(zenith) and azimuth such that their range runs from 0 to 1. For a given seed we add its rescaled cos(zenith) and azimuth to all grid points of the sequence, and if necessary they are moved back into the unit square by taking their coordinates modulo 1.0. These translated grid points are then scaled back to regular cos(zenith) and azimuth values.
A “Niederreiter2” sequence another quasirandom number generator, just like Sobol. Both sequences are implemented in GSL (
gsl_qrng_sobol
andgsl_qrng_niederreiter_2
).See also
Public Functions
-
~I3IterativeFitter()¶
-
void Configure()¶
-
void Geometry(I3FramePtr frame)¶
-
void Physics(I3FramePtr frame)¶
-
void Finish()¶
Private Functions
-
I3IterativeFitter()¶
-
I3IterativeFitter(const I3IterativeFitter &source)¶
-
I3IterativeFitter &operator=(const I3IterativeFitter &source)¶
-
void InitGridPoints(bool sobol)¶
Initialize the Sobol/Niederreiter grid with
nIterations_
grid points- Parameters:
sobol – [in] Use Sobol (true) or Niederreiter2 (false)
-
std::pair<double, double> GetNormDir(const I3Direction &dir)¶
Rescale a direction to a point in the unit square (applied to seed direction)
-
void SetGridDir(I3ParticlePtr particle, const std::pair<double, double> &offnorm, unsigned int index)¶
Set the direction based on a Sobol/Niederreiter grid point and an offset (based on seed direction).
-
I3LogLikelihoodFitPtr BestFit(unsigned int nseeds)¶
-
void Fit(int seed, int iteration, I3LogLikelihoodFitPtr fit, std::vector<I3LogLikelihoodFit> &goodFits)¶
Convenience method: do fit, and store it in a vector (if successful)
- SET_LOGGER ("I3IterativeFitter")
Private Members
-
I3GulliverPtr fitterCore_¶
The core Gulliver object for basic tracks.
-
I3RandomServicePtr randomService_¶
Random service to generate random directions.
-
I3EventLogLikelihoodBasePtr likelihood_¶
Event loglikelihood calcuation service.
-
I3ParametrizationBasePtr parametrization_¶
Track/shower/anything parametrization service.
-
I3MinimizerBasePtr minimizer_¶
Minimizer service.
-
I3SeedServiceBasePtr seedService_¶
Seed collection and preparation service.
-
I3SeedServiceBasePtr tweakService_¶
Seed preparation service for the (pseudo) random iteration seeds.
-
unsigned int nIterations_¶
Number of random iterations on each seed.
-
std::vector<double> cosZenithRange_¶
Range of cosine zenith angles within which to generate directions.
-
std::vector<std::pair<double, double>> gridPoints_¶
Sobol or Niederreiter grid (used to generate pseudo-random directions)
-
unsigned int eventNr_¶
Number of Physics calls, used in log messages.
-
unsigned int nSeeds_¶
Number of Physics calls with at least one good seed.
-
unsigned int nSuccessFits_¶
Number of seeds resulting in a successful fit.
-
unsigned int nSuccessEvents_¶
Number of Physics calls with a successful fit.
-
~I3IterativeFitter()¶
-
class I3LogLikelihoodCalculator : public I3ConditionalModule¶
- #include <I3LogLikelihoodCalculator.h>
Gulliver-based module to compute the log-likelihood for a given track.
This module does not do any reconstruction, it just computes the loglikelihood for a given I3Particle and a likelihood service.
It stores the same set of loglikelihood fitparams that I3SimpleFitter or I3IterativeFitter would store when run with a hypothetical minimizer which would evaluate the likelihood function only once (for the seed track) and then declares convergence.
This module is only for the simplest cases of an I3Particle object, no composite or otherwise complicated tracks which require a seed service and/or additional non-standard objects to specify the physics variables of the event hypothesis.
In a later stage I might add that extra functionality, for now this simple functionality is what I need.
By default the result is an object of type I3LogLikelihoodFitParams (defined in the gulliver project), which is stored in the frame under the modules own name (that is, the second argument of the AddModule(module,name) command in a typical icetray python script). If you set the KeepFitName flag and your input fit has the name “FooBar” then the result will be stored as “FooBar_FitParams” instead, as if the “FooBar” was the result of a I3SimpleFitter or I3IterativeFitter.
There are at least two use cases for this module:
A loglikelihood fit was performed in some previous processing level (e.g. online filtering, for the Moon) but the loglikelihood fitparams were not kept in the subsequent reading & writing from/to i3 files. Only the I3Particle was kept. With this module you can re-compute the the likelihood without doing a possibly time-costly fit. The number of minimizer steps (or rather: the number of llh function evaluations) cannot be retrieved anymore.
If you would like to know the likelihood of some track without doing a full reconstruction, e.g. to see the difference between the seed and the full fit.
- Todo:
Extend this module such that you can also use it with nonstandard event hypotheses (which need an additional object besides the I3Particle in order to store the relevant physics variables; and/or need a seed service to tweak and such).
Public Functions
-
inline ~I3LogLikelihoodCalculator()¶
destructor
-
void Configure()¶
configure (get & check configurables)
-
void Geometry(I3FramePtr frame)¶
handle geometry
-
void Physics(I3FramePtr frame)¶
do a reconstruction
-
void Finish()¶
say bye
Private Functions
-
I3LogLikelihoodCalculator()¶
-
I3LogLikelihoodCalculator(const I3LogLikelihoodCalculator &source)¶
inhibited
-
I3LogLikelihoodCalculator &operator=(const I3LogLikelihoodCalculator &source)¶
inhibited
-
bool CheckParticle(I3ParticleConstPtr p)¶
inhibited
paranoid check that input particle is healthy
- SET_LOGGER ("I3LogLikelihoodCalculator")
counts nr of events in which a non-NAN llh was computed
Private Members
-
I3EventLogLikelihoodBasePtr llhService_¶
log-likelihood service
-
I3ParametrizationBasePtr paramService_¶
-
int nPar_¶
fit name
-
unsigned int nEvent_¶
pretent there were nPar fitting parameters
-
unsigned int nFitOK_¶
counts physics frames
-
unsigned int nFitBad_¶
counts I3Particle with OK fitstatus
-
unsigned int nSuccess_¶
counts I3Particle that are absent or with non-OK fitstatus
-
class I3SimpleFitter : public I3ConditionalModule¶
- #include <I3SimpleFitter.h>
Gulliver-based module to perform simple generic log-likelihood reconstructions.
This module obtains seeds from a seed service (specified with the “SeedService” option), then performs reconstruction using a minimizer service (specified with the “Minimizer” option), a loglikelihood service (specified with the “LogLikelihood” option) to find a best fit, with the fittable variables determined by a parametrization service (specified with the “Parametrization” option).
This module is very generic, it does not make any assumptions about what kind of event you are trying to reconstruct. The actual physics is taken of by the Gulliver services.
For more general information about Gulliver services, read the documentation of of the
gulliver
project, which defines the interfaces for gulliver services. Thelilliput
project contains a collection of actual implementations.The result of the fit consists of two or three objects. First: an I3Particle object, stored under the module name. Second: if the kind of events to be reconstructed required “nonstandard” variables, organized in some class deriving from I3FrameObject, then an object of that class stored under the name composed as the module name plus the “Params” suffix. Third: the details of the log-likelihood reconstruction, such as the final (proper and reduced) likelihood, the number of times the likelihood function was evaluated and the number of degrees of freedom; these are stored as a I3LogLikelihoodFitParams object, under the name composed as the module name plus the “FitParams” suffix.
Public Functions
-
inline ~I3SimpleFitter()¶
-
void Configure()¶
-
void Geometry(I3FramePtr frame)¶
-
void Physics(I3FramePtr frame)¶
-
void Finish()¶
Private Types
Private Functions
-
I3SimpleFitter()¶
-
I3SimpleFitter(const I3SimpleFitter &source)¶
-
I3SimpleFitter &operator=(const I3SimpleFitter &source)¶
-
I3LogLikelihoodFitPtr Fit(I3FramePtr frame, const I3EventHypothesis &seed)¶
-
I3LogLikelihoodFitPtr Fit(I3FramePtr frame, unsigned int nseeds, I3VectorI3ParticlePtr allFits, I3LogLikelihoodFitParamsVectPtr params, std::vector<I3VectorDoublePtr> &traces)¶
- SET_LOGGER ("I3SimpleFitter")
Private Members
-
I3GulliverPtr fitterCore_¶
The core Gulliver object for basic tracks.
-
I3SeedServiceBasePtr seedService_¶
Seed preparation service.
-
I3EventLogLikelihoodBasePtr likelihood_¶
Event loglikelihood calcuation service.
-
I3ParametrizationBasePtr parametrization_¶
Track/shower/anything parametrization service.
-
I3MinimizerBasePtr minimizer_¶
Minimizer service.
-
StoragePolicyType storagePolicy_¶
Option to store single or multiple results.
-
TraceModeType traceMode_¶
Option to store fit tracing information (for debugging)
-
unsigned int eventNr_¶
-
unsigned int nSeeds_¶
-
unsigned int nSuccessFits_¶
-
unsigned int nSuccessEvents_¶
-
inline ~I3SimpleFitter()¶
-
namespace std
STL namespace.
- file I3IterativeFitter.cxx
- #include <algorithm>#include <cassert>#include <cmath>#include <string>#include <utility>#include <vector>#include <boost/make_shared.hpp>#include <gsl/gsl_qrng.h>#include “gulliver-modules/I3IterativeFitter.h”#include “dataclasses/I3Direction.h”#include “dataclasses/geometry/I3Geometry.h”#include “dataclasses/physics/I3Particle.h”#include “icetray/I3ConditionalModule.h”#include “icetray/I3Context.h”#include “icetray/I3Frame.h”#include “icetray/I3Logging.h”#include “icetray/I3Units.h”#include “gulliver/I3EventHypothesis.h”#include “gulliver/I3Gulliver.h”#include “gulliver/I3LogLikelihoodFit.h”#include “gulliver/utilities/ordinal.h”#include “phys-services/I3RandomService.h”
copyright (C) 2005 the icecube collaboration $Id$
- Version
$Revision$
- Date
$Date$
- Author
boersma
Functions
-
I3_MODULE(I3IterativeFitter)¶
Variables
-
static const char *minimizer_optionname = "Minimizer"¶
-
static const char *param_optionname = "Parametrization"¶
-
static const char *llh_optionname = "LogLikelihood"¶
-
static const char *seedservice_optionname = "SeedService"¶
-
static const char *tweakservice_optionname = "IterationTweakService"¶
-
static const char *randomservice_optionname = "RandomService"¶
-
static const char *coszenithrange_optionname = "CosZenithRange"¶
-
static const char *niter_optionname = "NIterations"¶
- file I3IterativeFitter.h
- #include <string>#include <utility>#include <vector>#include “dataclasses/physics/I3Particle.h”#include “gulliver/I3EventLogLikelihoodBase.h”#include “gulliver/I3Gulliver.h”#include “gulliver/I3LogLikelihoodFit.h”#include “gulliver/I3MinimizerBase.h”#include “gulliver/I3ParametrizationBase.h”#include “gulliver/I3SeedServiceBase.h”#include “icetray/I3ConditionalModule.h”#include “icetray/IcetrayFwd.h”#include “phys-services/I3RandomService.h”
copyright (C) 2005 the icecube collaboration $Id$
- Version
$Revision$
- Date
$Date$
- Author
boersma
- file I3LogLikelihoodCalculator.cxx
- #include “icetray/IcetrayFwd.h”#include <string>#include “gulliver-modules/I3LogLikelihoodCalculator.h”#include “icetray/I3Context.h”#include “icetray/I3Frame.h”#include “dataclasses/geometry/I3Geometry.h”#include “gulliver/I3EventLogLikelihoodBase.h”#include “gulliver/I3EventHypothesis.h”#include “gulliver/utilities/ordinal.h”
copyright (C) 2005 the icecube collaboration $Id$
- Version
$Revision$
- Date
$Date$
- Author
boersma
Functions
-
I3_MODULE(I3LogLikelihoodCalculator)¶
- file I3LogLikelihoodCalculator.h
- #include “gulliver/I3Gulliver.h”#include “gulliver/I3SeedServiceBase.h”#include “icetray/I3ConditionalModule.h”#include “icetray/IcetrayFwd.h”#include “dataclasses/physics/I3Particle.h”
copyright (C) 2007 the icecube collaboration $Id$
- Version
$Revision$
- Date
$Date$
- Author
boersma
- file I3SimpleFitter.cxx
- #include <algorithm>#include <cassert>#include <iomanip>#include <sstream>#include <string>#include <vector>#include <boost/make_shared.hpp>#include “gulliver-modules/I3SimpleFitter.h”#include “dataclasses/I3Vector.h”#include “dataclasses/geometry/I3Geometry.h”#include “dataclasses/physics/I3Particle.h”#include “gulliver/I3EventHypothesis.h”#include “gulliver/I3Gulliver.h”#include “gulliver/I3LogLikelihoodFit.h”#include “gulliver/utilities/ordinal.h”#include “icetray/I3ConditionalModule.h”#include “icetray/I3Context.h”#include “icetray/I3Frame.h”#include “icetray/I3Logging.h”#include “icetray/I3Units.h”
copyright (C) 2005 the icecube collaboration $Id$
- Version
$Revision$
- Date
$Date$
- Author
boersma
Functions
-
I3_MODULE(I3SimpleFitter)¶
- file I3SimpleFitter.h
- #include <string>#include <vector>#include “dataclasses/I3Vector.h”#include “dataclasses/physics/I3Particle.h”#include “gulliver/I3Gulliver.h”#include “gulliver/I3SeedServiceBase.h”#include “gulliver/I3EventLogLikelihoodBase.h”#include “gulliver/I3ParametrizationBase.h”#include “gulliver/I3MinimizerBase.h”#include “gulliver/I3LogLikelihoodFit.h”#include “icetray/I3ConditionalModule.h”#include “icetray/IcetrayFwd.h”
copyright (C) 2005 the icecube collaboration $Id$
- Version
$Revision$
- Date
$Date$
- Author
boersma
- page todo
- Class I3LogLikelihoodCalculator
Extend this module such that you can also use it with nonstandard event hypotheses (which need an additional object besides the I3Particle in order to store the relevant physics variables; and/or need a seed service to tweak and such).
- dir gulliver-modules
- dir gulliver-modules
- dir gulliver-modules
- dir icetray
- dir private
- dir public