Belle II Software development
KLMDQM2Module.cc
1/**************************************************************************
2 * basf2 (Belle II Analysis Software Framework) *
3 * Author: The Belle II Collaboration *
4 * *
5 * See git log for contributors and copyright holders. *
6 * This file is licensed under LGPL-3.0, see LICENSE.md. *
7 **************************************************************************/
8
9/* Own header. */
10#include <klm/modules/KLMDQM2/KLMDQM2Module.h>
11
12/* Basf2 headers. */
13#include <mdst/dataobjects/Track.h>
14
15/* ROOT headers. */
16#include <TDirectory.h>
17
18/* CLHEP headers. */
19#include <CLHEP/Vector/ThreeVector.h>
20
21/* C++ headers. */
22#include <string>
23
24using namespace Belle2;
25
26//-----------------------------------------------------------------
27// Register the Module
28//-----------------------------------------------------------------
29REG_MODULE(KLMDQM2);
30
31//-----------------------------------------------------------------
32// Implementation
33//-----------------------------------------------------------------
34
37 m_ElementNumbers(&(KLMElementNumbers::Instance())),
40 m_GeometryBKLM{nullptr},
41 m_MatchedHitsBKLM{nullptr},
42 m_AllExtHitsBKLM{nullptr},
43 m_MatchedHitsEKLM{nullptr},
44 m_AllExtHitsEKLM{nullptr},
49{
50 // Set module properties
51 setDescription(R"DOC(Additional Module for KLMDQM plots after HLT filters
52
53 An additional module developed to display plane efficiencies for the KLM during runs (i.e. for online analyses).
54 This module would be called after HLT filter in order to use mumu-tight skim to select reasonable events.
55 The output histograms would be plane efficiencies = MatchedDigits/AllExtits.
56 )DOC");
57
58 // Parameter definitions
59 addParam("MuonListName", m_MuonListName, "Muon list name.",
60 std::string("mu+:all"));
61 addParam("AllowedDistance1D", m_AllowedDistance1D,
62 "Maximal distance in the units of strip number from ExtHit to "
63 "matching KLMDigit (not for multi-strip hits).", double(8));
64 addParam("MinimalMatchingDigits", m_MinimalMatchingDigits,
65 "Minimal number of matching digits.", 0);
66 addParam("MinimalMatchingDigitsOuterLayers",
68 "Minimal number of matching digits in outer layers.", 0);
69 addParam("MinimalMomentumNoOuterLayers", m_MinimalMomentumNoOuterLayers,
70 "Minimal momentum in case there are no hits in outer layers.", 0.0);
71 addParam("RemoveUnusedMuons", m_RemoveUnusedMuons,
72 "Whether to remove unused muons.", false);
73 addParam("IgnoreBackwardPropagation", m_IgnoreBackwardPropagation,
74 "Whether to ignore ExtHits with backward propagation.", false);
76 addParam("histogramDirectoryName", m_HistogramDirectoryName,
77 "Directory for KLM DQM histograms in ROOT file.",
78 std::string("KLMEfficiencyDQM"));
79 addParam("SoftwareTriggerName", m_SoftwareTriggerName,
80 "Software Trigger for event selection",
81 std::string("software_trigger_cut&skim&accept_mumutight"));
82
83}
84
88
90{
91
92 TDirectory* newDirectory{gDirectory->mkdir(m_HistogramDirectoryName.c_str())};
93 TDirectory::TContext context{gDirectory, newDirectory};
94
95
98
99
100 /* Number of hits per channel. */
101 /* KLM General Related. */
102 m_MatchedHitsBKLM = new TH1F("matched_hitsBKLM",
103 "Matched Hits in BKLM Plane",
104 m_PlaneNumBKLM, 0.5, 0.5 + m_PlaneNumBKLM);
105 m_MatchedHitsBKLM->GetXaxis()->SetTitle("Layer Number");
106
107 m_AllExtHitsBKLM = new TH1F("all_ext_hitsBKLM",
108 "All ExtHits in BKLM Plane",
109 m_PlaneNumBKLM, 0.5, 0.5 + m_PlaneNumBKLM);
110 m_AllExtHitsBKLM->GetXaxis()->SetTitle("Layer number");
111
112
113
114 m_MatchedHitsEKLM = new TH1F("matched_hitsEKLM",
115 "Matched Hits in EKLM Plane",
116 m_PlaneNumEKLM, 0.5, m_PlaneNumEKLM + 0.5);
117 m_MatchedHitsEKLM->GetXaxis()->SetTitle("Plane number");
118
119 m_AllExtHitsEKLM = new TH1F("all_ext_hitsEKLM",
120 "All ExtHits in EKLM Plane",
121 m_PlaneNumEKLM, 0.5, m_PlaneNumEKLM + 0.5);
122 m_AllExtHitsEKLM->GetXaxis()->SetTitle("Plane number");
123
124
125
126 /********************/
127 /********************/
128 /* binned by sector */
129 /********************/
130 /********************/
131
132 m_MatchedHitsBKLMSector = new TH1F("matched_hitsBKLMSector",
133 "Matched Hits in BKLM Sector",
134 BKLMMaxSectors, 0.5, 0.5 + BKLMMaxSectors);
135 m_MatchedHitsBKLMSector->GetXaxis()->SetTitle("Sector Number");
136
137 m_AllExtHitsBKLMSector = new TH1F("all_ext_hitsBKLMSector",
138 "All ExtHits in BKLM Sector",
139 BKLMMaxSectors, 0.5, 0.5 + BKLMMaxSectors);
140 m_AllExtHitsBKLMSector->GetXaxis()->SetTitle("Sector number");
141
142
143
144 m_MatchedHitsEKLMSector = new TH1F("matched_hitsEKLMSector",
145 "Matched Hits in EKLM Sector",
146 EKLMMaxSectors, 0.5, EKLMMaxSectors + 0.5);
147 m_MatchedHitsEKLMSector->GetXaxis()->SetTitle("Sector number");
148
149 m_AllExtHitsEKLMSector = new TH1F("all_ext_hitsEKLMSector",
150 "All ExtHits in EKLM Sector",
151 EKLMMaxSectors, 0.5, EKLMMaxSectors + 0.5);
152 m_AllExtHitsEKLMSector->GetXaxis()->SetTitle("Sector number");
153
154}//end of defineHisto
155
157{
158 REG_HISTOGRAM;
159 //inputs
160 m_softwareTriggerResult.isOptional();
161 m_MuonList.isRequired(m_MuonListName);
162 m_Digits.isOptional();
164}
165
167{
168 //start by restarting histograms
169
170 /* KLM General Related. */
171 m_MatchedHitsBKLM->Reset();
172 m_AllExtHitsBKLM->Reset();
173 m_MatchedHitsEKLM->Reset();
174 m_AllExtHitsEKLM->Reset();
175
177 m_AllExtHitsBKLMSector->Reset();
179 m_AllExtHitsEKLMSector->Reset();
180}
181
183{
184 if (triggerFlag() || m_SoftwareTriggerName == "") {
185 unsigned int nMuons = m_MuonList->getListSize();
186 for (unsigned int i = 0; i < nMuons; ++i) {
187 const Particle* muon = m_MuonList->getParticle(i);
192
193 }
194
195 }
196}
197
199{
200}
201
205
207{
208
209 bool passed = false;
211 try {
213 true : false;
214 } catch (const std::out_of_range&) {
215 passed = false;
216 }
217 }
218 return passed;
219
220}
221
223 struct HitData* hitData)
224{
225 for (const KLMDigit& digit : m_Digits) {
226 if (!(digit.getSubdetector() == hitData->subdetector &&
227 digit.getSection() == hitData->section &&
228 digit.getLayer() == hitData->layer &&
229 digit.getSector() == hitData->sector &&
230 digit.getPlane() == hitData->plane) ||
231 !digit.isGood())
232 continue;
233
234 // Defining quantities for distance cut
235 double stripPosition = digit.getStrip();
236 double allowedDistance1D = m_AllowedDistance1D;
237
238 if (digit.isMultiStrip()) {
239 // Due to a firmware bug, we have to be wary with the allowed distance...
240 stripPosition = 0.5 * (digit.getLastStrip() + digit.getStrip());
241 allowedDistance1D *= (digit.getLastStrip() - digit.getStrip() + 1);
242 }
243 if (fabs(stripPosition - hitData->strip) < allowedDistance1D) {
244 hitData->digit = &digit;
245 return;
246 }
247 }
248}
249
251 std::map<KLMPlaneNumber, struct HitData>& hitMap,
252 KLMPlaneNumber planeGlobal, struct HitData* hitData)
253{
254 std::map<KLMPlaneNumber, struct HitData>::iterator it;
255 it = hitMap.find(planeGlobal);
256 /*
257 * There may be more than one such hit e.g. if track crosses the edge
258 * of the strips or WLS fiber groove. Select only one hit per plane.
259 */
260 if (it == hitMap.end()) {
261 hitMap.insert(std::pair<KLMPlaneNumber, struct HitData>(
262 planeGlobal, *hitData));
263 }
264}
265
267 const Particle* muon, TH1F* matchedHitsBKLM, TH1F* allExtHitsBKLM,
268 TH1F* matchedHitsEKLM, TH1F* allExtHitsEKLM, TH1F* matchedHitsBKLMSec, TH1F* allExtHitsBKLMSec,
269 TH1F* matchedHitsEKLMSec, TH1F* allExtHitsEKLMSec)
270{
271 const int nExtrapolationLayers =
273 const Track* track = muon->getTrack();
274 RelationVector<ExtHit> extHits = track->getRelationsTo<ExtHit>();
275 std::map<KLMPlaneNumber, struct HitData> selectedHits;
276 std::map<KLMPlaneNumber, struct HitData>::iterator it;
277 KLMChannelNumber channel;
279 struct HitData hitData, hitDataPrevious;
280 ROOT::Math::XYZVector extHitPosition;
281 CLHEP::Hep3Vector extHitPositionCLHEP, localPosition;
282 int layer;
283 int extHitLayer[nExtrapolationLayers] = {0};
284 int digitLayer[nExtrapolationLayers] = {0};
285 // initialize hitDataPrevious components
286 hitDataPrevious.subdetector = -1;
287 hitDataPrevious.section = -1;
288 hitDataPrevious.sector = -1;
289 hitDataPrevious.layer = -1;
290 for (const ExtHit& hit : extHits) {
291 /*
292 * Choose hits that exit the sensitive volume.
293 * It is not possible to use entry hits because of a bug in Geant4:
294 * the step does not always have a correct status (fGeomBoundary),
295 * and, consequently, ExtHits are not created.
296 */
297 if (hit.getStatus() != EXT_EXIT)
298 continue;
299 /*
300 * Ignore ExtHits with backward propagation. This affects cosmic events
301 * only. The removal of hits with backward propagation is normally
302 * not needed, however, it is added because of backward error propagation
303 * bug in Geant4 10.6.
304 */
306 if (hit.isBackwardPropagated())
307 continue;
308 }
309 KLMPlaneNumber planeGlobal = 0;
310 hitData.hit = &hit;
311 hitData.digit = nullptr;
312 if (hit.getDetectorID() == Const::EDetector::EKLM) {
313 int stripGlobal = hit.getCopyID();
315 m_eklmElementNumbers->stripNumberToElementNumbers(
316 stripGlobal, &hitData.section, &hitData.layer, &hitData.sector,
317 &hitData.plane, &hitData.strip);
318 channel = m_ElementNumbers->channelNumberEKLM(
319 hitData.section, hitData.sector, hitData.layer,
320 hitData.plane, hitData.strip);
321 status = m_ChannelStatus->getChannelStatus(channel);
322 if (status == KLMChannelStatus::c_Unknown)
323 B2FATAL("Incomplete KLM channel status data.");
324 if (status == KLMChannelStatus::c_Normal) {
325 layer = m_ElementNumbers->getExtrapolationLayer(
326 hitData.subdetector, hitData.layer);
327 extHitLayer[layer - 1]++;
328 planeGlobal = m_ElementNumbers->planeNumberEKLM(
329 hitData.section, hitData.sector, hitData.layer,
330 hitData.plane);
331 addHit(selectedHits, planeGlobal, &hitData);
332 }
333 } else if (hit.getDetectorID() == Const::EDetector::BKLM) {
334 int moduleNumber = hit.getCopyID();
337 moduleNumber, &hitData.section, &hitData.sector, &hitData.layer);
339 /*
340 * For scintillators, the plane and strip numbers are recorded
341 * in the copy number.
342 */
344 moduleNumber, &hitData.section, &hitData.sector, &hitData.layer,
345 &hitData.plane, &hitData.strip);
346 channel = m_ElementNumbers->channelNumberBKLM(
347 hitData.section, hitData.sector, hitData.layer,
348 hitData.plane, hitData.strip);
349 status = m_ChannelStatus->getChannelStatus(channel);
350 if (status == KLMChannelStatus::c_Unknown)
351 B2FATAL("Incomplete KLM channel status data.");
352 if (status == KLMChannelStatus::c_Normal) {
353 layer = m_ElementNumbers->getExtrapolationLayer(
354 hitData.subdetector, hitData.layer);
355 extHitLayer[layer - 1]++;
356 planeGlobal = m_ElementNumbers->planeNumberBKLM(
357 hitData.section, hitData.sector, hitData.layer,
358 hitData.plane);
359 addHit(selectedHits, planeGlobal, &hitData);
360 }
361 } else {
362 /* For RPCs, the sensitive volume corresponds to both readout planes. */
363 extHitPosition = hit.getPosition();
364 extHitPositionCLHEP.setX(extHitPosition.X());
365 extHitPositionCLHEP.setY(extHitPosition.Y());
366 extHitPositionCLHEP.setZ(extHitPosition.Z());
367 const bklm::Module* module =
368 m_GeometryBKLM->findModule(hitData.section, hitData.sector,
369 hitData.layer);
370 localPosition = module->globalToLocal(extHitPositionCLHEP);
372 hitData.strip = module->getZStrip(localPosition);
373 /*
374 * FIXME:
375 * There are 2 hits per module in RPCs, but the plane information is
376 * not available in ExtHit. For now, 2 entries are created (one for
377 * each plane) for the first hit, and the second one is removed.
378 */
379 if ((hitData.subdetector == hitDataPrevious.subdetector) &&
380 (hitData.section == hitDataPrevious.section) &&
381 (hitData.sector == hitDataPrevious.sector) &&
382 (hitData.layer == hitDataPrevious.layer))
383 continue;
384 std::memcpy(&hitDataPrevious, &hitData, sizeof(struct HitData));
385 /* The returned strip may be out of the valid range. */
387 hitData.section, hitData.sector, hitData.layer, hitData.plane,
388 hitData.strip, false)) {
389 channel = m_ElementNumbers->channelNumberBKLM(
390 hitData.section, hitData.sector, hitData.layer,
391 hitData.plane, hitData.strip);
392 status = m_ChannelStatus->getChannelStatus(channel);
393 if (status == KLMChannelStatus::c_Unknown)
394 B2FATAL("Incomplete KLM channel status data.");
395 if (status == KLMChannelStatus::c_Normal) {
396 layer = m_ElementNumbers->getExtrapolationLayer(
397 hitData.subdetector, hitData.layer);
398 extHitLayer[layer - 1]++;
399 hitData.localPosition = localPosition.z();
400 planeGlobal = m_ElementNumbers->planeNumberBKLM(
401 hitData.section, hitData.sector, hitData.layer,
402 hitData.plane);
403 addHit(selectedHits, planeGlobal, &hitData);
404 }
405 }
407 hitData.strip = module->getPhiStrip(localPosition);
408 /* The returned strip may be out of the valid range. */
410 hitData.section, hitData.sector, hitData.layer, hitData.plane,
411 hitData.strip, false)) {
412 channel = m_ElementNumbers->channelNumberBKLM(
413 hitData.section, hitData.sector, hitData.layer,
414 hitData.plane, hitData.strip);
415 status = m_ChannelStatus->getChannelStatus(channel);
416 if (status == KLMChannelStatus::c_Unknown)
417 B2FATAL("Incomplete KLM channel status data.");
418 if (status == KLMChannelStatus::c_Normal) {
419 layer = m_ElementNumbers->getExtrapolationLayer(
420 hitData.subdetector, hitData.layer);
421 extHitLayer[layer - 1]++;
422 hitData.localPosition = localPosition.y();
423 planeGlobal = m_ElementNumbers->planeNumberBKLM(
424 hitData.section, hitData.sector, hitData.layer,
425 hitData.plane);
426 addHit(selectedHits, planeGlobal, &hitData);
427 } //end of channel status check
428 } //end of channel number check
429 }//end of detector condition
430 } else
431 continue;
432 }
433 /* Find matching digits. */
434 int nDigits = 0;
435 for (it = selectedHits.begin(); it != selectedHits.end(); ++it) {
436 findMatchingDigit(&(it->second));
437 if (it->second.digit != nullptr) {
438 nDigits++;
439 layer = m_ElementNumbers->getExtrapolationLayer(
440 it->second.subdetector, it->second.layer);
441 digitLayer[layer - 1]++;
442 }
443 }
444 if (nDigits < m_MinimalMatchingDigits)
445 return false;
446 /* Write efficiency histograms */
447 for (it = selectedHits.begin(); it != selectedHits.end(); ++it) {
448 int matchingDigits = 0;
449 int matchingDigitsOuterLayers = 0;
450 int extHitsOuterLayers = 0;
451 layer = m_ElementNumbers->getExtrapolationLayer(
452 it->second.subdetector, it->second.layer) - 1;
453 for (int i = 0; i < nExtrapolationLayers; ++i) {
454 if (i != layer)
455 matchingDigits += digitLayer[i];
456 if (i > layer) {
457 matchingDigitsOuterLayers += digitLayer[i];
458 extHitsOuterLayers += extHitLayer[i];
459 }
460 }
461 /* Check the number of matching digits in other layers. */
462 if (matchingDigits < m_MinimalMatchingDigits)
463 continue;
464 /*
465 * Check the number of matching digits in outer layers relatively to
466 * this hit.
467 */
468 if (matchingDigitsOuterLayers < m_MinimalMatchingDigitsOuterLayers) {
474 if (extHitsOuterLayers >= m_MinimalMatchingDigitsOuterLayers)
475 continue;
476 /*
477 * If the number of ExtHits is insufficient, then check the momentum.
478 * The muons with sufficiently large momentum have a very small
479 * probability to get absorbed in the detector.
480 */
482 continue;
483 }
484 //Filling AddExtHits and MatchedHits histograms
485 if (it->second.subdetector == KLMElementNumbers::c_EKLM) {
486 int planeNum = m_eklmElementNumbers->planeNumber(it->second.section, it->second.layer, it->second.sector, it->second.plane);
487 int sectorNum = (it->second.section - 1) * EKLMElementNumbers::getMaximalSectorNumber() + it->second.sector;
488 allExtHitsEKLM->Fill(planeNum);
489 allExtHitsEKLMSec->Fill(sectorNum);
490 if (it->second.digit) {
491 matchedHitsEKLM->Fill(planeNum);
492 matchedHitsEKLMSec->Fill(sectorNum);
493 }
494 }//end of if loop
495
496
497 else if (it->second.subdetector == KLMElementNumbers::c_BKLM) {
499 it->second.section, it->second.sector, it->second.layer);
500 int sectorGlobal = it->second.section * BKLMElementNumbers::getMaximalSectorNumber() + (it->second.sector);
501 allExtHitsBKLM->Fill(layerGlobal);
502 allExtHitsBKLMSec->Fill(sectorGlobal);
503 if (it->second.digit) {
504 matchedHitsBKLM->Fill(layerGlobal);
505 matchedHitsBKLMSec->Fill(sectorGlobal);
506 }
507 } else {
508 B2FATAL("Had a hit that didn't come from BKLM or EKLM?");
509 }
510
511 } //end of selectedHits for loop
512 return true;
513} //end of collectTrackData
static void channelNumberToElementNumbers(KLMChannelNumber channel, int *section, int *sector, int *layer, int *plane, int *strip)
Get element numbers by channel number.
static bool checkChannelNumber(int section, int sector, int layer, int plane, int strip, bool fatalError=true)
Check channel number.
@ c_FirstRPCLayer
First RPC layer.
static constexpr int getMaximalSectorNumber()
Get maximal sector number (1-based).
static void moduleNumberToElementNumbers(KLMModuleNumber module, int *section, int *sector, int *layer)
Get element numbers by module number.
static constexpr int getMaximalSectorGlobalNumber()
Get maximal sector global number.
static int layerGlobalNumber(int section, int sector, int layer)
Get layer global number.
static constexpr int getMaximalSectorNumber()
Get maximal sector number.
static constexpr int getMaximalSectorGlobalNumberKLMOrder()
Get maximal sector global number with KLM ordering (section, sector).
Store one Ext hit as a ROOT object.
Definition ExtHit.h:31
HistoModule()
Constructor.
Definition HistoModule.h:32
@ c_Normal
Normally operating channel.
@ c_Unknown
Unknown status (no data).
TH1F * m_MatchedHitsBKLMSector
Matched hits in sector for BKLM.
TH1F * m_MatchedHitsEKLM
Matched hits in plane for EKLM.
const bklm::GeometryPar * m_GeometryBKLM
BKLM geometry.
void addHit(std::map< KLMPlaneNumber, struct HitData > &hitMap, KLMPlaneNumber planeGlobal, struct HitData *hitData)
Add hit to map.
StoreArray< KLMDigit > m_Digits
KLM digits.
DBObjPtr< KLMChannelStatus > m_ChannelStatus
Channel status.
TH1F * m_AllExtHitsBKLMSector
Extrapolated hits in sector for BKLM.
virtual void initialize() override
Register input and output data.
double m_MinimalMomentumNoOuterLayers
Minimal momentum in case there are no hits in outer layers.
bool m_RemoveUnusedMuons
Whether to remove unused muons.
virtual void event() override
Selection for mumu_tight_skim, then DQM plot filling.
const KLMElementNumbers * m_ElementNumbers
KLM element numbers.
const KLMPlaneArrayIndex * m_PlaneArrayIndex
Plane array index.
virtual void endRun() override
Called if the current run ends.
const EKLMElementNumbers * m_eklmElementNumbers
EKLM element numbers.
virtual void terminate() override
Called at the end of the event processing.
void findMatchingDigit(struct HitData *hitData)
Find matching digit.
std::string m_SoftwareTriggerName
Software Trigger Name.
TH1F * m_AllExtHitsEKLMSector
Extrapolated hits in sector for EKLM.
TH1F * m_MatchedHitsBKLM
Matched hits in plane for BKLM.
StoreObjPtr< ParticleList > m_MuonList
Muons.
virtual void beginRun() override
Called when entering a new run.
bool m_IgnoreBackwardPropagation
Whether to ignore ExtHits with backward propagation.
double m_AllowedDistance1D
Maximal distance in the units of strip number from ExtHit to matching KLMDigit.
KLMDQM2Module()
Constructor: Sets the description, the properties and the parameters of the module.
~KLMDQM2Module()
Destructor.
TH1F * m_AllExtHitsBKLM
Extrapolated hits in plane for BKLM.
TH1F * m_MatchedHitsEKLMSector
Matched hits in sector for EKLM.
const int m_PlaneNumBKLM
Number of layers/planes for BKLM.
int m_MinimalMatchingDigitsOuterLayers
Minimal number of matching digits in outer layers.
std::string m_MuonListName
Muon list name.
int m_MinimalMatchingDigits
Minimal number of matching digits.
StoreObjPtr< SoftwareTriggerResult > m_softwareTriggerResult
Trigger Information.
const int m_PlaneNumEKLM
Number of layers/planes for EKLM.
TH1F * m_AllExtHitsEKLM
Extrapolated hits in plane for EKLM.
bool collectDataTrack(const Particle *muon, TH1F *matchedHitsBKLM, TH1F *allExtHitsBKLM, TH1F *matchedHitsEKLM, TH1F *allExtHitsEKLM, TH1F *matchedHitsBKLMSec, TH1F *allExtHitsBKLMSec, TH1F *matchedHitsEKLMSec, TH1F *allExtHitsEKLMSec)
Collect the data for one muon.
bool triggerFlag()
Uses TrigResult along with desired software cut to determine whether histograms are filled or not for...
std::string m_HistogramDirectoryName
Directory for KLM DQM histograms in ROOT file.
void defineHisto() override
Definition of the histograms.
KLM digit (class representing a digitized hit in RPCs or scintillators).
Definition KLMDigit.h:29
static constexpr int getMaximalExtrapolationLayer()
Get maximal extrapolation layer.
KLM plane array index.
void setDescription(const std::string &description)
Sets the description of the module.
Definition Module.cc:214
void setPropertyFlags(unsigned int propertyFlags)
Sets the flags for the module properties.
Definition Module.cc:208
@ c_ParallelProcessingCertified
This module can be run in parallel processing mode safely (All I/O must be done through the data stor...
Definition Module.h:80
Class to store reconstructed particles.
Definition Particle.h:76
const Track * getTrack() const
Returns the pointer to the Track object that was used to create this Particle (ParticleType == c_Trac...
Definition Particle.cc:916
double getP() const
Returns momentum magnitude (same as getMomentumMagnitude but with shorter name)
Definition Particle.h:598
Class for type safe access to objects that are referred to in relations.
Class that bundles various TrackFitResults.
Definition Track.h:25
static GeometryPar * instance(void)
Static method to get a reference to the singleton GeometryPar instance.
Define the geometry of a BKLM module Each sector [octant] contains Modules.
Definition Module.h:76
void addParam(const std::string &name, T &paramVariable, const std::string &description, const T &defaultValue)
Adds a new parameter to the module.
Definition Module.h:559
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition Module.h:649
uint16_t KLMChannelNumber
Channel number.
uint16_t KLMPlaneNumber
Plane number.
Abstract base class for different kinds of events.
double localPosition
Local coordinate.
const KLMDigit * digit
Digit.
const ExtHit * hit
Extrapolation hit.