Belle II Software  release-05-01-25
KLMDigitizerModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2010 Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Timofey Uglov, Kirill Chilikin *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 /* Own header. */
12 #include <klm/modules/KLMDigitizer/KLMDigitizerModule.h>
13 
14 /* KLM headers. */
15 #include <klm/dataobjects/KLMChannelIndex.h>
16 #include <klm/dataobjects/KLMScintillatorFirmwareFitResult.h>
17 #include <klm/simulation/ScintillatorSimulator.h>
18 
19 /* Belle 2 headers. */
20 #include <framework/dataobjects/BackgroundMetaData.h>
21 #include <mdst/dataobjects/MCParticle.h>
22 
23 /* ROOT headers. */
24 #include <TRandom.h>
25 
26 using namespace Belle2;
27 
28 REG_MODULE(KLMDigitizer)
29 
31  Module(),
32  m_ElementNumbers(&(KLMElementNumbers::Instance())),
33  m_eklmElementNumbers(&(EKLMElementNumbers::Instance())),
34  m_ChannelSpecificSimulation(false),
35  m_EfficiencyMode(c_Plane),
36  m_Fitter(nullptr)
37 {
38  setDescription("KLM digitization module: create KLMDigits from BKLMSimHits and EKLMSimHits.");
39  setPropertyFlags(c_ParallelProcessingCertified);
40  addParam("SimulationMode", m_SimulationMode,
41  "Simulation mode (\"Generic\" or \"ChannelSpecific\").",
42  std::string("Generic"));
43  addParam("DigitizationInitialTime", m_DigitizationInitialTime,
44  "Initial digitization time in TDC periods.", -40);
45  addParam("SaveFPGAFit", m_SaveFPGAFit, "Save FPGA fit data and set a relation with KLMDigits.", false);
46  addParam("Efficiency", m_Efficiency,
47  "Efficiency determination mode (\"Strip\" or \"Plane\").",
48  std::string("Plane"));
49  addParam("Debug", m_Debug,
50  "Debug mode (generates additional output files with histograms).",
51  false);
52 }
53 
55 {
56 }
57 
59 {
60  m_bklmSimHits.isRequired();
61  m_eklmSimHits.isRequired();
62  m_Digits.registerInDataStore();
63  m_Digits.registerRelationTo(m_bklmSimHits);
64  m_Digits.registerRelationTo(m_eklmSimHits);
65  if (m_SaveFPGAFit) {
66  m_FPGAFits.registerInDataStore();
67  m_Digits.registerRelationTo(m_FPGAFits);
68  }
69  if (m_SimulationMode == "Generic") {
70  /* Nothing to do. */
71  } else if (m_SimulationMode == "ChannelSpecific") {
73  } else {
74  B2FATAL("Unknown simulation mode: " << m_SimulationMode);
75  }
76  if (m_Efficiency == "Strip")
78  else if (m_Efficiency == "Plane")
80  else
81  B2FATAL("Unknown efficiency mode: " << m_EfficiencyMode);
82 }
83 
85 {
86  KLMChannelIndex klmChannels;
87  for (KLMChannelIndex eklmChannel = klmChannels.beginEKLM();
88  eklmChannel != klmChannels.endEKLM(); ++eklmChannel) {
89  int stripGlobal = m_eklmElementNumbers->stripNumber(
90  eklmChannel.getSection(), eklmChannel.getLayer(),
91  eklmChannel.getSector(), eklmChannel.getPlane(),
92  eklmChannel.getStrip());
93  const EKLMChannelData* channel = m_Channels->getChannelData(stripGlobal);
94  if (channel == nullptr)
95  B2FATAL("Incomplete channel data.");
96  if (channel->getPhotoelectronAmplitude() <= 0) {
97  B2ERROR("Non-positive photoelectron amplitude. The requested "
98  "channel-specific simulation is impossible. "
99  "KLMDigitizer is switched to the generic mode."
100  << LogVar("Section", eklmChannel.getSection())
101  << LogVar("Layer", eklmChannel.getLayer())
102  << LogVar("Sector", eklmChannel.getSector())
103  << LogVar("Plane", eklmChannel.getPlane())
104  << LogVar("Strip", eklmChannel.getStrip()));
106  return;
107  }
108  }
109 }
110 
112 {
113  if (!m_DigPar.isValid())
114  B2FATAL("KLM digitization parameters are not available.");
115  if (!m_TimeConversion.isValid())
116  B2FATAL("KLM time conversion parameters are not available.");
117  if (!m_Channels.isValid())
118  B2FATAL("EKLM channel data are not available.");
119  if (!m_ChannelStatus.isValid())
120  B2FATAL("KLM channel status data are not available.");
121  if (!m_StripEfficiency.isValid())
122  B2FATAL("KLM strip efficiency data are not available.");
125  m_Fitter = new KLM::ScintillatorFirmware(m_DigPar->getNDigitizations());
126 }
127 
128 /*
129  * Digitization. Light propagation into the fiber, SiPM and electronics effects
130  * are simulated in KLM::ScintillatorSimulator class.
131  */
132 
133 bool KLMDigitizerModule::checkActive(uint16_t channel)
134 {
135  enum KLMChannelStatus::ChannelStatus status =
136  m_ChannelStatus->getChannelStatus(channel);
137  if (status == KLMChannelStatus::c_Unknown)
138  B2FATAL("Incomplete KLM channel status data.");
139  return (status != KLMChannelStatus::c_Dead);
140 }
141 
143 {
144  if (isnan(efficiency))
145  B2FATAL("Incomplete KLM efficiency data.");
146  double selection = gRandom->Rndm();
147  return (selection < efficiency);
148 }
149 
151 {
152  int tdc;
153  KLM::ScintillatorSimulator simulator(&(*m_DigPar), m_Fitter, 0, false);
154  std::multimap<uint16_t, const BKLMSimHit*>::iterator it, it2, ub;
155  for (it = m_bklmSimHitChannelMap.begin(); it != m_bklmSimHitChannelMap.end();
156  it = m_bklmSimHitChannelMap.upper_bound(it->first)) {
157  const BKLMSimHit* simHit = it->second;
158  ub = m_bklmSimHitChannelMap.upper_bound(it->first);
159  bool rpc = simHit->inRPC();
160  if (m_EfficiencyMode == c_Strip) {
161  float efficiency = m_StripEfficiency->getBarrelEfficiency(
162  simHit->getSection(), simHit->getSector(),
163  simHit->getLayer(), simHit->getPlane(),
164  simHit->getStrip());
165  if (!efficiencyCorrection(efficiency))
166  continue;
167  }
168  if (rpc) {
171  /* Select hit that has the smallest time. */
172  it2 = it;
173  const BKLMSimHit* hit = it->second;
174  double time = hit->getTime();
175  ++it2;
176  while (it2 != ub) {
177  if (it2->second->getTime() < time) {
178  time = it2->second->getTime();
179  hit = it2->second;
180  }
181  ++it2;
182  }
183  KLMDigit* bklmDigit = m_Digits.appendNew(hit, strip);
184  bklmDigit->addRelationTo(simHit);
185  } else {
186  simulator.simulate(it, ub);
187  if (simulator.getNGeneratedPhotoelectrons() == 0)
188  continue;
189  KLMDigit* bklmDigit = m_Digits.appendNew(simHit);
190  bklmDigit->addRelationTo(simHit);
191  bklmDigit->setMCTime(simulator.getMCTime());
192  bklmDigit->setSiPMMCTime(simulator.getSiPMMCTime());
194  if (simulator.getFitStatus() ==
195  KLM::c_ScintillatorFirmwareSuccessfulFit) {
196  tdc = simulator.getFPGAFit()->getStartTime();
197  bklmDigit->setCharge(simulator.getFPGAFit()->getMinimalAmplitude());
198  } else {
199  tdc = 0;
200  bklmDigit->setCharge(m_DigPar->getADCRange() - 1);
201  }
202  bklmDigit->setTDC(tdc);
203  bklmDigit->setTime(m_TimeConversion->getTimeSimulation(tdc, true));
204  bklmDigit->setFitStatus(simulator.getFitStatus());
205  bklmDigit->setNPhotoelectrons(simulator.getNPhotoelectrons());
206  bklmDigit->setEnergyDeposit(simulator.getEnergy());
207  if (simulator.getFitStatus() == KLM::c_ScintillatorFirmwareSuccessfulFit &&
208  m_SaveFPGAFit) {
210  m_FPGAFits.appendNew(*simulator.getFPGAFit());
211  bklmDigit->addRelationTo(fit);
212  }
213  }
214  }
215 }
216 
218 {
219  uint16_t tdc;
220  int strip;
221  KLM::ScintillatorSimulator simulator(
222  &(*m_DigPar), m_Fitter,
224  const EKLMChannelData* channelData;
225  std::multimap<uint16_t, const EKLMSimHit*>::iterator it, ub;
226  for (it = m_eklmSimHitChannelMap.begin(); it != m_eklmSimHitChannelMap.end();
227  it = m_eklmSimHitChannelMap.upper_bound(it->first)) {
228  const EKLMSimHit* simHit = it->second;
229  ub = m_eklmSimHitChannelMap.upper_bound(it->first);
230  if (m_EfficiencyMode == c_Strip) {
231  float efficiency = m_StripEfficiency->getEndcapEfficiency(
232  simHit->getSection(), simHit->getSector(),
233  simHit->getLayer(), simHit->getPlane(),
234  simHit->getStrip());
235  if (!efficiencyCorrection(efficiency))
236  continue;
237  }
240  simHit->getSection(), simHit->getLayer(), simHit->getSector(),
241  simHit->getPlane(), simHit->getStrip());
242  channelData = m_Channels->getChannelData(strip);
243  if (channelData == nullptr)
244  B2FATAL("Incomplete EKLM channel data.");
245  simulator.setChannelData(channelData);
246  }
247  /* Simulation for a strip. */
248  simulator.simulate(it, ub);
249  if (simulator.getNGeneratedPhotoelectrons() == 0)
250  continue;
251  KLMDigit* eklmDigit = m_Digits.appendNew(simHit);
252  eklmDigit->addRelationTo(simHit);
253  eklmDigit->setMCTime(simulator.getMCTime());
254  eklmDigit->setSiPMMCTime(simulator.getSiPMMCTime());
256  if (simulator.getFitStatus() == KLM::c_ScintillatorFirmwareSuccessfulFit) {
257  tdc = simulator.getFPGAFit()->getStartTime();
258  eklmDigit->setCharge(simulator.getFPGAFit()->getMinimalAmplitude());
259  } else {
260  tdc = 0;
261  eklmDigit->setCharge(m_DigPar->getADCRange() - 1);
262  }
263  eklmDigit->setTDC(tdc);
264  eklmDigit->setTime(m_TimeConversion->getTimeSimulation(tdc, true));
265  eklmDigit->setFitStatus(simulator.getFitStatus());
266  eklmDigit->setNPhotoelectrons(simulator.getNPhotoelectrons());
267  eklmDigit->setEnergyDeposit(simulator.getEnergy());
268  if (simulator.getFitStatus() == KLM::c_ScintillatorFirmwareSuccessfulFit &&
269  m_SaveFPGAFit) {
271  m_FPGAFits.appendNew(*simulator.getFPGAFit());
272  eklmDigit->addRelationTo(fit);
273  }
274  }
275 }
276 
278 {
279  int i;
280  uint16_t channel;
281  m_bklmSimHitChannelMap.clear();
282  m_eklmSimHitChannelMap.clear();
283  if (m_EfficiencyMode == c_Plane) {
284  /* BKLM. */
285  {
286  m_bklmSimHitPlaneMap.clear();
287  for (i = 0; i < m_bklmSimHits.getEntries(); i++) {
288  const BKLMSimHit* hit = m_bklmSimHits[i];
289  if (hit->getStripMin() <= 0)
290  continue;
291  const MCParticle* particle = hit->getRelatedFrom<MCParticle>();
292  /*
293  * We do not simulate the plane efficiency for BKLMSimHits
294  * from beam background because there are no MCParticles associated
295  * to them.
296  */
297  if (particle != nullptr) {
298  uint16_t plane = m_ElementNumbers->planeNumberBKLM(
299  hit->getSection(), hit->getSector(),
300  hit->getLayer(), hit->getPlane());
301  m_bklmSimHitPlaneMap.insert(
302  std::pair<uint16_t, const BKLMSimHit*>(plane, hit));
303  } else {
304  B2ASSERT("The BKLMSimHit is not related to any MCParticle and "
305  "it is also not a beam background hit.",
306  hit->getBackgroundTag() != BackgroundMetaData::bg_none);
308  hit->getSection(), hit->getSector(), hit->getLayer(),
309  hit->getPlane(), hit->getStrip());
310  if (checkActive(channel))
311  m_bklmSimHitChannelMap.insert(
312  std::pair<uint16_t, const BKLMSimHit*>(channel, hit));
313  }
314  }
315  std::multimap<uint16_t, const BKLMSimHit*>::iterator it, it2;
316  std::multimap<const MCParticle*, const BKLMSimHit*> particleHitMap;
317  std::multimap<const MCParticle*, const BKLMSimHit*>::iterator
318  itParticle, it2Particle;
319  it = m_bklmSimHitPlaneMap.begin();
320  while (it != m_bklmSimHitPlaneMap.end()) {
321  particleHitMap.clear();
322  it2 = it;
323  while (true) {
324  const BKLMSimHit* hit = it2->second;
325  const MCParticle* particle = hit->getRelatedFrom<MCParticle>();
326  particleHitMap.insert(
327  std::pair<const MCParticle*, const BKLMSimHit*>(particle, hit));
328  ++it2;
329  if (it2 == m_bklmSimHitPlaneMap.end())
330  break;
331  if (it2->first != it->first)
332  break;
333  }
334  itParticle = particleHitMap.begin();
335  while (itParticle != particleHitMap.end()) {
336  it2Particle = itParticle;
337  const BKLMSimHit* hit = it2Particle->second;
338  float efficiency = m_StripEfficiency->getBarrelEfficiency(
339  hit->getSection(), hit->getSector(),
340  hit->getLayer(), hit->getPlane(),
341  hit->getStripMin());
342  bool hitSelected = efficiencyCorrection(efficiency);
343  while (true) {
344  hit = it2Particle->second;
345  if (hitSelected) {
346  for (int s = hit->getStripMin(); s <= hit->getStripMax(); ++s) {
348  hit->getSection(), hit->getSector(), hit->getLayer(),
349  hit->getPlane(), s);
350  if (checkActive(channel)) {
351  m_bklmSimHitChannelMap.insert(
352  std::pair<uint16_t, const BKLMSimHit*>(channel, hit));
353  }
354  }
355  }
356  ++it2Particle;
357  if (it2Particle == particleHitMap.end())
358  break;
359  if (it2Particle->first != itParticle->first)
360  break;
361  }
362  itParticle = it2Particle;
363  }
364  it = it2;
365  }
366  }
367  /* EKLM. */
368  {
369  m_eklmSimHitPlaneMap.clear();
370  for (i = 0; i < m_eklmSimHits.getEntries(); i++) {
371  const EKLMSimHit* hit = m_eklmSimHits[i];
372  const MCParticle* particle = hit->getRelatedFrom<MCParticle>();
373  /*
374  * We do not simulate the plane efficiency for EKLMSimHits
375  * from beam background because there are no MCParticles
376  * associated to them.
377  */
378  if (particle != nullptr) {
379  uint16_t plane = m_ElementNumbers->planeNumberEKLM(
380  hit->getSection(), hit->getSector(),
381  hit->getLayer(), hit->getPlane());
382  m_eklmSimHitPlaneMap.insert(
383  std::pair<uint16_t, const EKLMSimHit*>(plane, hit));
384  } else {
385  B2ASSERT("The EKLMSimHit is not related to any MCParticle and "
386  "it is also not a beam background hit.",
387  hit->getBackgroundTag() != BackgroundMetaData::bg_none);
389  hit->getSection(), hit->getSector(), hit->getLayer(),
390  hit->getPlane(), hit->getStrip());
391  if (checkActive(channel))
392  m_eklmSimHitChannelMap.insert(
393  std::pair<uint16_t, const EKLMSimHit*>(channel, hit));
394  }
395  }
396  std::multimap<uint16_t, const EKLMSimHit*>::iterator it, it2;
397  std::multimap<const MCParticle*, const EKLMSimHit*> particleHitMap;
398  std::multimap<const MCParticle*, const EKLMSimHit*>::iterator
399  itParticle, it2Particle;
400  it = m_eklmSimHitPlaneMap.begin();
401  while (it != m_eklmSimHitPlaneMap.end()) {
402  particleHitMap.clear();
403  it2 = it;
404  while (true) {
405  const EKLMSimHit* hit = it2->second;
406  const MCParticle* particle = hit->getRelatedFrom<MCParticle>();
407  particleHitMap.insert(
408  std::pair<const MCParticle*, const EKLMSimHit*>(particle, hit));
409  ++it2;
410  if (it2 == m_eklmSimHitPlaneMap.end())
411  break;
412  if (it2->first != it->first)
413  break;
414  }
415  itParticle = particleHitMap.begin();
416  while (itParticle != particleHitMap.end()) {
417  it2Particle = itParticle;
418  const EKLMSimHit* hit = it2Particle->second;
419  float efficiency = m_StripEfficiency->getEndcapEfficiency(
420  hit->getSection(), hit->getSector(),
421  hit->getLayer(), hit->getPlane(),
422  hit->getStrip());
423  bool hitSelected = efficiencyCorrection(efficiency);
424  while (true) {
425  hit = it2Particle->second;
426  if (hitSelected) {
428  hit->getSection(), hit->getSector(), hit->getLayer(),
429  hit->getPlane(), hit->getStrip());
430  if (checkActive(channel)) {
431  m_eklmSimHitChannelMap.insert(
432  std::pair<uint16_t, const EKLMSimHit*>(channel, hit));
433  }
434  }
435  ++it2Particle;
436  if (it2Particle == particleHitMap.end())
437  break;
438  if (it2Particle->first != itParticle->first)
439  break;
440  }
441  itParticle = it2Particle;
442  }
443  it = it2;
444  }
445  }
446  } else {
447  for (i = 0; i < m_bklmSimHits.getEntries(); i++) {
448  const BKLMSimHit* hit = m_bklmSimHits[i];
449  if (hit->inRPC()) {
450  if (hit->getStripMin() <= 0)
451  continue;
452  for (int s = hit->getStripMin(); s <= hit->getStripMax(); ++s) {
454  hit->getSection(), hit->getSector(), hit->getLayer(),
455  hit->getPlane(), s);
456  if (checkActive(channel)) {
457  m_bklmSimHitChannelMap.insert(
458  std::pair<uint16_t, const BKLMSimHit*>(channel, hit));
459  }
460  }
461  } else {
463  hit->getSection(), hit->getSector(), hit->getLayer(),
464  hit->getPlane(), hit->getStrip());
465  if (checkActive(channel)) {
466  m_bklmSimHitChannelMap.insert(
467  std::pair<uint16_t, const BKLMSimHit*>(channel, hit));
468  }
469  }
470  }
471  for (i = 0; i < m_eklmSimHits.getEntries(); i++) {
472  const EKLMSimHit* hit = m_eklmSimHits[i];
474  hit->getSection(), hit->getSector(), hit->getLayer(),
475  hit->getPlane(), hit->getStrip());
476  if (checkActive(channel)) {
477  m_eklmSimHitChannelMap.insert(
478  std::pair<uint16_t, const EKLMSimHit*>(channel, hit));
479  }
480  }
481  }
482  digitizeBKLM();
483  digitizeEKLM();
484 }
485 
487 {
488  delete m_Fitter;
489 }
490 
492 {
493 }
Belle2::KLMChannelIndex::endEKLM
KLMChannelIndex & endEKLM()
Last channel for EKLM.
Definition: KLMChannelIndex.cc:211
Belle2::EKLMHitBase::getSector
int getSector() const
Get sector number.
Definition: EKLMHitBase.h:92
Belle2::EKLMHitBase::getSection
int getSection() const
Get section number.
Definition: EKLMHitBase.h:56
Belle2::EKLMChannelData
EKLM channel data.
Definition: EKLMChannelData.h:33
Belle2::KLMDigitizerModule::initialize
virtual void initialize() override
Initializer.
Definition: KLMDigitizerModule.cc:58
Belle2::EKLMHitBase::getLayer
int getLayer() const
Get layer number.
Definition: EKLMHitBase.h:74
Belle2::KLMDigit::setMCTime
void setMCTime(float time)
Set MC time.
Definition: KLMDigit.h:401
Belle2::KLMDigitizerModule::m_SaveFPGAFit
bool m_SaveFPGAFit
Save FPGA fit data (KLMScintillatorFirmwareFitResult).
Definition: KLMDigitizerModule.h:164
Belle2::KLMDigitizerModule::digitizeBKLM
void digitizeBKLM()
Digitization in BKLM.
Definition: KLMDigitizerModule.cc:150
Belle2::BKLMElementNumbers::getStripByModule
static int getStripByModule(int module)
Get strip number by module identifier.
Definition: BKLMElementNumbers.h:324
Belle2::KLMDigitizerModule::m_SimulationMode
std::string m_SimulationMode
Simulation mode.
Definition: KLMDigitizerModule.h:155
Belle2::KLMElementNumbers::localChannelNumberBKLM
int localChannelNumberBKLM(uint16_t channel) const
Get local BKLM channel number.
Definition: KLMElementNumbers.cc:92
Belle2::EKLMElementNumbers
EKLM element numbers.
Definition: EKLMElementNumbers.h:34
Belle2::KLMDigitizerModule::m_eklmSimHitPlaneMap
std::multimap< uint16_t, const EKLMSimHit * > m_eklmSimHitPlaneMap
Simulation hit map for EKLM (by plane).
Definition: KLMDigitizerModule.h:185
Belle2::KLM::ScintillatorSimulator
Digitize EKLMSim2Hits to get EKLM StripHits.
Definition: ScintillatorSimulator.h:39
Belle2::BKLMSimHit::getSection
int getSection() const
Get section number.
Definition: BKLMSimHit.h:70
Belle2::KLMDigit::setSiPMMCTime
void setSiPMMCTime(float time)
Set SiPM MC time.
Definition: KLMDigit.h:419
Belle2::BKLMSimHit::getLayer
int getLayer() const
Get layer number.
Definition: BKLMSimHit.h:84
Belle2::KLMDigitizerModule::c_Plane
@ c_Plane
Plane.
Definition: KLMDigitizerModule.h:101
Belle2::BKLMSimHit::getStrip
int getStrip() const
Get strip number of this hit.
Definition: BKLMSimHit.h:106
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::KLMChannelStatus::ChannelStatus
ChannelStatus
Channel status.
Definition: KLMChannelStatus.h:44
Belle2::KLMDigitizerModule::m_Digits
StoreArray< KLMDigit > m_Digits
KLM digits.
Definition: KLMDigitizerModule.h:197
Belle2::KLMElementNumbers::planeNumberBKLM
uint16_t planeNumberBKLM(int section, int sector, int layer, int plane) const
Get plane number for BKLM.
Definition: KLMElementNumbers.cc:128
Belle2::KLMDigitizerModule::efficiencyCorrection
bool efficiencyCorrection(float efficiency)
Efficiency correction.
Definition: KLMDigitizerModule.cc:142
Belle2::KLMChannelIndex::beginEKLM
KLMChannelIndex beginEKLM()
First channel for EKLM.
Definition: KLMChannelIndex.cc:205
Belle2::KLMDigit::setNGeneratedPhotoelectrons
void setNGeneratedPhotoelectrons(int nPhotoelectrons)
Set generated number of photoelectrons.
Definition: KLMDigit.h:356
Belle2::KLMDigit::setCharge
void setCharge(uint16_t charge)
Set charge.
Definition: KLMDigit.h:248
Belle2::KLMDigitizerModule::m_bklmSimHitChannelMap
std::multimap< uint16_t, const BKLMSimHit * > m_bklmSimHitChannelMap
Simulation hit map for BKLM (by channel).
Definition: KLMDigitizerModule.h:176
Belle2::KLM::ScintillatorSimulator::getNGeneratedPhotoelectrons
int getNGeneratedPhotoelectrons()
Get generated number of photoelectrons.
Definition: ScintillatorSimulator.cc:531
Belle2::KLMDigitizerModule::m_DigitizationInitialTime
int m_DigitizationInitialTime
Initial digitization time in TDC periods.
Definition: KLMDigitizerModule.h:161
Belle2::RelationsInterface::addRelationTo
void addRelationTo(const RelationsInterface< BASE > *object, float weight=1.0, const std::string &namedRelation="") const
Add a relation from this object to another object (with caching).
Definition: RelationsObject.h:144
Belle2::KLMDigitizerModule::m_Efficiency
std::string m_Efficiency
Efficiency determination mode ("Strip" or "Plane").
Definition: KLMDigitizerModule.h:167
Belle2::KLMDigitizerModule::event
virtual void event() override
This method is called for each event.
Definition: KLMDigitizerModule.cc:277
Belle2::BKLMSimHit::getSector
int getSector() const
Get sector number.
Definition: BKLMSimHit.h:77
Belle2::KLMElementNumbers::channelNumberBKLM
uint16_t channelNumberBKLM(int section, int sector, int layer, int plane, int strip) const
Get channel number for BKLM.
Definition: KLMElementNumbers.cc:50
Belle2::KLM::ScintillatorSimulator::getEnergy
double getEnergy()
Get total energy deposited in the strip (sum over ssimulation hits).
Definition: ScintillatorSimulator.cc:536
Belle2::KLM::ScintillatorSimulator::getFPGAFit
KLMScintillatorFirmwareFitResult * getFPGAFit()
Get fit data.
Definition: ScintillatorSimulator.cc:513
Belle2::KLMDigitizerModule
KLM digitization module.
Definition: KLMDigitizerModule.h:51
Belle2::KLMDigitizerModule::m_StripEfficiency
DBObjPtr< KLMStripEfficiency > m_StripEfficiency
Strip efficiency.
Definition: KLMDigitizerModule.h:146
Belle2::KLMDigit
KLM digit (class representing a digitized hit in RPCs or scintillators).
Definition: KLMDigit.h:40
Belle2::KLMDigit::setTime
void setTime(float time)
Set hit time.
Definition: KLMDigit.h:302
Belle2::KLM::ScintillatorSimulator::setChannelData
void setChannelData(const EKLMChannelData *channelData)
Set channel data.
Definition: ScintillatorSimulator.cc:146
Belle2::KLM::ScintillatorSimulator::getSiPMMCTime
float getSiPMMCTime() const
Get SiPM MC time.
Definition: ScintillatorSimulator.h:159
Belle2::KLMChannelStatus::c_Dead
@ c_Dead
Dead channel (no signal).
Definition: KLMChannelStatus.h:53
Belle2::KLMDigitizerModule::m_eklmSimHitChannelMap
std::multimap< uint16_t, const EKLMSimHit * > m_eklmSimHitChannelMap
Simulation hit map for EKLM (by channel).
Definition: KLMDigitizerModule.h:182
Belle2::KLMDigitizerModule::m_ElementNumbers
const KLMElementNumbers * m_ElementNumbers
Element numbers.
Definition: KLMDigitizerModule.h:149
Belle2::Module
Base class for Modules.
Definition: Module.h:74
Belle2::KLMScintillatorFirmwareFitResult::getStartTime
int getStartTime() const
Get signal start time (in TDC counts).
Definition: KLMScintillatorFirmwareFitResult.h:67
Belle2::KLMDigitizerModule::m_bklmSimHits
StoreArray< BKLMSimHit > m_bklmSimHits
BKLM simulation hits.
Definition: KLMDigitizerModule.h:191
Belle2::EKLMSimHit
Class EKLMSimHit stores information on particular Geant step; using information from TrackID and Pare...
Definition: EKLMSimHit.h:41
Belle2::BKLMSimHit::inRPC
bool inRPC() const
Determine whether this hit is in an RPC or scintillator.
Definition: BKLMSimHit.h:63
Belle2::KLMDigitizerModule::beginRun
virtual void beginRun() override
Called when entering a new run.
Definition: KLMDigitizerModule.cc:111
Belle2::KLMDigitizerModule::m_Fitter
KLM::ScintillatorFirmware * m_Fitter
FPGA fitter.
Definition: KLMDigitizerModule.h:188
Belle2::KLMDigit::setTDC
void setTDC(uint16_t tdc)
Set TDC.
Definition: KLMDigit.h:284
Belle2::KLMDigitizerModule::m_bklmSimHitPlaneMap
std::multimap< uint16_t, const BKLMSimHit * > m_bklmSimHitPlaneMap
Simulation hit map for BKLM (by plane).
Definition: KLMDigitizerModule.h:179
Belle2::KLM::ScintillatorSimulator::getNPhotoelectrons
double getNPhotoelectrons()
Get number of photoelectrons (fit result).
Definition: ScintillatorSimulator.cc:523
Belle2::KLMDigit::setEnergyDeposit
void setEnergyDeposit(float eDep)
Set EnergyDeposit.
Definition: KLMDigit.h:320
Belle2::KLMDigitizerModule::m_FPGAFits
StoreArray< KLMScintillatorFirmwareFitResult > m_FPGAFits
FPGA fits.
Definition: KLMDigitizerModule.h:200
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::KLMDigitizerModule::~KLMDigitizerModule
virtual ~KLMDigitizerModule()
Destructor.
Definition: KLMDigitizerModule.cc:54
Belle2::EKLMSimHit::getPlane
int getPlane() const
Get plane number.
Definition: EKLMSimHit.h:113
Belle2::KLMDigitizerModule::terminate
virtual void terminate() override
This method is called at the end of the event processing.
Definition: KLMDigitizerModule.cc:491
Belle2::KLMDigitizerModule::m_DigPar
DBObjPtr< KLMScintillatorDigitizationParameters > m_DigPar
Digitization parameters.
Definition: KLMDigitizerModule.h:134
Belle2::KLM::ScintillatorSimulator::simulate
void simulate(const std::multimap< uint16_t, const BKLMSimHit * >::iterator &firstHit, const std::multimap< uint16_t, const BKLMSimHit * >::iterator &end)
Simulate BKLM strip.
Definition: ScintillatorSimulator.cc:169
LogVar
Class to store variables with their name which were sent to the logging service.
Definition: LogVariableStream.h:24
Belle2::EKLMElementNumbers::stripNumber
int stripNumber(int section, int layer, int sector, int plane, int strip) const
Get strip number.
Definition: EKLMElementNumbers.cc:212
Belle2::KLMDigitizerModule::m_EfficiencyMode
EfficiencyMode m_EfficiencyMode
Efficiency determination mode (converted from the string parameter).
Definition: KLMDigitizerModule.h:170
Belle2::KLMDigitizerModule::m_eklmSimHits
StoreArray< EKLMSimHit > m_eklmSimHits
EKLM simulation hits.
Definition: KLMDigitizerModule.h:194
Belle2::KLMDigitizerModule::digitizeEKLM
void digitizeEKLM()
Digitization in EKLM.
Definition: KLMDigitizerModule.cc:217
Belle2::BackgroundMetaData::bg_none
@ bg_none
No background.
Definition: BackgroundMetaData.h:41
Belle2::KLMDigitizerModule::checkActive
bool checkActive(uint16_t channel)
Check if channel is active (status is not KLMChannelStatus::c_Dead).
Definition: KLMDigitizerModule.cc:133
Belle2::KLMDigitizerModule::m_TimeConversion
DBObjPtr< KLMTimeConversion > m_TimeConversion
Time conversion.
Definition: KLMDigitizerModule.h:137
Belle2::KLMElementNumbers::planeNumberEKLM
uint16_t planeNumberEKLM(int section, int sector, int layer, int plane) const
Get channel number for EKLM.
Definition: KLMElementNumbers.cc:136
Belle2::EKLMSimHit::getStrip
int getStrip() const
Get strip number.
Definition: EKLMSimHit.h:131
Belle2::KLMScintillatorFirmwareFitResult::getMinimalAmplitude
int getMinimalAmplitude() const
Get minimal amplitude (ADC output).
Definition: KLMScintillatorFirmwareFitResult.h:118
Belle2::KLMDigitizerModule::c_Strip
@ c_Strip
Strip.
Definition: KLMDigitizerModule.h:98
Belle2::KLMDigit::setNPhotoelectrons
void setNPhotoelectrons(float nPhotoelectrons)
Set number of photoelectrons.
Definition: KLMDigit.h:338
Belle2::KLMDigit::setFitStatus
void setFitStatus(int s)
Set fit status.
Definition: KLMDigit.h:383
Belle2::KLMElementNumbers::channelNumberEKLM
uint16_t channelNumberEKLM(int section, int sector, int layer, int plane, int strip) const
Get channel number for EKLM.
Definition: KLMElementNumbers.cc:64
Belle2::KLMChannelStatus::c_Unknown
@ c_Unknown
Unknown status (no data).
Definition: KLMChannelStatus.h:47
Belle2::KLMChannelIndex
KLM channel index.
Definition: KLMChannelIndex.h:33
Belle2::KLM::ScintillatorSimulator::getFitStatus
enum ScintillatorFirmwareFitStatus getFitStatus() const
Get fit status.
Definition: ScintillatorSimulator.cc:518
Belle2::KLM::ScintillatorFirmware
FPGA fitter class.
Definition: ScintillatorFirmware.h:33
Belle2::KLMDigitizerModule::m_eklmElementNumbers
const EKLMElementNumbers * m_eklmElementNumbers
EKLM element numbers.
Definition: KLMDigitizerModule.h:152
Belle2::BKLMSimHit
Store one simulation hit as a ROOT object.
Definition: BKLMSimHit.h:35
Belle2::MCParticle
A Class to store the Monte Carlo particle information.
Definition: MCParticle.h:43
Belle2::KLMElementNumbers
KLM element numbers.
Definition: KLMElementNumbers.h:37
Belle2::KLMDigitizerModule::m_Debug
bool m_Debug
Use debug mode in EKLM::ScintillatorSimulator or not.
Definition: KLMDigitizerModule.h:173
Belle2::KLMScintillatorFirmwareFitResult
FPGA fit simulation data.
Definition: KLMScintillatorFirmwareFitResult.h:50
Belle2::KLM::ScintillatorSimulator::getMCTime
float getMCTime() const
Get MC time.
Definition: ScintillatorSimulator.h:150
Belle2::KLMDigitizerModule::checkChannelParameters
void checkChannelParameters()
Check channel parameters for channel-specific simulation.
Definition: KLMDigitizerModule.cc:84
Belle2::KLMDigitizerModule::m_ChannelSpecificSimulation
bool m_ChannelSpecificSimulation
Whether the simulation is channel-specific.
Definition: KLMDigitizerModule.h:158
Belle2::KLMDigitizerModule::m_ChannelStatus
DBObjPtr< KLMChannelStatus > m_ChannelStatus
Channel status.
Definition: KLMDigitizerModule.h:143
Belle2::KLMDigitizerModule::m_Channels
DBObjPtr< EKLMChannels > m_Channels
Channel data.
Definition: KLMDigitizerModule.h:140
Belle2::KLMDigitizerModule::endRun
virtual void endRun() override
This method is called if the current run ends.
Definition: KLMDigitizerModule.cc:486
Belle2::BKLMSimHit::getPlane
int getPlane() const
Get plane number.
Definition: BKLMSimHit.h:91