paraboloid C++ API Reference¶
-
class GridStar¶
- #include <ParaboloidImpl.h>
Create ellipses in a global spherical coordinate system centered around a given point.
- Todo:
finished
See also
Public Functions
-
inline GridStar()¶
-
inline ~GridStar()¶
-
int Init(int n_Pts_on_Circle, int n_Number_Circles, double n_Azi_Reach, double n_Zen_Reach, I3Direction &track)¶
-
int GridNext(double &Azi_Local, double &Zen_Local, double &Azi_Global, double &Zen_Global)¶
Public Members
-
RotateLocalNet localnet¶
-
class I3ParaboloidFitParams : public I3LogLikelihoodFitParams¶
- #include <I3ParaboloidFitParams.h>
The I3LogLikelihoodFitParams extended with paraboloid-specific information.
copyright (C) 2004 the icecube collaboration $Id$
- Version
$Revision$
- Date
$Date$
- Author
David Boersma boersma@icecube.wisc.edu
Public Types
-
enum ParaboloidFitStatus¶
ParaboloidFitStatus: Initialized to PBF_UNDEFINED = -1000 The general scheme is:
Negative values occur when paraboloid cannot successfully complete the fit.
Positive values occur after additional checks, when the fit is deemed “not good”. (Somewhat subjective?)
Identically Zero means everything was okay = PBF_SUCCESS
Values:
-
enumerator PBF_UNDEFINED¶
undefined (initvalue before fit)
-
enumerator PBF_NO_SEED¶
missing input track
-
enumerator PBF_INCOMPLETE_GRID¶
on too many gridpoints vertex refit failed
-
enumerator PBF_FAILED_PARABOLOID_FIT¶
no parabola could be fitted to grid
-
enumerator PBF_SINGULAR_CURVATURE_MATRIX¶
curvature matrix singular
-
enumerator PBF_SUCCESS¶
successful fit
-
enumerator PBF_NON_POSITIVE_ERRS¶
both errors are negative (llh landscape looks like *maximum instead of minimum)
-
enumerator PBF_NON_POSITIVE_ERR_1¶
only error2>0 (saddle point)
-
enumerator PBF_NON_POSITIVE_ERR_2¶
only error1>0 (saddle point)
-
enumerator PBF_TOO_SMALL_ERRS¶
error is infinitesimally small, probably a numerical problem happened
Public Functions
-
inline I3ParaboloidFitParams(double l = NAN, double r = NAN, int dof = -1, int mini = -1)¶
default constructor
-
inline virtual ~I3ParaboloidFitParams()¶
destructor
-
void Reset()¶
put all datamembers to default values
Public Members
-
ParaboloidFitStatus pbfStatus_¶
Paraboloid fit status: initialized to PBF_UNDEFINED = -1000
See also
-
double pbfZen_¶
theta-angle of paraboloid fit (global coordinates)
-
double pbfAzi_¶
phi-angle of paraboloid fit (global coordinates)
-
double pbfErr1_¶
parameter 1 of error ellipse (using diagonalized hesse matrix)
-
double pbfErr2_¶
parameter 2 of error ellipse (using diagonalized hesse matrix)
-
double pbfRotAng_¶
rotation angle of error ellipse (using diagonalized hesse matri
-
double pbfCenterLlh_¶
llh-value of grid center
-
double pbfLlh_¶
llh-value of paraboloid fit
-
double pbfZenOff_¶
theta offset of paraboloid fit (local coordinates)
-
double pbfAziOff_¶
phi-offset of paraboloid fit (local coordinates)
-
double pbfCurv11_¶
(theta) curvature-value of paraboloid fit
-
double pbfCurv12_¶
(1/cov) curvature-value of paraboloid fit
-
double pbfCurv22_¶
(phi) curvature-value of paraboloid fit
-
double pbfChi2_¶
chi^2 of parabola fit
-
double pbfDetCurvM_¶
determinant of curvature matrix of parabola fit
-
double pbfSigmaZen_¶
theta error of fit using inverted curvature matrix
-
double pbfSigmaAzi_¶
phi error of fit using inverted curvature matrix
-
double pbfCovar_¶
covariance from inverted curvature matrix
-
double pbfTrOffZen_¶
true zenith (of highest energy muon) in local paraboloid coordinate system
-
double pbfTrOffAzi_¶
true azimuth (of highest energy muon) in local paraboloid coordinate system
-
double pbfTrZen_¶
true zenith (of highest energy muon)
-
double pbfTrAzi_¶
true azimuth (of highest energy muon)
Friends
- friend class icecube::serialization::access
-
class I3ParaboloidFitParamsConverter : public I3ConverterImplementation<I3ParaboloidFitParams>¶
- #include <I3ParaboloidFitParamsConverter.h>
copyright (C) 2010 The Icecube Collaboration
$Id$
- Version
$Revision$
- Date
$LastChangedDate$
- Author
Fabian Kislat fabian.kislat@desy.de $LastChangedBy$
Private Functions
-
I3TableRowDescriptionPtr CreateDescription(const I3ParaboloidFitParams &p)¶
copyright (C) 2010 The Icecube Collaboration
$Id$
- Version
$Revision$
- Date
$LastChangedDate$
- Author
Fabian Kislat fabian.kislat@desy.de, last changed by $LastChangedBy$
-
size_t FillRows(const I3ParaboloidFitParams &p, I3TableRowPtr rows)¶
-
class I3ParaboloidFitter : public I3ConditionalModule¶
- #include <I3ParaboloidFitter.h>
Evaluate the likelihood on a zenith-azimuth grid near a given track, fit a paraboloid.
The direction given input track is transformed such that it’s on the equator of the direction sphere and the zenith and azimuth are almost Cartesian. For each direction in a grid of directions close to the input directions, the vertex of the track is refit (reoptimized), preferrably using the same minimizer and likelihood function as used to obtain the input fit. The vertex-optimized likelihood values of the grid points define a “2D likelihood landscape” that can be approximated with a paraboloid surface. The minimum of this paraboloid can be used as a (slightly) improved fit, while the geometrical properties of this paraboloid can be used to define quality cuts. For example “error ellips” can defined as the curve in zenith-azimuth space on which the (reoptimized) log-likelihood is 0.5 higher than that of the paraboloid minimum.
Paraboloid is a Gulliver module, so it uses likelihood functions, minimizers, parametrizations and seed services. Since we do not expect/want the vertex at the grid points to be dramatically different from the vertex of the input track, it’s important to configure the parametrization with a small stepsize (default 5m).
(A seed service can do seed tweaking before a fit, e.g. setting the vertex time such that the time residuals of all hits are positive. It can do that both with an external first guess event hypotheses and with an event hypothesis generated by the calling code, the module.)
You want minimal or no tweaking of the input fit (that fit should already be very good), but for the grid points tweaking of the vertex time is useful. For this reason paraboloid can be configured to use the two different seed services: one to provide the input track, and one to doctor the grid seeds.
When using the results of Paraboloid is always very important to use the paraboloid status. In case of a (partial) failure the track variables are not necessarily set to NAN. If you only want to use the results with status 0 (success). The non-zero value of the paraboloid status tells you exactly what went wrong; this can help you to improve the configuration of the fitter modules or to diagnose problems.
The original ideas and implementation for the paraboloid fitter were by Till Neunhoeffer.
See also
https://dx.doi.org/10.25358/openscience-2482 (Till’s PhD thesis, in German)
See also
Astroparticle Physics 25 (2006) 220-225, April 2006: Estimating the Angular Resolution of Tracks in Neutrino Telescopes Based on a Likelihood Analysis
Public Functions
-
I3ParaboloidFitter(const I3Context &ctx)¶
constructor for I3Tray::AddModule (defines configuration parameters)
-
~I3ParaboloidFitter()¶
destructor
-
void Configure()¶
configure (retrieves and checks configuration parameters, initializes the module)
-
void Geometry(I3FramePtr frame)¶
get the geoemtry from the frame and pass it to the LLH
-
void Physics(I3FramePtr frame)¶
do the paraboloid fit for a particular physics event
Private Functions
-
I3ParaboloidFitter()¶
-
I3ParaboloidFitter(const I3ParaboloidFitter &source)¶
inhibited default constructor
-
I3ParaboloidFitter &operator=(const I3ParaboloidFitter &source)¶
inhibited copy constructor
-
int GetParaboloidDatapoints(const I3Frame &f, const I3EventHypothesis &hypothesis)¶
generate grid of directions
-
int GetErrorsFromCurvature(I3ParaboloidFitParams ¶bParam)¶
store direction error estimate
-
int StoreResults(I3ParaboloidFitParams ¶bParam)¶
store other results
-
int GetTruth(I3ParaboloidFitParams ¶bParam, I3FramePtr f)¶
compare with MC track (if given)
-
int GetLogLHFitParams(const I3Frame &f, I3LogLikelihoodFitPtr fitptr)¶
define own fit with the direction from the minimum of the parabola, and a refitted vertex
-
bool IsInDepthRange(const I3Particle &p, double zmin, double zmax)¶
helper method to check if a particle traverses a given depth range
-
bool NDFOK(const I3Frame &f, I3EventHypothesis &h)¶
helper method to check NDF
- SET_LOGGER ("I3ParaboloidFitter")
macro, setting the “MyClass” name for log4cplus.conf
Private Members
-
I3GulliverPtr fitterCore_¶
inhibited assignment operator
the core Gulliver object for basic tracks (a pointer, so that we can postpone object creation to Configure())
-
I3SeedServiceBasePtr seedService_¶
seed service for input track
-
I3SeedServiceBasePtr gridSeedService_¶
seed service for tweaking the grid point seed
-
I3EventLogLikelihoodBasePtr eventLLH_¶
event likelihood service (should be same as the one used for the input fit)
-
I3MinimizerBasePtr minimizer_¶
-
int nMaxMissingGridPoints_¶
name of minimizer service
-
double vertexStepSize_¶
max number of grid points for which vertex refit may fail, before we give up
-
double azimuthReach_¶
initial stepsize for vertex refit (not too big)
-
double zenithReach_¶
grid scale for azimuth
-
unsigned int nSteps_¶
grid scale for zenith
-
unsigned int nSamplingPoints_¶
number of grid rings around input direction
-
int ndf_¶
number of grid points on each grid ring
NDF: compute once, use anywhere during physics event
-
unsigned int eventNr_¶
event counter, for log_messages
-
unsigned int zollCounter_¶
-
I3SimpleParametrizationPtr vertexParametrization_¶
Vertex parametrization In contrast to most other likelihood fitters, paraboloid does not fit arbitrary parameters. For each grid point, the vertex position is reoptimized. In a later stage we might want to rethink this. Maybe in some cases also the energy needs to be reoptimized, for example.
-
ParaboloidImpl::GridStar grid_¶
object generating the grid of directions
-
ParaboloidImpl::ParabolaFit paraFit_¶
object containing paraboloid fit details
-
ParaboloidImpl::ParabolaTypePol_2 paraPol_¶
object containing paraboloid fit details
-
ParaboloidImpl::ParabolaTypeStd_2 paraStd_¶
object containing paraboloid fit details
-
double chi2_¶
“chi squared” (likelihood)
-
double seedLlh_¶
likelihood of input fit (after vertex refit)
-
I3ParaboloidFitter(const I3Context &ctx)¶
-
class ParabolaFit¶
- #include <ParaboloidImpl.h>
Fits a two dimensional parabola to a symmetric set of sampling points.
- Todo:
finished
Public Functions
-
inline ParabolaFit(size_t n = Standard_Vector_Size)¶
-
inline ~ParabolaFit()¶
-
inline void Clear()¶
-
inline size_t size() const¶
-
int lin_reg_parabola_2_sym(ParabolaTypePol_2 ¶, double *chi2 = NULL, double *detM = NULL)¶
-
class ParabolaTypePol_2¶
- #include <ParaboloidImpl.h>
Describes a 2 dimensional Parabola using alpha, beta, gamma as parameters.
- Todo:
finished
See also
Public Functions
-
inline ParabolaTypePol_2()¶
-
inline ParabolaTypePol_2(ParabolaTypeStd_2 &std)¶
-
inline ~ParabolaTypePol_2()¶
-
ParabolaTypePol_2 &operator=(ParabolaTypeStd_2 &std)¶
-
double get_y(double x[2])¶
-
class ParabolaTypeStd_2¶
- #include <ParaboloidImpl.h>
Describes a 2 dimensional Parabola using A, b, c as parameters.
- Todo:
finished
See also
Public Functions
-
inline ParabolaTypeStd_2()¶
-
inline ParabolaTypeStd_2(ParabolaTypePol_2 &pol)¶
-
inline ~ParabolaTypeStd_2()¶
-
ParabolaTypeStd_2 &operator=(ParabolaTypePol_2 &pol)¶
-
double get_y(double x[2])¶
-
class RotateLocalNet¶
- #include <ParaboloidImpl.h>
Rotate coordinates in a local coordinate system (theta,phi) to a global coordinate system and vice versa.
- Todo:
finished
See also
-
class SumType¶
Public Functions
-
inline SumType()¶
-
inline ~SumType()¶
-
void BuildManySums2(const std::vector<ParaboloidImpl::XYTupel> &data_points)¶
-
bool IsSymmetric()¶
Private Functions
-
void Clear()¶
-
inline SumType()¶
-
class XYTupel¶
- #include <ParaboloidImpl.h>
Basic data type for ParabolaFit.
- Todo:
finished
See also
-
namespace ParaboloidImpl¶
Functions
-
double determinant(const bnu::matrix<double> &input)¶
Surprisingly, boost::numeric::ublas does not have its own determinant calculation.
- Parameters:
input – [in] a square matrix
-
bool invert_matrix(const bnu::matrix<double> &input, bnu::matrix<double> &inverse)¶
Surprisingly, boost::numeric::ublas does not have its own matrix inverse calculation. But it’s relatively easy to get that from the available ublas functions.
- Parameters:
input – [in] a square matrix
inverse – [out] a matrix with the same shape as the input
- Returns:
true = success false = failure (LU facturization failed)
-
int solve_linear_equation(bnu::matrix<double> input, const bnu::vector<double> &inhomo, bnu::vector<double> &solution)¶
input is passed by value, because we need a modifiable copy
-
int Kramer(const bnu::matrix<double> &tmatrixd, const bnu::vector<double> &inhomo, double &det, bnu::vector<double> &solution)¶
solve a inhomogeneous linear equation using kramer’s method (not very efficient for n>3)
Kramer
- Todo:
finished
- Parameters:
tmatrixd – [in] inhomogenous equation system
inhomo – [in] inhomogenous equation system
det – [out] Determinant of matrix
solution – [out] solution of equation system saved in vector
- Returns:
-1 = something went wrong (e.g. dimensions of matrix/vector not compatible) 0 = determinant of matrix is zero -> no unique solution 1 = success
-
int diagonalize_sym_2_2_matrix(bnu::matrix<double> &tmatrixd, double &eigen1, double &eigen2, double &angle)¶
diagonalize a symmetric 2x2 matrix, calculate the eigenvalues and the angle between the eigen vectors
diagonalize_sym_2_2_matrix
- Todo:
finished
- Parameters:
tmatrixd – [in]
eigen1 – [out] eigenvalue of the matrix
eigen2 – [out] eigenvalues of the matrix
angle – [out] = angle between eigen vectors in radians
- Returns:
-1 = something went wrong (e.g. matrix not 2x2) 0 = determinant of matrix < 0 -> no eigenvalues 1 = success
-
double determinant(const bnu::matrix<double> &input)¶
-
namespace std
STL namespace.
- file I3ParaboloidFitParams.cxx
- #include “icetray/serialization.h”#include “paraboloid/I3ParaboloidFitParams.h”
copyright (C) 2006 the icecube collaboration $Id$
- Version
$Revision$
- Date
$Date$
- Author
David Boersma boersma@icecube.wisc.edu
Functions
-
std::ostream &operator<<(std::ostream &os, const I3ParaboloidFitParams &p)¶
-
I3_SERIALIZABLE(I3ParaboloidFitParams)¶
- file I3ParaboloidFitParams.h
- #include “gulliver/I3LogLikelihoodFitParams.h”
Functions
-
std::ostream &operator<<(std::ostream&, const I3ParaboloidFitParams&)
- I3_CLASS_VERSION (I3ParaboloidFitParams, 2)
-
I3_POINTER_TYPEDEFS(I3ParaboloidFitParams)¶
-
std::ostream &operator<<(std::ostream&, const I3ParaboloidFitParams&)
- file I3ParaboloidFitParamsConverter.cxx
- #include “I3ParaboloidFitParamsConverter.h”
- file I3ParaboloidFitParamsConverter.h
- #include <tableio/I3Converter.h>#include <paraboloid/I3ParaboloidFitParams.h>
- file I3ParaboloidFitter.cxx
- #include “icetray/IcetrayFwd.h”#include <string>#include <sstream>#include <cmath>#include “paraboloid/I3ParaboloidFitter.h”#include “icetray/I3Context.h”#include “icetray/I3Frame.h”#include “dataclasses/physics/I3MCTreeUtils.h”#include “dataclasses/physics/I3EventHeader.h”#include “dataclasses/geometry/I3Geometry.h”#include “gulliver/I3EventLogLikelihoodBase.h”#include “lilliput/parametrization/I3SimpleParametrization.h”#include “gulliver/I3Gulliver.h”#include “gulliver/I3MinimizerBase.h”#include “gulliver/I3EventHypothesis.h”#include “gulliver/utilities/ordinal.h”#include “paraboloid/I3ParaboloidFitParams.h”
copyright (C) 2005 the icecube collaboration $Id$
- Version
$Revision$
- Date
$Date$
- Author
boersma
Functions
-
I3_MODULE(I3ParaboloidFitter)¶
- file I3ParaboloidFitter.h
- #include <string>#include “gulliver/I3Gulliver.h”#include “gulliver/I3SeedServiceBase.h”#include “lilliput/parametrization/I3SimpleParametrization.h”#include “icetray/I3ConditionalModule.h”#include “icetray/IcetrayFwd.h”#include “dataclasses/physics/I3Particle.h”#include “paraboloid/ParaboloidImpl.h”#include “paraboloid/I3ParaboloidFitParams.h”
copyright (C) 2005 the icecube collaboration $Id$
- Version
$Revision$
- Date
$Date$
- Author
David Boersma boersma@icecube.wisc.edu
- Author
Chad Finley cfinley@icecube.wisc.edu
- Author
Jon Dumm jdumm@icecube.wisc.edu
- Author
Timo Griesel timo.griesel@uni-mainz.de
- file ParaboloidImpl.cxx
- #include <cmath>#include <boost/numeric/ublas/matrix.hpp>#include <boost/numeric/ublas/lu.hpp>#include <boost/numeric/ublas/io.hpp>#include <dataclasses/I3Direction.h>#include <dataclasses/I3Quaternion.h>#include <icetray/I3Logging.h>#include “paraboloid/ParaboloidImpl.h”
(C) copyright 2006 the icecube collaboration $Id$
This code was almost literally copied from sieglinde/reconstruction/SLParaboloidAux.(hh|cc). The code in SLParaboloidAux.(hh|cc) was written by Marc Hellwig (with minor edits by DJB), as a thorough rewrite of the original parabola fit by Till Neunhoeffer in siegmund/recoos.
Imported/Copied/Adapted to IceTray/Gulliver by David Boersma Furthere eveloped and maintained by Timo Griesel (and DJB)
- Version
$Revision$
- Date
$Date$
- Author
David Boersma boersma@icecube.wisc.edu
- Author
Timo Griesel timo.griesel@uni-mainz.de
- file ParaboloidImpl.h
- #include <vector>#include <dataclasses/I3Direction.h>#include <dataclasses/I3Quaternion.h>#include <boost/numeric/ublas/matrix.hpp>#include <boost/numeric/ublas/lu.hpp>#include <boost/numeric/ublas/io.hpp>
(C) 2006 the icecube collaboration $Id$
This code was copied from sieglinde/reconstruction/SLParaboloidAux.hh Written by Marc Hellwig Imported to IceTray/Gulliver by David Boersma Further developed and maintained by Timo Griesel, Chad Finley and Jon Dumm
- Version
$Revision$
- Date
$Date$
- Author
David Boersma boersma@icecube.wisc.edu
- Author
Timo Griesel griesel@uni-mainz.de
- Author
Chad Finley cfinley@icecube.wisc.edu
- Author
Jon Dumm jdumm@icecube.wisc.edu
- page todo
- Member ParaboloidImpl::diagonalize_sym_2_2_matrix (bnu::matrix< double > &tmatrixd, double &eigen1, double &eigen2, double &angle)
finished
- Class ParaboloidImpl::GridStar
finished
- Member ParaboloidImpl::Kramer (const bnu::matrix< double > &tmatrixd, const bnu::vector< double > &inhomo, double &det, bnu::vector< double > &solution)
finished
- Class ParaboloidImpl::ParabolaFit
finished
- Class ParaboloidImpl::ParabolaTypePol_2
finished
- Class ParaboloidImpl::ParabolaTypeStd_2
finished
- Class ParaboloidImpl::RotateLocalNet
finished
- Class ParaboloidImpl::XYTupel
finished
- dir icetray
- dir paraboloid
- dir paraboloid
- dir paraboloid
- dir private
- dir public