Python API¶
IceTray segments¶
- icecube.lilliput.segments.I3SinglePandelFitter(tray, name, fitname='', pulses='OfflinePulses', seeds=[], minimizer=None, parametrization=None, domllh='SPE1st', tstype='TFirst', noiserate=1e-08, If=None)
Run a single Pandel fit.
Use standard settings as for IC40-IC79 and bulk ice properties.
Log-likelihood functions:
SPE1st: use time of first pulse in each DOM.
SPEAll: use times of all pulses in each DOM.
SPEqAll: use charge-weighted times of all pulses in each DOM.
MPE: use time of first pulse and total charge in each DOM.
- Parameters:
fitname (str) – Fit name
pulses (str) – Name of pulse map
minimizer (I3MinimizerBase or str, optional) – Minimizer service
parametrization (I3ParametrizationBase or str, optional) – Track parametrization service
seeds (list) – Seed names
domllh ({"SPE1st", "SPEAll", "SPEqAll", "MPE"}, optional) – Likelihood function
tstype (str, optional) – Use
"TNone"
if seed has a reliable vertex time.noiserate (float, optional) – Noise rate
- Returns:
Names of services and fit for re-usage convenience
- Return type:
Warning
Although SPEAll and SPEqAll were thought to be fundamentally correct, they work surprisingly worse than SPE1st.
Notes
The default noise rate of 10Hz does not have any physical basis. It was set at the time when we used only HLC pulses and based on flawed logic (the HLC dark noise rate was estimated to be of order 1Hz and then 10Hz seemed sort of conservative). Studies on low filtering level indicate that higher noise level reduce the rate of mis-reconstructed background, while at high cut level the effect of this setting on well reconstructed neutrino events seems to be minor. New studies and developments are underway (as of November 2011).
- icecube.lilliput.segments.I3IterativePandelFitter(tray, name, fitname='', pulses='OfflinePulses', n_iterations=4, minimizer=None, parametrization=None, seeds=[], domllh='SPE1st', tstype='TFirst', noiserate=1e-08, If=None)
Run an iterative Pandel Fit.
Use standard settings as for IC40-IC79 and bulk ice properties.
- Parameters:
fitname (str) – Fit name
pulses (str) – Name of pulse map
n_iterations (int, optional) – Number of iterations
minimizer (I3MinimizerBase or str, optional) – Minimizer service
parametrization (I3ParametrizationBase or str, optional) – Track parametrization service
seeds (list) – Seed names
domllh ({"SPE1st", "SPEAll", "SPEqAll", "MPE"}, optional) – Likelihood function
tstype (str, optional) – Use
"TNone"
if seed has a reliable vertex time.noiserate (float, optional) – Noise rate
- Returns:
Names of services and fit for re-usage convenience
- Return type:
- icecube.lilliput.segments.I3ParaboloidPandelFitter(tray, name, fitname='', pulses='OfflinePulses', minimizer=None, parametrization=None, inputtrack=None, domllh='SPE1st', input_tstype='TFirst', grid_tstype='TFirst', noiserate=1e-08, If=None)
Run a Paraboloid fit with a Pandel likelihood.
Use standard settings as used for OnlineL2 and bulk ice properties This tray segment only works for gulliver Pandel track fits.
- Parameters:
fitname (str) – Fit name
pulses (str) – Name of pulse map as for input track
minimizer (I3MinimizerBase or str, optional) – Minimizer service
parametrization (I3ParametrizationBase or str, optional) – Track parametrization service
inputtrack (str) – Name of input track
domllh ({"SPE1st", "SPEAll", "SPEqAll", "MPE"}, optional) – Likelihood as for input track
input_tstype (str, optional) – Vertex time strategy as for input track
grid_tstype (str, optional) – Use
"TFirst"
or"TChargeFraction"
for grid.noiserate (float) – Noise rate as for input track
- Returns:
Names of services and fit for re-usage convenience
- Return type:
Notes
The default vertex time strategy of
"TFirst"
for the input track was chosen because it is the standard for online level 2 and level 3. Supposing the input track has already an excellent timing,"TNone"
should be the logical and optimal choice. According to some incomplete testing there should not be any dramatic differences between the two strategies.
- icecube.lilliput.segments.add_minuit_simplex_minimizer_service(tray, minimizer=None)
Add Minuit (Simplex) minimizer service to tray.
We need only one Minuit minimizer service in the tray, this can be used by all gulliver fitters that need it.
- icecube.lilliput.segments.add_simple_track_parametrization_service(tray, parametrization=None)
Add track parametrization service to tray.
We need only one track parametrization service in the tray, this can be used by all gulliver fitters that need it.
- icecube.lilliput.segments.add_pandel_likelihood_service(tray, pulses, domllh='SPE1st', noiserate=1e-08)
Add Pandel service to tray.
For every desired Pandel configuration we need at most one instance of the Pandel service in the tray. This service can be used by all gulliver fitters that need it.
With this function you can do some configuration, but not all. If you want to do more fancy configuration, then you’ll need to add the service yourself instead of using this convenience function.
- icecube.lilliput.segments.add_seed_service(tray, pulses, seeds=[], tstype='TFirst')
Add seed service to tray.
For every desired seeding configuration we need at most one instance of the basic seed service in the tray. This service can be used by all gulliver modules that need it.
Minimizer¶
- class icecube.lilliput.scipymin.SciPyMinimizer(name, method='L-BFGS-B', tolerance=10000000.0, max_iterations=100, options=None)
Bases:
I3Minimizer
Multi-dimensional minimizers in SciPy
- Parameters:
- Minimize((I3Minimizer)arg1, (I3GulliverBase)arg2, (I3FitParameterInitSpecsSeries)arg3) I3MinimizerResult :
- C++ signature :
I3MinimizerResult Minimize(I3MinimizerBase {lvalue},I3GulliverBase {lvalue},std::__1::vector<I3FitParameterInitSpecs, std::__1::allocator<I3FitParameterInitSpecs>>)
Minimize( (I3Minimizer)arg1, (I3GulliverBase)arg2, (I3FitParameterInitSpecsSeries)arg3) -> None :
- C++ signature :
void Minimize(I3MinimizerWrapper {lvalue},I3GulliverBase {lvalue},std::__1::vector<I3FitParameterInitSpecs, std::__1::allocator<I3FitParameterInitSpecs>>)
- UsesGradient((I3Minimizer)arg1) bool :
- C++ signature :
bool UsesGradient(I3MinimizerBase {lvalue})
- class icecube.lilliput.i3minuit.IMinuitMinimizer(name='iminuit', Tolerance=0.1, MaxIterations=10000, MinuitPrintLevel=0, MinuitStrategy=2, Algorithm='SIMPLEX', FlatnessCheck=False, WithGradients=False, CheckGradient=False, IgnoreEDM=False, MinuitPrecision=None)
Bases:
I3Minimizer
Wrapper for the iminuit minimization algorithm for use with the Gulliver suite.
It is intended to be a drop-in replacement for the ROOT-based minimizers
I3GulliverMinuitFactory
andI3GulliverMinuit2Factory
found in the Gulliver suite. iminuit is a stand alone Python package which contains Minuit2. It is relativly well maintained and unlike the version which comes with ROOT, it does not print anything to stdout when you set the print level to zero.If you do not have iminuit, it can be insalled very easily using the following command:
pip install iminuit
- Parameters:
name (string, optional) – String for the Gulliver module to identify the minimizer
tolerance (float, optional) – Tolerance for finding the minimum
MaxIterations (int, optional) – Maximum number of iterations to perform before giving up
MinuitPrintLevel (int, optional) – Set print level: 0 = quiet, 1 = normal, 2 = paranoid, 3 = really paranoid
MinuitStrategy (int,optional) – 0 = fast, 1 = default, 2 = slow but valid
Algorithm (string,optional) – Currently SIMPLEX, MIGRAD
FlatnessCheck (noop, not implemented) –
WithGradients (bool, optional) – Whether to use analytic gradients
CheckGradients (noop, iminuit does not check the first gradient numerically) –
IgnoreEDM (bool, optional) – Ignore EDM when checking for convergence. If True requires that fmin be both valid and accurate.
See also
Further reading: Wikipedia, iminuit documentation, Minuit CERN page, Minuit paper.
- Minimize((I3Minimizer)arg1, (I3GulliverBase)arg2, (I3FitParameterInitSpecsSeries)arg3) I3MinimizerResult :
- C++ signature :
I3MinimizerResult Minimize(I3MinimizerBase {lvalue},I3GulliverBase {lvalue},std::__1::vector<I3FitParameterInitSpecs, std::__1::allocator<I3FitParameterInitSpecs>>)
Minimize( (I3Minimizer)arg1, (I3GulliverBase)arg2, (I3FitParameterInitSpecsSeries)arg3) -> None :
- C++ signature :
void Minimize(I3MinimizerWrapper {lvalue},I3GulliverBase {lvalue},std::__1::vector<I3FitParameterInitSpecs, std::__1::allocator<I3FitParameterInitSpecs>>)
- UsesGradient((I3Minimizer)arg1) bool :
- C++ signature :
bool UsesGradient(I3MinimizerBase {lvalue})
Example¶
Example how to use the minimizers in an IceTray script:
# Add the standard gulliver service factories.
tray.AddService("I3SimpleParametrizationFactory", "param", ...)
tray.AddService("I3RecoLLHFactory", "pandel", ...)
tray.AddService("I3BasicSeedServiceFactory", "seed", ...)
# Add the python class directly into the tray's context
# (no need for a service factory).
tray.context["iminuit"] = lilliput.IMinuitMinimizer(tolerance=0.001)
# Now call the gulliver module in the normal way.
tray.AddModule("I3SimpleFitter", "singlefit",
SeedService="seed",
Parametrization="param",
LogLikelihood="pandel",
Minimizer="iminuit")