spline-reco C++ API Reference¶
-
class I3SplineRecoLikelihood : public I3EventLogLikelihoodBase¶
- #include <I3SplineRecoLikelihood.h>
Public Types
Public Functions
-
inline I3SplineRecoLikelihood()¶
-
inline virtual ~I3SplineRecoLikelihood()¶
-
void SetGeometry(const I3Geometry &f)¶
-
double GetLogLikelihood(const I3EventHypothesis&)¶
-
double GetLogLikelihood(const I3EventHypothesis&, I3EventHypothesis*, bool solve_for_energies, double weight = 1)¶
-
double GetPerDomLogLikelihood(const I3RecoPulseSeries &pulses, const OMKey &omkey, const I3Particle &p, PhotonicsSource&)¶
-
unsigned int GetMultiplicity()¶
Public Members
-
I3PhotonicsServicePtr ps¶
basic spline, probably bare muon
-
I3PhotonicsServicePtr stochastics_ps¶
stochastics spline
-
I3PhotonicsServicePtr random_noise_ps¶
noise spline = bare spline with jitter = 1000ns (I really should change that, and use a parametrisation - Kai Schatto)
-
double floorWeight¶
-
bool modelStochastics¶
-
int confidence_level¶
Kolmogorov Smirnov options.
-
int chargeCalcStep¶
-
bool EnergyDependentMPE¶
-
bool EnergyDependentJitter¶
-
double preJitter¶
jitter options
-
double postJitter¶
-
double noiseRate¶
-
bool chargeWeight¶
Private Functions
-
double GetMPEModifier()¶
-
double GetEnergyDependentPostJitter()¶
-
unsigned int GetMultiplicityNCh(I3RecoPulseSeriesMapConstPtr pulses)¶
-
unsigned int GetMultiplicityNPulses(I3RecoPulseSeriesMapConstPtr pulses)¶
-
unsigned int GetMultiplicityNPulsesChargeWeight(I3RecoPulseSeriesMapConstPtr pulses)¶
-
double CalculateTimeWindow(I3RecoPulseSeriesMapConstPtr pulses)¶
-
double BuildWeightedDensityFunction(double bare_pes, double stoch_pes, double bareDF, double rDF, double floorDF, double afterDF, double stochDF = 0.)¶
-
double BuildWeightedPdf(double tres, double bare_pes, double stoch_pes)¶
Those two functions construct pdf/cdf by merging different pdfs/cdfs depending on options like “modelStochastics” and “noiseModel”. Bare muon pes and stochastic spline pes are needed for weighting. the result is the normalized pdf/cdf value at tres
-
double BuildWeightedCdf(double tres, double bare_pes, double stoch_pes, double hitTime)¶
-
double GetCdf(double tres, I3PhotonicsServicePtr &spline)¶
returns cdf/pdf from given spline at tres
-
double GetPdf(double tres, I3PhotonicsServicePtr &spline)¶
-
double ConvolvedCdf(double tres, I3PhotonicsServicePtr &spline, double jitter)¶
wrapper to return approximative fast gauss convoluted cdf/pdf (makes no sense for cdf, though): Prejitter.
-
double ConvolvedPdf(double tres, I3PhotonicsServicePtr &spline, double jitter)¶
-
double BuildMPE(double npe, double tres, int photon_counter, double mean_pes, double stochMean_pes, double noiseprob, double hitTime)¶
returns MPE/SPE. Use photon_counter (k) in loop for MPEAll
-
double BuildSPE(double tres, double mean_pes, double stochMean_pes, double noiseprob, double charge)¶
-
double ConvolvedMPE(double npe, double tres, int photon_counter, double mean_pes, double stochMean_pes, double noiseprob, double jitter, double hitTime)¶
wrapper to return approximative fast gauss convoluted MPE/SPE: Postjitter
-
double ConvolvedSPE(double tres, double mean_pes, double stochMean_pes, double noiseprob, double charge, double jitter)¶
-
void GaussianIIR1D(float *data, long length, float sigma, int numsteps)¶
fast approximative gaussian convolution
-
void GetNoiseWeights(double &physicsWeight, double &randomWeight, double &afterWeight, double &preLateWeight, double expected_pes)¶
fitted noise weights for noise modelling. See https://butler.physik.uni-mainz.de/icecube/thesis/bachelor_joschua.pdf
-
double KS0Limit(double npulses, int level)¶
Kolmogorov-Smirnov limits.
-
double GetApprovedCharge(const I3RecoPulseSeries &pulses, const OMKey &omkey, const I3Particle p, const double geotime, const double mean_pes, const double stochMean_pes)¶
returns only the charge of a dom that fits to pdf by performing Kolmogorov-Smirnov-Test
-
long double MPECoeff(int n, int k)¶
calculates n!/((k-1)!*(n-k)!) avoiding rounding errors due to huge factorials. Could be faster, too, because it’s cancelling down before calculation. Loops and picks one number from denominator and one number from numerator at a time and multiplies the quotient to the result.
- SET_LOGGER ("I3SplineRecoLikelihood")
Private Members
-
I3RecoPulseSeriesMapConstPtr pulses_¶
-
I3GeometryConstPtr geometry_¶
-
std::map<OMKey, double> domCharge_¶
stores the dom charge calculated by Kolmogorov-Smirnov Test or in case of KS-test switched off just the dom charge avoids recalculation every iteration
-
const double jitterWidth_¶
reasonable gauss conv settings
-
const double posResLimit_¶
calculation width of jitter in multiples of sigma(preJitter/postJitter)
-
const int convDataLength_¶
convolving only if tres < posResLimit;
-
const int halfSupport_¶
must be odd. Number of sample points
-
const int numSteps_¶
-
double enEstimate_¶
accuracy of approximation
-
double timeWindow_¶
-
double startTime_¶
-
double endTime_¶
-
int chargeCalcCount_¶
-
inline I3SplineRecoLikelihood()¶
-
class I3SplineRecoLikelihoodFactory : public I3ServiceFactory¶
- #include <I3SplineRecoLikelihoodFactory.h>
Public Functions
-
void Configure()¶
Private Members
-
I3SplineRecoLikelihoodPtr llh_¶
-
void Configure()¶
-
namespace std
STL namespace.
- file gaussianiir1d.cxx
- #include <math.h>#include “spline-reco/I3SplineRecoLikelihood.h”
Fast 1D Gaussian convolution IIR approximation.
Fast 1D Gaussian convolution
$Id$
Copyright (c) 2011, Pascal Getreuer All rights reserved.
(gaussianiir1d.c)
This program is free software: you can redistribute it and/or modify it under, at your option, the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version, or the terms of the simplified BSD license.
- Date
- Rcs
- Author
Pascal Getreuer getreuer@gmail.com
You should have received a copy of these licenses along with this program. If not, see http://www.gnu.org/licenses/ and http://www.opensource.org/licenses/bsd-license.html.
Implements the fast Gaussian convolution algorithm of Alvarez and Mazorra, where the Gaussian is approximated by a cascade of first-order infinite impulsive response (IIR) filters. Boundaries are handled with half-sample symmetric extension.
Gaussian convolution is approached as approximating the heat equation and each timestep is performed with an efficient recursive computation. Using more steps yields a more accurate approximation of the Gaussian. A reasonable default value for
numsteps
is 4.Reference: Alvarez, Mazorra, “Signal and Image Restoration using Shock Filters and
Anisotropic Diffusion,” SIAM J. on Numerical Analysis, vol. 31, no. 2, pp. 590-605, 1994.
- param data:
the data to be convolved, modified in-place
- param length:
number of elements
- param sigma:
the standard deviation of the Gaussian in pixels
- param numsteps:
number of timesteps, more steps implies better accuracy
- file I3SplineRecoLikelihood.cxx
- #include “spline-reco/I3SplineRecoLikelihood.h”#include <boost/shared_ptr.hpp>#include <boost/python.hpp>#include “dataclasses/physics/I3MCTreePhysicsLibrary.hh”#include “gulliver/I3Gulliver.h”
This file contains the implementation of the I3SplineRecoLikelihood class, which is a modification of Jake Feintzeig’s splineReco module to perform basic spline reconstructions.
Implementation of I3SplineRecoLikelihood
$Id$
(c) 2012 The IceCube Collaboration http://www.icecube.wisc.edu
This file is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation; either version 3 of the License, or (at your option) any later version.
- Version
- Rcs
93523
- Date
- Rcs
2012-09-26 11:16:32 +0200 (Wed, 26 Sep 2012)
- Author
Kai Schatto KaiSchatto@gmx.de See https://wiki.icecube.wisc.edu/index.php/Improved_likelihood_reconstruction for more information.
This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.
You should have received a copy of the GNU General Public License along with this program. If not, see http://www.gnu.org/licenses/
- file I3SplineRecoLikelihood.h
- #include <string>#include <vector>#include “gulliver/I3EventLogLikelihoodBase.h”#include “gulliver/I3EventHypothesis.h”#include “dataclasses/geometry/I3Geometry.h”#include “dataclasses/physics/I3Particle.h”#include “dataclasses/physics/I3RecoPulse.h”#include “photonics-service/I3PhotonicsService.h”#include “icetray/OMKey.h”
This file contains the definition of the I3SplineRecoLikelihood class, which is a modification of Jake Feintzeig’s splineReco module to perform basic spline reconstructions.
Definition of I3SplineRecoLikelihood
$Id$
Copyright (C) 2012 The IceCube Collaboration http://www.icecube.wisc.edu
This file is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation; either version 3 of the License, or (at your option) any later version.
- Version
$Revision$
- Date
$Date$
- Author
Kai Schatto KaiSchatto@gmx.de See https://wiki.icecube.wisc.edu/index.php/Improved_likelihood_reconstruction for more information.
This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.
You should have received a copy of the GNU General Public License along with this program. If not, see http://www.gnu.org/licenses/
Defines
-
I3SPLINERECOLIKELIHOOD_LlhChoice¶
Functions
-
I3_POINTER_TYPEDEFS(I3SplineRecoLikelihood)¶
- file I3SplineRecoLikelihoodFactory.cxx
- #include “spline-reco/I3SplineRecoLikelihoodFactory.h”#include “icetray/I3Units.h”#include “photonics-service/I3PhotonicsService.h”
This file contains the implementation of the I3SplineRecoLikelihoodFactory service factory, which is a modification of Jake Feintzeig’s splineReco module to perform basic spline reconstructions.
Implementation of I3SplineRecoLikelihoodFactory
$Id$
(c) 2012 The IceCube Collaboration http://www.icecube.wisc.edu
This class, I3SplineRecoLikelihoodFactory, is the one with the interface to IceTray. Use this class as service in your IceTray scripts, with the name of an I3PhotoSplineServiceFactory in the parameter `PhotonicsService’. Give the name of this service to your fitter, e.g. SimpleFitter, in parameter `LogLikelihood’. Look at the example scripts in the resources directory.
- Version
- Rcs
93523
- Date
- Rcs
2012-09-26 11:16:32 +0200 (Wed, 26 Sep 2012)
- Author
Kai Schatto KaiSchatto@gmx.de See https://wiki.icecube.wisc.edu/index.php/Improved_likelihood_reconstruction for more information.
This file is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation; either version 3 of the License, or (at your option) any later version.
This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.
You should have received a copy of the GNU General Public License along with this program. If not, see http://www.gnu.org/licenses/
Functions
-
I3_SERVICE_FACTORY(I3SplineRecoLikelihoodFactory)¶
- file I3SplineRecoLikelihoodFactory.h
- #include <math.h>#include “icetray/I3ServiceFactory.h”#include “spline-reco/I3SplineRecoLikelihood.h”
This file contains the definition of the I3SplineRecoLikelihoodFactory module, which is a modification of Jake Feintzeig’s splineReco module to perform basic spline reconstructions.
Definition of I3SplineRecoLikelihoodFactory
$Id$
Copyright (C) 2012 The IceCube Collaboration http://www.icecube.wisc.edu
This class, I3SplineRecoLikelihoodFactory, is the one with the interface to IceTray. Use this class as service in your IceTray scripts, with the name of an I3PhotoSplineServiceFactory in the parameter `PhotonicsService’. Give the name of this service to your fitter, e.g. SimpleFitter, in parameter `LogLikelihood’. Look at the example scripts in the resources directory.
- Version
$Revision$
- Date
$Date$
- Author
Kai Schatto KaiSchatto@gmx.de See https://wiki.icecube.wisc.edu/index.php/Improved_likelihood_reconstruction for more information.
This file is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation; either version 3 of the License, or (at your option) any later version.
This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.
You should have received a copy of the GNU General Public License along with this program. If not, see http://www.gnu.org/licenses/
- dir icetray
- dir private
- dir public
- dir spline-reco
- dir spline-reco
- dir spline-reco