segmented-spline-reco C++ API Reference¶
-
class I3EquatorParametrization : public I3ServiceBase, public I3ParametrizationBase¶
- #include <I3EquatorParametrization.h>
A parametrization which restricts the direction phase space to a the hemisphere centered around the seed track direction.
Public Functions
-
I3EquatorParametrization(const std::string &name, double ddir, double dxyz, double dLogE, double dt, bool use_rotated_vertex_para, bool use_differential_eloss_para, bool fit_effective_domeff)¶
constructor with full initialization, for unit tests
-
virtual ~I3EquatorParametrization()¶
destructor
-
void Configure() override¶
configure
-
bool InitChainRule(bool wantgradient, bool wanthessian) override¶
-
void ApplyChainRule() override¶
-
void UpdatePhysicsVariables() override¶
this should calculate datamembers of the I3Particle from the values in par
this should calculate datamembers of EmissionHypothesis from the values in par_
-
void UpdateParameters() override¶
This should calculate the values in par_ from datamembers of the I3Particle If relevant it should also update stepsizes.
This should calculate the values in par_ from datamembers of EmissionHypothesis If relevant it should also update stepsizes.
-
inline bool CheckEnergyFit() override¶
- SET_LOGGER ("I3EquatorParametrization")
macro which sets the name to use to configure log4cplus.conf
Private Functions
-
I3EquatorParametrization()¶
inhibit the default constructor
-
I3EquatorParametrization(const I3EquatorParametrization&)¶
inhibit the copy constructor
-
I3EquatorParametrization operator=(const I3EquatorParametrization &rhs)¶
inhibit the assignment operator
-
void InitializeFitSpecs()¶
initialize the FitParameterInitSpecs and parameter vectors
Private Members
-
I3Particle seed_track_¶
-
std::vector<I3Particle> seed_cascades_¶
-
bool var_dir_¶
-
double step_dir_¶
-
bool var_loge_¶
-
double step_loge_¶
-
int var_vertex_pos_¶
-
double step_vertex_pos_¶
-
bool use_rotated_vertex_para_¶
-
bool differential_eloss_para_¶
-
bool fit_effective_domeff_¶
-
bool var_dt_¶
-
double step_dt_¶
-
I3EquatorParametrization(const std::string &name, double ddir, double dxyz, double dLogE, double dt, bool use_rotated_vertex_para, bool use_differential_eloss_para, bool fit_effective_domeff)¶
-
class I3SegmentedRecoSeedService : public I3SeedServiceBase¶
- #include <I3SegmentedRecoSeedService.h>
Base class seed track collection & preparation services.
Basic seed services collects first guess tracks and can do a simple vertex correction. It needs recopulses to do this vertex corrections. Only useful for pure I3Particle seeds, any non-standard topologies won’t work with this.
See also
I3Gulliver
See also
I3GulliverBase
See also
I3SeedServiceBase
See also
I3EventHypothesis
See also
Public Types
-
enum I3TimeShiftType¶
enumerate the different ways to fix the vertex time (of an infinite track)
See also
Values:
-
enumerator TUndefined¶
-
enumerator TNone¶
-
enumerator TMean¶
keep vertex time from first guess
-
enumerator TFirst¶
set vertex time such that the mean time residual is zero
-
enumerator TChargeFraction¶
set vertex time such that the smallest time residual is zero
-
enumerator TDirectChargeFraction¶
set vertex time such that a fixed fraction (e.g. 90%) of the charge has positive time residual
-
enumerator TUndefined¶
-
enum I3PositionShiftType¶
enumerate the different ways to fix the vertex position (of an infinite track)
See also
Values:
-
enumerator PUndefined¶
-
enumerator PNone¶
-
enumerator PCOG¶
keep vertex position from first guess
-
enumerator PUndefined¶
-
enum I3SeedAlternatives¶
The basic behavior of the seed service is to generate only one seed per first guess. We could however also add more seeds, with some trivial transformations. This enum serves to enumerate various ways to generate 1 or more extra seeds.
Values:
-
enumerator SeedAlt_None¶
-
enumerator SeedAlt_Reverse¶
no extra seeds
-
enumerator SeedAlt_Tetrahedron¶
for each first guess, add 1 extra seed, with opposite direction
-
enumerator SeedAlt_Cube¶
add 3 first guesses, with tetrahedric angles
-
enumerator SeedAlt_None¶
Public Functions
-
I3SegmentedRecoSeedService(const std::string &name, const std::vector<std::string> &fgnames, const std::string &inputreadout, double fixedE, std::vector<double> eguesspar, I3PositionShiftType pstype, I3TimeShiftType tstype, I3TimeShiftType alt_tstype, bool speedpolice, double tresmeanmax, double frac, int altMode, bool onlyAlt, bool use_differential_energy_loss_, double shower_spacing, I3Surfaces::SurfacePtr bounding_surface)¶
- Parameters:
name – [in] name of seed service factory, for log messages
fgnames – [in] Labels of first guess fit results
inputreadout – [in] Name of I3RecoPulseMap (for vertex correction)
fixedE – [in] fixed energy value (overrides FG energy)
eguesspar – [in] coefficients of polynomial for E guess (overrides FG energy)
pstype – [in] position shift type
tstype – [in] time shift type
alt_tstype – [in]
speedpolice – [in] flag to enable/disable enforcement of v=c for infinite tracks
tresmeanmax – [in] max mean time residual (protection with TFirst time shift type)
frac – [in] charge fraction value for “TChargeFraction” vertex time correction
altMode – [in] specifies if alternative seeds should be created, with which algorithm
onlyAlt – [in] specifies if nonalternative seeds (based directly on input fg) should be skipped
use_differential_energy_loss_ – [in]
shower_spacing – [in]
bounding_surface – [in]
-
virtual ~I3SegmentedRecoSeedService()¶
base class destructor (idle)
-
virtual unsigned int SetEvent(const I3Frame &f)¶
provide event data
purge old seed tracks
get the first guess tracks, prepare new seeds
do time/space corrections (if wanted)
- Parameters:
f – [in] Frame with event data
- Returns:
number of available seeds
-
virtual I3EventHypothesis GetSeed(unsigned int iseed) const¶
get a seed
-
virtual I3EventHypothesis GetDummy() const¶
get a dummy seed (useful in case of fg failure)
-
virtual void Tweak(I3EventHypothesis &eh) const¶
Space and time coordinates of the vertex are tweaked, for numerical convenience (minimumizer algorithms). If no hits/pulses are configured, then no corrections are applied. Vertex time is corrected such that the timeresiduals are either positive or centered around 0, depending on configuration flags; xyz is chosen such that it is closest to the COG of the hits.
-
void FillInTheBlanks(I3EventHypothesis &eh) const¶
If “FixedEnergy” or “NChEnergyGuessPolynomial” is configured for this service, then this service will use that energy instead of the FG energy.
-
virtual I3EventHypothesis GetCopy(const I3EventHypothesis &eh) const¶
This will just make a (deep) copy of the I3Particle datamember of the I3EventHypothesis. The “nonstd” datamember is not copied.
-
virtual I3EventHypothesis GetCopyNonStd(const I3EventHypothesis &eh) const¶
This will just make a (deep) copy of the I3ParticlePtr and I3FrameObjectPtr datamember of the I3EventHypothesis. The “nonstd” datamember IS now copied.
Public Static Functions
-
static std::vector<I3ParticlePtr> Reverse(I3ParticleConstPtr p)¶
Make a new first guess with direction opposite to input track.
This static method was made public method purely for unit test convenience. If other people would find it useful in other classes, then maybe it would better move to some utility namespace, e.g. I3Calculator.
- Parameters:
p – [in] input particle
- Returns:
a vector with just 1 new particle (pointer) with opposite direction
-
static std::vector<I3ParticlePtr> Tetrahedron(I3ParticleConstPtr p)¶
Make three new first guesses with directions tetrahedric w.r.t. input track (i.e. the original direction and the three new directions are arranged like the points of a tetrahedron).
This static method was made public method purely for unit test convenience. If other people would find it useful in other classes, then maybe it would better move to some utility namespace e.g. I3Calculator.
- Parameters:
p – [in] input particle
- Returns:
a vector with 3 new particles (pointers)
-
static std::vector<I3ParticlePtr> Cube(I3ParticleConstPtr p)¶
Make five new first guesses with directions “cubic” w.r.t. input track (1 opposite + 4 perpendicular)
You should think of the faces of the cube, rather than its points, otherwise this method should have been called “Octahedron”. Are you confused yet?
This static method was made public method purely for unit test convenience. If other people would find it useful in other classes, then maybe it would better move to some utility namespace, e.g. I3Calculator.
- Parameters:
p – [in] input particle
- Returns:
a vector with 5 new particles (pointers)
Private Functions
-
bool SeedOK(const I3Particle &seed) const¶
basic tests: OK fit status, non-NAN datamembers
-
void GetFirstGuesses(const I3Frame &f)¶
get FG tracks from the frame, using the configured track labels
-
void GetNonStdHypothesis(const I3Frame &f)¶
get track and cascades FG from the frame, using the configured track labels, and construct nonstd hypothesis
-
void ShiftVertex(I3Particle &p) const¶
shift the vertex position as close as possible to the COG of the hits (only for inf. tracks)
-
void ShiftTime(I3Particle &p, const I3RecoPulseSeriesMap &hitmap) const¶
adjust the vertex time (such that timeresiduals are more or less OK)
-
double PolynomialEnergyGuess() const¶
guess energy with a polynomial logE = p0+p1*x+p2*x*x+… with x=log(NCh)
-
void GetAlternatives()¶
extend seed set with reverse/tetrahedric/cubic alternatives (if wanted)
- SET_LOGGER ("I3SegmentedRecoSeedService")
Private Members
-
double cascade_spacing_¶
-
I3Surfaces::SurfacePtr bounding_surface_¶
-
bool NonStdHypothesis_¶
-
std::vector<I3EventHypothesis> seeds_¶
-
I3EventHypothesis dummy_¶
-
double fixedEnergy_¶
-
I3PositionShiftType posShiftType_¶
-
I3TimeShiftType timeShiftType_¶
-
I3TimeShiftType altTimeShiftType_¶
-
bool speedPolice_¶
-
double maxTResMean_¶
-
double chargeFraction_¶
-
int altMode_¶
-
bool onlyAlternatives_¶
-
bool use_differential_energy_loss_¶
-
unsigned int nMissingFG_¶
-
unsigned int nBadFG_¶
number of events without FG
-
unsigned int nReadoutMissing_¶
number of bad first guesses (fit status not OK), can be larger than nSetEvent_ in case of multiple FG per event
-
unsigned int nSetEvent_¶
number of events without pulses (while needed)
-
I3Position cog_¶
number of events
-
I3GeometryConstPtr geometry_¶
-
I3RecoPulseSeriesMapConstPtr pulses_¶
-
enum I3TimeShiftType¶
-
class I3SegmentedRecoSeedServiceFactory : public I3ServiceFactory¶
- #include <I3SegmentedRecoSeedServiceFactory.h>
This class installs a I3SegmentedRecoSeedService.
You can configure the I3SegmentedRecoSeedService via eight parameters:
FirstGuess name first guess track/shower (I3Particle or I3Vector<I3Particle>)
FirstGuesses list of names of first guess tracks/showers (I3Particle or I3Vector<I3Particle>)
InputReadout name of hit/pulse series map
FixedEnergy energy value to use in seed (overrides FG energy)
NChEnergyGuessPolynomial smarter energy guess (overrides FG energy), specify coefficients
[c0,c1,c2…] of polynomial log10(E/GeV) = c0 + c1 * x + c2 *x *x + … with x=log10(NCh)
TimeShiftType how to compute the vertex time
AltTimeShiftType how to compute the vertex time for the “alternative” seeds
SpeedPolice flag, whether or not to enforce speed=1 for tracks and speed=0 for cascades.
MaxMeanTimeResidual bound on time shift when using timeshift type TFirst
AddAlternatives Add simple alternative seeds for each first guess; argument is a a string; possibilities: “None” (default, no alternatives), “Reverse” (add a track with the same vertex in the opposite direction), “Cubic” (add 5 more tracks: the reverse track and four perpendicular tracks).
The names of first guesses can refer to either single (I3Particle object in the frame) or multiple results (I3Vector<I3Particle>), both are accepted. By default no energy guess/fix is performed. In case of an NCh-based energy guess you need also to provide a hit/pulse seriesmap (NCh will be set to the size of that map).
If you specify a hit/pulse seriesmap, then the vertex of each seed with type “InfiniteTrack” will always be shifted to the point closest to the COG of the hits/pulses, and the time will be adjusted correspondingly.
Timeshifts are intended to improve the vertex time. If the input track already comes from a log-likelihood reconstruction, then you probably want to set it to “TNone”, no correction. If you set it to “TMean”/”TFirst” then the vertex time is adjusted such that the mean/smallest time residual is zero. For AMANDA analyses in the past TMean was most common (its the default now), for IceCube data it seems that TFirst works better. In order to protect the TFirst shift against anomalous very early hits, there is a maximum on the mean of time residual, which can be configured with MaxMeanTimeResidual.
- Version
$Id: I3SegmentedRecoSeedServiceFactory.h 165272 2018-09-07 10:57:06Z fbradascio $
- Author
boersma
Public Functions
-
virtual ~I3SegmentedRecoSeedServiceFactory()¶
destructor
-
virtual bool InstallService(I3Context &ctx)¶
Installed this objects service into the specified services object.
- Parameters:
ctx – the I3Context into which the service should be installed.
- Returns:
true if the service is successfully installed.
-
virtual void Configure()¶
Configure service prior to installing it.
Private Functions
-
I3SegmentedRecoSeedServiceFactory(const I3SegmentedRecoSeedServiceFactory &rhs)¶
-
I3SegmentedRecoSeedServiceFactory operator=(const I3SegmentedRecoSeedServiceFactory &rhs)¶
no copy constructor
- SET_LOGGER ("I3SegmentedRecoSeedServiceFactory")
Private Members
-
double cascade_spacing_¶
list of first guess names (new option)
-
I3Surfaces::SurfacePtr bounding_surface_¶
-
bool NonStdHypothesis_¶
hits, pulses needed for vertex correction
-
std::string name_¶
use a track (particle hypothesis) and a series of cascades (nonstd hypothesis) as seed
-
unsigned int nContext_¶
name of this service
-
I3SeedServiceBasePtr seedService_¶
number contexts into which a service pointer was installed
-
double fixedEnergy_¶
coefficients of E guess polynomial (if given: overrides FG energy)
-
I3SegmentedRecoSeedService::I3PositionShiftType posShiftType_¶
name of time shift type for the alternative seeds
-
I3SegmentedRecoSeedService::I3TimeShiftType timeShiftType_¶
actual position shift type
-
I3SegmentedRecoSeedService::I3TimeShiftType altTimeShiftType_¶
actual time shift type
-
bool speedPolice_¶
actual time shift type for alternative seeds
-
double maxTResMean_¶
whether or not to force the speed of inf tracks to be c
-
double chargeFraction_¶
max value of average time residual (protection against very early hits in case of TFirst shift type
-
std::string addAlternativesString_¶
charge fraction, for vertex time correction strategy based on tres-ordered hit charges
-
int addAlternatives_¶
name of recipe for generating more than one seed per first guess
-
bool onlyAlternatives_¶
number of the recipe for generating more than one seed per first guess
-
bool use_differential_energy_loss_¶
do not use the seed based on the input track at all
-
struct I3SegmentedSplineRecoDOMCache¶
- #include <I3SegmentedSplineRecoDOMCacheMap.h>
-
class I3SegmentedSplineRecoDOMCacheMap : public std::map<OMKey, struct I3SegmentedSplineRecoDOMCache>¶
- #include <I3SegmentedSplineRecoDOMCacheMap.h>
Public Functions
-
~I3SegmentedSplineRecoDOMCacheMap()¶
-
void UpdateParams(I3GeometryConstPtr geometry, I3CalibrationConstPtr calib, I3DetectorStatusConstPtr det_status, I3RecoPulseSeriesMapConstPtr pulses, double effective_dom_efficiency, double spe_charge_correction, I3TimeWindowSeriesMap excluded_doms)¶
-
~I3SegmentedSplineRecoDOMCacheMap()¶
-
class I3SegmentedSplineRecoLikelihood : public I3EventLogLikelihoodBase, public I3ServiceBase¶
- #include <I3SegmentedSplineRecoLikelihood.h>
Public Functions
-
inline virtual ~I3SegmentedSplineRecoLikelihood()¶
-
virtual void Configure() override¶
I3ServiceBase interface.
-
void SetGeometry(const I3Geometry &f) override¶
-
double GetLogLikelihood(const I3EventHypothesis&) override¶
-
double GetLogLikelihood(const I3EventHypothesis&, I3EventHypothesis*, bool solve_for_energies, double weight = 1) override¶
-
double GetLogLikelihoodCascades(const I3EventHypothesis&, I3EventHypothesis*, I3EventHypothesis*)¶
-
double GetLogLikelihoodHess(const I3EventHypothesis&, I3EventHypothesis&, I3EventHypothesis&) override¶
-
unsigned int GetMultiplicityNPulsesChargeWeight(I3RecoPulseSeriesMapConstPtr pulses)¶
-
unsigned int GetMultiplicityNPulses(I3RecoPulseSeriesMapConstPtr pulses)¶
-
unsigned int GetMultiplicityNCh(I3RecoPulseSeriesMapConstPtr pulses)¶
-
unsigned int GetMultiplicity() override¶
-
void ChainHess(bu::symmetric_matrix<double> &result_hess, bu::symmetric_matrix<double> &hess_upper, std::vector<double> &grad_upper, std::vector<bu::symmetric_matrix<double>> &inner_hess_vec, bu::matrix<double> &inner_jacobian)¶
-
void ChainHessSecond(unsigned int oi1, unsigned int oi2, unsigned int jacobian_row_dim, bu::symmetric_matrix<double> &result_hess, bu::symmetric_matrix<double> &hess_upper, bu::matrix<double> &inner_jacobian)¶
-
inline bool HasGradient() override¶
-
inline bool HasHessian() override¶
Public Members
-
I3PhotonicsServicePtr muon_ps_¶
basic spline, probably bare muon
-
I3PhotonicsServicePtr cascade_ps_¶
basic spline cascade
-
double preJitter_¶
jitter options
-
double postJitter_¶
-
bool PrePulsesReco_¶
Prepulse reconstruction option, set to true by default (true = prepulses are reconstructed)
-
bool MuonReco_¶
Muon reconstruction option, without cascade’s information.
-
double cutoff_distance_¶
cutoff_distance: distance DOM-cascade after which cascade contributions are not considered anymore
-
double noiseRate_¶
-
bool chargeWeight_¶
-
double noise_scale_factor_¶
-
double energy_regularization_¶
-
bool use_cache_¶
-
bool update_eval_cache_¶
-
double effective_dom_efficiency_¶
-
double spe_charge_correction_¶
-
bool EnergyDependentPostJitter_¶
Private Functions
-
double CalculateTimeWindow(I3RecoPulseSeriesMapConstPtr pulses)¶
-
void GaussianIIR1D(double *data, long length, double sigma, int numsteps)¶
-
void GaussianIIR1D_logsumexp(double *data, long length, double sigma, int numsteps)¶
-
double CheckEnDependentJitterSettings()¶
-
double llh_per_pulse(I3OMGeoMap::const_iterator omGeo, std::vector<photo_source> photo_sources, std::vector<I3RecoPulse>::const_iterator hit, int hit_id, double hitTime, std::vector<double> sources_energy, std::vector<double> src_vertex_shifts, double noiseprob, double RDE, double totcharge, std::vector<double> &om_src_geotimes, std::vector<double> &om_src_dists, std::vector<double> &om_src_meanamps, std::vector<std::vector<double>> &om_src_amp_grads, std::vector<bu::symmetric_matrix<double>> &om_src_amp_hesses, std::vector<bu::matrix<double>> &jacobian_casc_to_track_incl_distance, std::vector<bu::symmetric_matrix<double>> &hess_vec_casc_to_track, I3EventHypothesis *gradient, double basic_grad[6], I3MatrixPtr hessian, TableEvalCache &rescale_cache, I3SegmentedSplineRecoDOMCacheMap::const_iterator domcache_iter)¶
Function for cascades
-
void InitSymmMatrix(bu::symmetric_matrix<double> &mat)¶
-
void InitMatrix(bu::matrix<double> &mat)¶
-
void InitMatrixIdentity(bu::matrix<double> &mat)¶
-
void InitBasics(const I3EventHypothesis&)¶
Private Members
-
I3RecoPulseSeriesMapConstPtr pulses_¶
-
I3RecoPulseSeriesMapPtr cleaned_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
-
I3CalibrationConstPtr calib_¶
-
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_¶
-
I3SegmentedSplineRecoDOMCacheMap domCache_¶
-
std::map<OMKey, TableEvalCache> meanamp_cache_¶
-
std::map<OMKey, TableEvalCacheList> prob_cache_¶
-
std::map<OMKey, TableEvalCacheList> cdf_cache_¶
-
std::map<OMKey, TableEvalCache> rescale_cache_¶
-
I3Particle last_optimization_particle_¶
-
std::vector<photo_source> photo_sources_¶
-
unsigned num_cascades_¶
-
const double probPP_ = 0.003¶
-
const double delayPP_ = 31.8 * pow(1345. / 1300, 1. / 2)¶
-
inline virtual ~I3SegmentedSplineRecoLikelihood()¶
-
class I3TimingEquatorFitter : public I3ConditionalModule¶
- #include <I3TimingEquatorFitter.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 ~I3TimingEquatorFitter()¶
-
void Configure()¶
-
void Geometry(I3FramePtr frame)¶
-
void Physics(I3FramePtr frame)¶
-
void Finish()¶
Private Types
Private Functions
-
I3TimingEquatorFitter()¶
-
I3TimingEquatorFitter(const I3TimingEquatorFitter &source)¶
-
I3TimingEquatorFitter &operator=(const I3TimingEquatorFitter &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 ("I3TimingEquatorFitter")
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 ~I3TimingEquatorFitter()¶
-
struct photo_source¶
- #include <I3SegmentedSplineRecoLikelihood.h>
-
struct TableEvalCache¶
- #include <I3SegmentedSplineRecoDOMCacheMap.h>
-
struct TableEvalCacheList¶
- #include <I3SegmentedSplineRecoDOMCacheMap.h>
Public Members
-
std::vector<TableEvalCache> hit_info¶
-
std::vector<TableEvalCache> hit_info¶
-
namespace std
STL namespace.
- file gaussianiir1d.cxx
- #include <math.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
$Date$
- 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
Functions
-
double log_sum_exp(double v1, double v2)¶
- file I3EquatorParametrization.cxx
-
#include “gulliver/I3EventHypothesis.h”#include “icetray/I3Logging.h”#include “icetray/I3SingleServiceFactory.h”#include “dataclasses/physics/I3Particle.h”#include “dataclasses/I3Double.h”#include “dataclasses/I3Vector.h”#include “dataclasses/I3Matrix.h”#include “dataclasses/I3Position.h”#include “dataclasses/I3Direction.h”#include “dataclasses/I3Constants.h”#include “phys-services/I3Calculator.h”#include “gulliver/I3Gulliver.h”#include <boost/numeric/ublas/matrix.hpp>#include <boost/numeric/ublas/matrix_proxy.hpp>#include <boost/numeric/ublas/io.hpp>
implementation of the I3EquatorParametrization class
(c) 2018 the IceCube Collaboration $Id$
- Version
$Revision$
- Date
$Date$
- Author
Federica Bradascio <federica.bradascio.desy.de
Typedefs
-
typedef I3SingleServiceFactory<I3EquatorParametrization, I3ParametrizationBase> I3EquatorParametrizationFactory¶
Functions
-
I3_SERVICE_FACTORY(I3EquatorParametrizationFactory)¶
- file I3EquatorParametrization.h
- #include <string>#include “icetray/IcetrayFwd.h”#include “icetray/I3ServiceBase.h”#include “gulliver/I3ParametrizationBase.h”#include “gulliver/I3FitParameterInitSpecs.h”
declaration of the I3HalfSphereParametrization class
(c) 2007 the IceCube Collaboration
$Id$
- Version
$Revision$
- Date
$Date$
- Author
David Boersma boersma@icecube.wisc.edu
- file I3SegmentedRecoSeedService.cxx
- #include <cmath>#include <vector>#include <string>#include <list>#include “icetray/I3Frame.h”#include “icetray/I3FrameObject.h”#include “icetray/I3Units.h”#include “phys-services/I3Cuts.h”#include “phys-services/I3Calculator.h”#include “dataclasses/physics/I3Particle.h”
copyright (C) 2004 the icecube collaboration $Id:$
- Version
$Revision$
- Date
$Date$
- Author
David Boersma boersma@icecube.wisc.edu
Functions
-
double MeanTimeResidual(const I3Particle &particle, const I3Geometry &geometry, const I3RecoPulseSeriesMap &hitmap)¶
-
double MinimumTimeResidual(const I3Particle &particle, const I3Geometry &geometry, const I3RecoPulseSeriesMap &hitmap, double tresmeanmax)¶
-
double ChargeFractionTimeResidual(const I3Particle &particle, const I3Geometry &geometry, const I3RecoPulseSeriesMap &hitmap, double fraction, bool all, const std::string &logname)¶
- file I3SegmentedRecoSeedService.h
- #include <vector>#include <string>#include “boost/shared_ptr.hpp”#include “icetray/IcetrayFwd.h”#include “icetray/I3DefaultName.h”#include “dataclasses/geometry/I3Geometry.h”#include “dataclasses/physics/I3RecoHit.h”#include “dataclasses/physics/I3RecoPulse.h”#include “phys-services/surfaces/Surface.h”#include “gulliver/I3SeedServiceBase.h”
copyright (C) 2004 the icecube collaboration $Id$
- Version
$Revision$
- Date
$Date$
- Author
David Boersma boersma@icecube.wisc.edu
Functions
-
I3_POINTER_TYPEDEFS(I3SegmentedRecoSeedService)¶
-
I3_DEFAULT_NAME(I3SegmentedRecoSeedService)¶
- file I3SegmentedRecoSeedServiceFactory.cxx
-
#include “icetray/IcetrayFwd.h”
Variables
-
static const char *fgname_optionname = "FirstGuess"
-
static const char *fgnames_optionname = "FirstGuesses"
-
static const char *nonstdhypothesis_optionname = "NonStdHypothesis"¶
-
static const char *inputreadout_optionname = "InputReadout"
-
static const char *fixedE_optionname = "FixedEnergy"
-
static const char *energynch_optionname = "NChEnergyGuessPolynomial"
-
static const char *posshift_optionname = "PositionShiftType"
-
static const char *timeshift_optionname = "TimeShiftType"
-
static const char *alttimeshift_optionname = "AltTimeShiftType"
-
static const char *speed_optionname = "SpeedPolice"
-
static const char *tresmax_optionname = "MaxMeanTimeResidual"
-
static const char *chargefraction_optionname = "ChargeFraction"
-
static const char *addalt_optionname = "AddAlternatives"
-
static const char *onlyalt_optionname = "OnlyAlternatives"
-
static const char *fgname_optionname = "FirstGuess"
- file I3SegmentedRecoSeedServiceFactory.h
- #include <string>#include “icetray/IcetrayFwd.h”#include “icetray/I3ServiceFactory.h”#include “gulliver/I3SeedServiceBase.h”
- file I3SegmentedSplineRecoDOMCacheMap.cxx
- #include “icetray/I3Units.h”#include “dataclasses/calibration/I3Calibration.h”#include “dataclasses/status/I3DetectorStatus.h”#include “dataclasses/geometry/I3Geometry.h”#include “dataclasses/physics/I3Particle.h”#include “dataclasses/physics/I3RecoPulse.h”#include “dataclasses/I3DOMFunctions.h”#include “boost/shared_ptr.hpp”#include “boost/foreach.hpp”#include <string>#include <vector>#include <stack>
- file I3SegmentedSplineRecoDOMCacheMap.h
- #include <string>#include <vector>#include “gulliver/I3EventLogLikelihoodBase.h”#include “gulliver/I3EventHypothesis.h”#include “dataclasses/geometry/I3Geometry.h”#include “dataclasses/calibration/I3Calibration.h”#include “dataclasses/status/I3DetectorStatus.h”#include “dataclasses/physics/I3Particle.h”#include “dataclasses/physics/I3RecoPulse.h”#include <dataclasses/I3TimeWindow.h>#include <boost/numeric/ublas/symmetric.hpp>
- file I3SegmentedSplineRecoLikelihood.cxx
-
#include “dataclasses/I3Vector.h”#include <boost/shared_ptr.hpp>#include <boost/python.hpp>#include “dataclasses/physics/I3MCTreePhysicsLibrary.hh”#include “gulliver/I3Gulliver.h”#include “dataclasses/I3Matrix.h”#include <dataclasses/I3TimeWindow.h>#include <boost/numeric/ublas/symmetric.hpp>#include <boost/numeric/ublas/matrix.hpp>#include <boost/numeric/ublas/matrix_proxy.hpp>#include <boost/numeric/ublas/io.hpp>
Implementation of I3SegmentedSplineRecoLikelihood
$Id$
(c) 2012 The IceCube Collaboration http://www.icecube.wisc.edu
- Version
$Revision$
- Date
$Date$
- file I3SegmentedSplineRecoLikelihood.h
- #include <string>#include <vector>#include “icetray/I3ServiceBase.h”#include “gulliver/I3EventLogLikelihoodBase.h”#include “gulliver/I3EventHypothesis.h”#include “dataclasses/geometry/I3Geometry.h”#include “dataclasses/calibration/I3Calibration.h”#include “dataclasses/physics/I3Particle.h”#include “dataclasses/physics/I3RecoPulse.h”#include “photonics-service/I3PhotonicsService.h”#include “dataclasses/I3Matrix.h”#include <boost/numeric/ublas/symmetric.hpp>#include <boost/numeric/ublas/io.hpp>
Definition of I3SegmentedSplineRecoLikelihood
$Id$
Copyright (C) 2012 The IceCube Collaboration http://www.icecube.wisc.edu
- Version
$Revision$
- Date
$Date$
- Author
Federica Bradascio federica.bradascio@desy.de
Functions
-
I3_POINTER_TYPEDEFS(I3SegmentedSplineRecoLikelihood)¶
- file I3SegmentedSplineRecoLikelihoodFactory.cxx
-
#include “icetray/I3SingleServiceFactory.h”
Implementation of I3SegmentedSplineRecoLikelihoodFactory
$Id$
(c) 2012 The IceCube Collaboration http://www.icecube.wisc.edu
- Version
$Revision$
- Date
- Rcs
2012-09-26 11:16:32 +0200 (Wed, 26 Sep 2012)
Typedefs
-
typedef I3SingleServiceFactory<I3SegmentedSplineRecoLikelihood, I3EventLogLikelihoodBase> I3SegmentedSplineRecoLikelihoodFactory¶
Functions
-
I3_SERVICE_FACTORY(I3SegmentedSplineRecoLikelihoodFactory)¶
- file I3TimingEquatorFitter.cxx
- #include <algorithm>#include <cassert>#include <iomanip>#include <sstream>#include <string>#include <vector>#include <boost/make_shared.hpp>#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”#include “icetray/I3PhysicsTimer.h”
copyright (C) 2005 the icecube collaboration $Id$
- Version
$Revision$
- Date
$Date$
- Author
boersma
Functions
-
I3_MODULE(I3TimingEquatorFitter)¶
- file I3TimingEquatorFitter.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
- dir icetray
- dir private
- dir public
- dir segmented-spline-reco
- dir segmented-spline-reco
- dir segmented-spline-reco