Belle II Software  release-05-01-25
ScintillatorSimulator.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/simulation/ScintillatorSimulator.h>
13 
14 /* KLM headers. */
15 #include <klm/bklm/geometry/GeometryPar.h>
16 #include <klm/eklm/geometry/GeometryData.h>
17 
18 /* Belle 2 headers. */
19 #include <framework/dataobjects/EventMetaData.h>
20 #include <framework/datastore/StoreObjPtr.h>
21 #include <framework/gearbox/Unit.h>
22 #include <framework/logging/Logger.h>
23 
24 /* ROOT headers. */
25 #include <TFile.h>
26 #include <TH1D.h>
27 #include <TRandom.h>
28 
29 /* C++ headers. */
30 #include <algorithm>
31 #include <string>
32 
33 using namespace Belle2;
34 
35 static const char MemErr[] = "Memory allocation error.";
36 
37 static bool compareBKLMSimHits(const BKLMSimHit* hit1, const BKLMSimHit* hit2)
38 {
39  return hit1->getEnergyDeposit() < hit2->getEnergyDeposit();
40 }
41 
42 static bool compareEKLMSimHits(const EKLMSimHit* hit1, const EKLMSimHit* hit2)
43 {
44  return hit1->getEnergyDeposit() < hit2->getEnergyDeposit();
45 }
46 
48 {
50  /*
51  * Here there is a memory leak in case of realloc() failure, but it does not
52  * matter because a fatal error is issued in this case.
53  */
54  /* cppcheck-suppress memleakOnRealloc */
56  (struct Photoelectron*)realloc(m_Photoelectrons,
57  size * sizeof(struct Photoelectron));
58  /* cppcheck-suppress memleakOnRealloc */
60  size * sizeof(int));
61  /* cppcheck-suppress memleakOnRealloc */
63  size * sizeof(int));
64  if (size != 0) {
65  if (m_Photoelectrons == nullptr || m_PhotoelectronIndex == nullptr ||
66  m_PhotoelectronIndex2 == nullptr)
67  B2FATAL(MemErr);
68  }
69 }
70 
73  ScintillatorFirmware* fitter,
74  double digitizationInitialTime, bool debug) :
75  m_DigPar(digPar),
76  m_fitter(fitter),
77  m_DigitizationInitialTime(digitizationInitialTime),
78  m_Debug(debug),
79  m_FPGAStat(c_ScintillatorFirmwareNoSignal),
80  m_npe(0),
81  m_Energy(0),
82  m_MCTime(-1),
83  m_SiPMMCTime(-1)
84 {
85  int i;
86  /* cppcheck-suppress variableScope */
87  double time, attenuationTime;
92  /* Amplitude arrays. */
93  m_amplitudeDirect = (float*)malloc(m_DigPar->getNDigitizations() *
94  sizeof(float));
95  if (m_amplitudeDirect == nullptr)
96  B2FATAL(MemErr);
97  m_amplitudeReflected = (float*)malloc(m_DigPar->getNDigitizations() *
98  sizeof(float));
99  if (m_amplitudeReflected == nullptr)
100  B2FATAL(MemErr);
101  m_amplitude = (float*)malloc(m_DigPar->getNDigitizations() * sizeof(float));
102  if (m_amplitude == nullptr)
103  B2FATAL(MemErr);
104  m_ADCAmplitude = (int*)malloc(m_DigPar->getNDigitizations() * sizeof(int));
105  if (m_ADCAmplitude == nullptr)
106  B2FATAL(MemErr);
107  m_SignalTimeDependence = (double*)malloc((m_DigPar->getNDigitizations() + 1) *
108  sizeof(double));
109  if (m_SignalTimeDependence == nullptr)
110  B2FATAL(MemErr);
112  sizeof(double));
113  if (m_SignalTimeDependenceDiff == nullptr)
114  B2FATAL(MemErr);
115  attenuationTime = 1.0 / m_DigPar->getPEAttenuationFrequency();
116  for (i = 0; i <= m_DigPar->getNDigitizations(); i++) {
117  time = digPar->getADCSamplingTime() * i;
119  exp(-digPar->getPEAttenuationFrequency() * time) * attenuationTime /
120  digPar->getADCSamplingTime();
121  if (i > 0) {
124  }
125  }
126  m_Photoelectrons = nullptr;
127  m_PhotoelectronIndex = nullptr;
128  m_PhotoelectronIndex2 = nullptr;
130 }
131 
132 
134 {
135  free(m_amplitudeDirect);
136  free(m_amplitudeReflected);
137  free(m_amplitude);
138  free(m_ADCAmplitude);
139  free(m_SignalTimeDependence);
140  free(m_SignalTimeDependenceDiff);
141  free(m_Photoelectrons);
142  free(m_PhotoelectronIndex);
143  free(m_PhotoelectronIndex2);
144 }
145 
147  const EKLMChannelData* channelData)
148 {
149  m_Pedestal = channelData->getPedestal();
150  m_PhotoelectronAmplitude = channelData->getPhotoelectronAmplitude();
151  m_Threshold = channelData->getThreshold();
152 }
153 
155 {
156  m_MCTime = -1;
157  m_SiPMMCTime = -1;
158  m_npe = 0;
159  m_Energy = 0;
160  for (int i = 0; i < m_DigPar->getNDigitizations(); i++) {
161  if (m_Debug) {
162  m_amplitudeDirect[i] = 0;
163  m_amplitudeReflected[i] = 0;
164  } else
165  m_amplitude[i] = 0;
166  }
167 }
168 
170  const std::multimap<uint16_t, const BKLMSimHit*>::iterator& firstHit,
171  const std::multimap<uint16_t, const BKLMSimHit*>::iterator& end)
172 {
173  m_stripName = "strip_" + std::to_string(firstHit->first);
174  prepareSimulation();
176  const BKLMSimHit* hit = firstHit->second;
177  const bklm::Module* module =
178  geoPar->findModule(hit->getSection(), hit->getSector(), hit->getLayer());
179  double stripLength =
180  2.0 * (hit->isPhiReadout() ?
181  module->getPhiScintHalfLength(hit->getStrip()) :
182  module->getZScintHalfLength(hit->getStrip()));
183  std::vector<const BKLMSimHit*> hits;
184  for (std::multimap<uint16_t, const BKLMSimHit*>::iterator it = firstHit;
185  it != end; ++it)
186  hits.push_back(it->second);
187  std::sort(hits.begin(), hits.end(), compareBKLMSimHits);
188  for (std::vector<const BKLMSimHit*>::iterator it = hits.begin();
189  it != hits.end(); ++it) {
190  hit = *it;
191  m_Energy = m_Energy + hit->getEnergyDeposit();
192  /* Poisson mean for number of photons. */
193  double nPhotons = hit->getEnergyDeposit() * m_DigPar->getNPEperMeV();
194  /* Fill histograms. */
195  double sipmDistance = hit->getPropagationTime() *
196  m_DigPar->getFiberLightSpeed();
197  double time = hit->getTime() + hit->getPropagationTime();
198  if (m_MCTime < 0) {
199  m_MCTime = hit->getTime();
200  m_SiPMMCTime = time;
201  } else {
202  if (hit->getTime() < m_MCTime)
203  m_MCTime = hit->getTime();
204  if (time < m_SiPMMCTime)
205  m_SiPMMCTime = time;
206  }
207  int generatedPhotons = gRandom->Poisson(nPhotons);
208  generatePhotoelectrons(stripLength, sipmDistance, generatedPhotons,
209  hit->getTime(), false);
210  if (m_DigPar->getMirrorReflectiveIndex() > 0) {
211  generatedPhotons = gRandom->Poisson(nPhotons);
212  generatePhotoelectrons(stripLength, sipmDistance, generatedPhotons,
213  hit->getTime(), true);
214  }
215  }
216  performSimulation();
217 }
218 
220  const std::multimap<uint16_t, const EKLMSimHit*>::iterator& firstHit,
221  const std::multimap<uint16_t, const EKLMSimHit*>::iterator& end)
222 {
223  m_stripName = "strip_" + std::to_string(firstHit->first);
224  prepareSimulation();
225  const EKLMSimHit* hit = firstHit->second;
226  double stripLength = EKLM::GeometryData::Instance().getStripLength(
227  hit->getStrip()) / CLHEP::mm * Unit::mm;
228  std::vector<const EKLMSimHit*> hits;
229  for (std::multimap<uint16_t, const EKLMSimHit*>::iterator it = firstHit;
230  it != end; ++it)
231  hits.push_back(it->second);
232  std::sort(hits.begin(), hits.end(), compareEKLMSimHits);
233  for (std::vector<const EKLMSimHit*>::iterator it = hits.begin();
234  it != hits.end(); ++it) {
235  hit = *it;
236  m_Energy = m_Energy + hit->getEnergyDeposit();
237  /* Poisson mean for number of photons. */
238  double nPhotons = hit->getEnergyDeposit() * m_DigPar->getNPEperMeV();
239  /* Fill histograms. */
240  double sipmDistance = 0.5 * stripLength - hit->getLocalPosition().x();
241  double time = hit->getTime() +
242  sipmDistance / m_DigPar->getFiberLightSpeed();
243  if (m_MCTime < 0)
244  m_MCTime = time;
245  else
246  m_MCTime = time < m_MCTime ? time : m_MCTime;
247  int generatedPhotons = gRandom->Poisson(nPhotons);
248  generatePhotoelectrons(stripLength, sipmDistance, generatedPhotons,
249  hit->getTime(), false);
250  if (m_DigPar->getMirrorReflectiveIndex() > 0) {
251  generatedPhotons = gRandom->Poisson(nPhotons);
252  generatePhotoelectrons(stripLength, sipmDistance, generatedPhotons,
253  hit->getTime(), true);
254  }
255  }
256  performSimulation();
257 }
258 
260 {
261  if (m_Debug) {
262  fillSiPMOutput(m_amplitudeDirect, true, false);
263  fillSiPMOutput(m_amplitudeReflected, false, true);
264  for (int i = 0; i < m_DigPar->getNDigitizations(); i++)
265  m_amplitude[i] = m_amplitudeDirect[i] + m_amplitudeReflected[i];
266  } else
267  fillSiPMOutput(m_amplitude, true, true);
268  /* SiPM noise and ADC. */
269  if (m_DigPar->getMeanSiPMNoise() > 0)
270  addRandomSiPMNoise();
271  simulateADC();
272  m_FPGAStat = m_fitter->fit(m_ADCAmplitude, m_Threshold, &m_FPGAFit);
273  if (m_FPGAStat != c_ScintillatorFirmwareSuccessfulFit)
274  return;
275  if (m_Debug) {
276  if (m_npe >= 10)
277  debugOutput();
278  }
279 }
280 
282 {
283  int i;
284  for (i = 0; i < m_DigPar->getNDigitizations(); i++)
285  m_amplitude[i] = m_amplitude[i] +
286  gRandom->Poisson(m_DigPar->getMeanSiPMNoise());
287 }
288 
290 {
291  int* currentIndexArray, *newIndexArray, *tmpIndexArray;
292  int i, i1, i2, i1Max, i2Max, j, mergeSize;
293  currentIndexArray = m_PhotoelectronIndex;
294  newIndexArray = m_PhotoelectronIndex2;
295  mergeSize = 1;
296  while (mergeSize < nPhotoelectrons) {
297  for (i = 0; i < nPhotoelectrons; i = i + 2 * mergeSize) {
298  i1 = i;
299  j = i;
300  i2 = i + mergeSize;
301  if (i2 > nPhotoelectrons)
302  i2 = nPhotoelectrons;
303  i1Max = i2;
304  i2Max = i2 + mergeSize;
305  if (i2Max > nPhotoelectrons)
306  i2Max = nPhotoelectrons;
307  while (i1 < i1Max || i2 < i2Max) {
308  if (i1 < i1Max) {
309  if (i2 < i2Max) {
310  if (m_Photoelectrons[currentIndexArray[i1]].bin <
311  m_Photoelectrons[currentIndexArray[i2]].bin) {
312  newIndexArray[j] = currentIndexArray[i1];
313  i1++;
314  } else {
315  newIndexArray[j] = currentIndexArray[i2];
316  i2++;
317  }
318  } else {
319  newIndexArray[j] = currentIndexArray[i1];
320  i1++;
321  }
322  } else {
323  newIndexArray[j] = currentIndexArray[i2];
324  i2++;
325  }
326  j++;
327  }
328  }
329  tmpIndexArray = currentIndexArray;
330  currentIndexArray = newIndexArray;
331  newIndexArray = tmpIndexArray;
332  mergeSize = mergeSize * 2;
333  }
334  return currentIndexArray;
335 }
336 
338  double stripLen, double distSiPM, int nPhotons, double timeShift,
339  bool isReflected)
340 {
341  const double maxHitTime = m_DigPar->getNDigitizations() *
342  m_DigPar->getADCSamplingTime();
343  int i;
344  /* cppcheck-suppress variableScope */
345  double hitTime, deExcitationTime, cosTheta, hitDist, selection;
346  double inverseLightSpeed, inverseAttenuationLength;
347  inverseLightSpeed = 1.0 / m_DigPar->getFiberLightSpeed();
348  inverseAttenuationLength = 1.0 / m_DigPar->getAttenuationLength();
349  /* Generation of photoelectrons. */
350  for (i = 0; i < nPhotons; i++) {
351  if (m_npe >= m_PhotoelectronBufferSize)
352  reallocPhotoElectronBuffers(m_npe + 100);
353  cosTheta = gRandom->Uniform(m_DigPar->getMinCosTheta(), 1);
354  if (!isReflected)
355  hitDist = distSiPM / cosTheta;
356  else
357  hitDist = (2.0 * stripLen - distSiPM) / cosTheta;
358  /* Fiber absorption. */
359  selection = gRandom->Uniform();
360  if (selection > exp(-hitDist * inverseAttenuationLength))
361  continue;
362  /* Account for mirror reflective index. */
363  if (isReflected) {
364  selection = gRandom->Uniform();
365  if (selection > m_DigPar->getMirrorReflectiveIndex())
366  continue;
367  }
368  deExcitationTime =
369  gRandom->Exp(m_DigPar->getScintillatorDeExcitationTime()) +
370  gRandom->Exp(m_DigPar->getFiberDeExcitationTime());
371  hitTime = hitDist * inverseLightSpeed + deExcitationTime +
372  timeShift - m_DigitizationInitialTime;
373  if (hitTime >= maxHitTime)
374  continue;
375  if (hitTime >= 0)
376  m_Photoelectrons[m_npe].bin =
377  floor(hitTime / m_DigPar->getADCSamplingTime());
378  else
379  m_Photoelectrons[m_npe].bin = -1;
380  m_Photoelectrons[m_npe].expTime =
381  exp(m_DigPar->getPEAttenuationFrequency() * hitTime);
382  m_Photoelectrons[m_npe].isReflected = isReflected;
383  m_PhotoelectronIndex[m_npe] = m_npe;
384  m_npe++;
385  }
386 }
387 
388 /*
389  * ADC signal corresponding to a photoelectron is
390  *
391  * exp(-(t - tau) / t0),
392  *
393  * where tau is the hit time and t0 = 1.0 / m_DigPar->PEAttenuationFreq. Its
394  * integral from t1 to t2 is
395  *
396  * t0 * exp(-(t1 - tau) / t0) - t0 * exp(-(t2 - tau) / t0).
397  *
398  * The integration is performed over digitization bins from (t_dig * i) to
399  * (t_dig * (i + 1)), where t_dig = m_DigPar->ADCSamplingTime and i is the bin
400  * number. The integrals are
401  *
402  * I1 = t0 - t0 * exp(-(t_dig * (i + 1) - tau) / t0)
403  *
404  * for the first signal bin (t_dig * i <= tau < t_dig * (i + 1)) and
405  *
406  * I2 = t0 * exp(-(t_dig * i - tau) / t0) -
407  * t0 * exp(-(t_dig * (i + 1) - tau) / t0)
408  *
409  * for the following bins. The integrals contain expressions that do not depend
410  * on the hit time and may be calculated preliminary:
411  *
412  * B1 = t0 * exp(-(t_dig * (i + 1) / t0)
413  *
414  * and
415  *
416  * B2 = t0 * exp(-(t_dig * i) / t0) - t0 * exp(-(t_dig * (i + 1)) / t0).
417  *
418  * In terms of B1 and B2 the integrals are equal to
419  *
420  * I1 = t0 - B1 * exp(tau / t0)
421  *
422  * and
423  *
424  * I2 = B2 * exp(tau / t0).
425  *
426  * The sum of integrals for all photoelectrons is
427  *
428  * I = \sum_{hits before this bin} B2 * exp(tau_j / t0) +
429  * \sum_{hits in this bin} t0 - B1 * exp(tau_j / t0),
430  *
431  * where tau_j is the time of j-th hit. The bin-dependent expresions B1 and B2
432  * are the same for all hits:
433  *
434  * I = B2 * \sum_{hits before this bin} exp(tau_j / t0) +
435  * N_i * t0 - B1 * \sum_{hits in this bin} exp(tau_j / t0).
436  *
437  * where N_i is the number of hits in this (i-th) bin.
438  *
439  * In order to get the approximation of the function value, the integrals
440  * over bins are normalized:
441  *
442  * I -> I / t_dig.
443  *
444  * The normalization (1 / t_dig) is included into the definitions of B1 and B2.
445  */
446 void KLM::ScintillatorSimulator::fillSiPMOutput(float* hist, bool useDirect,
447  bool useReflected)
448 {
449  /* cppcheck-suppress variableScope */
450  int i, bin, maxBin;
451  double attenuationTime, sig, expSum;
452  /* cppcheck-suppress variableScope */
453  int ind1, ind2, ind3;
454  int* indexArray;
455  if (m_npe == 0)
456  return;
457  attenuationTime = 1.0 / m_DigPar->getPEAttenuationFrequency();
458  indexArray = sortPhotoelectrons(m_npe);
459  ind1 = 0;
460  expSum = 0;
461  while (1) {
462  ind2 = ind1;
463  bin = m_Photoelectrons[indexArray[ind1]].bin;
464  while (1) {
465  ind2++;
466  if (ind2 == m_npe)
467  break;
468  if (bin != m_Photoelectrons[indexArray[ind2]].bin)
469  break;
470  }
471  /* Now ind1 .. ind2 - photoelectrons in current bin. */
472  for (ind3 = ind1; ind3 != ind2; ind3++) {
473  if (m_Photoelectrons[indexArray[ind3]].isReflected && !useReflected)
474  continue;
475  if (!m_Photoelectrons[indexArray[ind3]].isReflected && !useDirect)
476  continue;
477  if (bin >= 0) {
478  sig = attenuationTime - m_Photoelectrons[indexArray[ind3]].expTime *
479  m_SignalTimeDependence[bin + 1];
480  hist[bin] = hist[bin] + sig;
481  }
482  expSum = expSum + m_Photoelectrons[indexArray[ind3]].expTime;
483  }
484  if (ind2 == m_npe)
485  maxBin = m_DigPar->getNDigitizations() - 1;
486  else
487  maxBin = m_Photoelectrons[indexArray[ind2]].bin;
488  for (i = bin + 1; i <= maxBin; i++) {
489  sig = m_SignalTimeDependenceDiff[i] * expSum;
490  hist[i] = hist[i] + sig;
491  }
492  if (ind2 == m_npe)
493  break;
494  ind1 = ind2;
495  }
496 }
497 
499 {
500  int i;
501  /* cppcheck-suppress variableScope */
502  double amp;
503  if (m_Pedestal == 0 || m_PhotoelectronAmplitude == 0)
504  B2FATAL("Incorrect EKLM ADC simulation parameters.");
505  for (i = 0; i < m_DigPar->getNDigitizations(); i++) {
506  amp = m_Pedestal - m_PhotoelectronAmplitude * m_amplitude[i];
507  if (amp < m_DigPar->getADCSaturation())
508  amp = m_DigPar->getADCSaturation();
509  m_ADCAmplitude[i] = floor(amp);
510  }
511 }
512 
514 {
515  return &m_FPGAFit;
516 }
517 
518 enum KLM::ScintillatorFirmwareFitStatus KLM::ScintillatorSimulator::getFitStatus() const
519 {
520  return m_FPGAStat;
521 }
522 
524 {
525  double intg;
526  intg = m_FPGAFit.getAmplitude();
527  return intg * m_DigPar->getPEAttenuationFrequency() /
528  m_PhotoelectronAmplitude;
529 }
530 
532 {
533  return m_npe;
534 }
535 
537 {
538  return m_Energy;
539 }
540 
542 {
543  int i;
544  std::string str;
546  TFile* hfile = nullptr;
547  TH1D* histAmplitudeDirect = nullptr;
548  TH1D* histAmplitudeReflected = nullptr;
549  TH1D* histAmplitude = nullptr;
550  TH1D* histADCAmplitude = nullptr;
551  try {
552  histAmplitudeDirect =
553  new TH1D("histAmplitudeDirect", m_stripName.c_str(),
554  m_DigPar->getNDigitizations(), 0, m_histRange);
555  histAmplitudeReflected =
556  new TH1D("histAmplitudeReflected", m_stripName.c_str(),
557  m_DigPar->getNDigitizations(), 0, m_histRange);
558  histAmplitude =
559  new TH1D("histAmplitude", m_stripName.c_str(),
560  m_DigPar->getNDigitizations(), 0, m_histRange);
561  histADCAmplitude =
562  new TH1D("histADCAmplitude", m_stripName.c_str(),
563  m_DigPar->getNDigitizations(), 0, m_histRange);
564  } catch (std::bad_alloc& ba) {
565  B2FATAL(MemErr);
566  }
567  for (i = 0; i < m_DigPar->getNDigitizations(); i++) {
568  histAmplitudeDirect->SetBinContent(i + 1, m_amplitudeDirect[i]);
569  histAmplitudeReflected->SetBinContent(i + 1, m_amplitudeReflected[i]);
570  histAmplitude->SetBinContent(i + 1, m_amplitude[i]);
571  histADCAmplitude->SetBinContent(i + 1, m_ADCAmplitude[i]);
572  }
573  str = std::string("experiment_") + std::to_string(event->getExperiment()) +
574  "_run_" + std::to_string(event->getRun()) + "_event_" +
575  std::to_string(event->getEvent()) + "_" + m_stripName + ".root";
576  try {
577  hfile = new TFile(str.c_str(), "NEW");
578  } catch (std::bad_alloc& ba) {
579  B2FATAL(MemErr);
580  }
581  hfile->Append(histAmplitudeDirect);
582  hfile->Append(histAmplitudeReflected);
583  hfile->Append(histAmplitude);
584  hfile->Append(histADCAmplitude);
585  hfile->Write();
586  hfile->Close();
587 }
588 
Belle2::EKLMHitBase::getEnergyDeposit
float getEnergyDeposit() const
Get energy deposit.
Definition: EKLMHitBase.h:110
Belle2::EKLMChannelData
EKLM channel data.
Definition: EKLMChannelData.h:33
Belle2::KLM::ScintillatorSimulator::m_PhotoelectronBufferSize
int m_PhotoelectronBufferSize
Size of photoelectron data buffer.
Definition: ScintillatorSimulator.h:247
Belle2::KLM::ScintillatorSimulator::ScintillatorSimulator
ScintillatorSimulator(const KLMScintillatorDigitizationParameters *digPar, ScintillatorFirmware *fitter, double digitizationInitialTime, bool debug)
Constructor.
Definition: ScintillatorSimulator.cc:71
Belle2::KLMScintillatorDigitizationParameters::getADCPedestal
float getADCPedestal() const
Get ADC pedestal.
Definition: KLMScintillatorDigitizationParameters.h:114
Belle2::KLMScintillatorDigitizationParameters::getPEAttenuationFrequency
float getPEAttenuationFrequency() const
Get attenuation frequency of a single photoelectron pulse.
Definition: KLMScintillatorDigitizationParameters.h:290
Belle2::KLM::ScintillatorSimulator::getNGeneratedPhotoelectrons
int getNGeneratedPhotoelectrons()
Get generated number of photoelectrons.
Definition: ScintillatorSimulator.cc:531
Belle2::KLM::ScintillatorSimulator::debugOutput
void debugOutput()
Debug output (signal and fit result histograms).
Definition: ScintillatorSimulator.cc:541
Belle2::KLMScintillatorDigitizationParameters::getNDigitizations
int getNDigitizations() const
Get number of digitizations (points) in one sample.
Definition: KLMScintillatorDigitizationParameters.h:98
Belle2::KLM::ScintillatorSimulator::m_SignalTimeDependence
double * m_SignalTimeDependence
Buffer for signal time dependence calculation.
Definition: ScintillatorSimulator.h:232
Belle2::KLM::ScintillatorSimulator::~ScintillatorSimulator
~ScintillatorSimulator()
Destructor.
Definition: ScintillatorSimulator.cc:133
Belle2::KLM::ScintillatorSimulator::m_ADCAmplitude
int * m_ADCAmplitude
Digital amplitude.
Definition: ScintillatorSimulator.h:229
Belle2::KLM::ScintillatorSimulator::getEnergy
double getEnergy()
Get total energy deposited in the strip (sum over ssimulation hits).
Definition: ScintillatorSimulator.cc:536
Belle2::KLM::ScintillatorSimulator::Photoelectron
Photoelectron data.
Definition: ScintillatorSimulator.h:44
Belle2::KLM::ScintillatorSimulator::getFPGAFit
KLMScintillatorFirmwareFitResult * getFPGAFit()
Get fit data.
Definition: ScintillatorSimulator.cc:513
Belle2::EKLMChannelData::getThreshold
int getThreshold() const
Get threshold.
Definition: EKLMChannelData.h:82
Belle2::bklm::Module
Define the geometry of a BKLM module Each sector [octant] contains Modules.
Definition: Module.h:86
Belle2::bklm::GeometryPar
Provides BKLM geometry parameters for simulation, reconstruction etc (from Gearbox or DataBase)
Definition: GeometryPar.h:48
Belle2::BKLMSimHit::getEnergyDeposit
double getEnergyDeposit() const
Get energy deposition.
Definition: BKLMSimHit.h:155
Belle2::KLM::ScintillatorSimulator::setChannelData
void setChannelData(const EKLMChannelData *channelData)
Set channel data.
Definition: ScintillatorSimulator.cc:146
Belle2::KLM::ScintillatorSimulator::m_amplitudeDirect
float * m_amplitudeDirect
Analog amplitude (direct).
Definition: ScintillatorSimulator.h:220
Belle2::EKLMChannelData::getPhotoelectronAmplitude
float getPhotoelectronAmplitude() const
Get photoelectron amplitude.
Definition: EKLMChannelData.h:66
Belle2::KLMScintillatorDigitizationParameters::getADCSamplingTime
float getADCSamplingTime() const
Get ADC sampling time in ns.
Definition: KLMScintillatorDigitizationParameters.h:82
Belle2::KLM::ScintillatorSimulator::m_Pedestal
double m_Pedestal
Pedestal.
Definition: ScintillatorSimulator.h:265
Belle2::EKLMSimHit
Class EKLMSimHit stores information on particular Geant step; using information from TrackID and Pare...
Definition: EKLMSimHit.h:41
Belle2::KLM::ScintillatorSimulator::m_Photoelectrons
struct Photoelectron * m_Photoelectrons
Buffer for photoelectron data.
Definition: ScintillatorSimulator.h:238
Belle2::KLM::ScintillatorSimulator::m_PhotoelectronIndex
int * m_PhotoelectronIndex
Buffer for photoelectron indices.
Definition: ScintillatorSimulator.h:241
Belle2::KLM::ScintillatorSimulator::m_Threshold
int m_Threshold
Threshold.
Definition: ScintillatorSimulator.h:271
Belle2::KLMScintillatorDigitizationParameters::getADCPEAmplitude
float getADCPEAmplitude() const
Get ADC photoelectron amplitude.
Definition: KLMScintillatorDigitizationParameters.h:130
Belle2::KLM::ScintillatorSimulator::m_SignalTimeDependenceDiff
double * m_SignalTimeDependenceDiff
Buffer for signal time dependence calculation.
Definition: ScintillatorSimulator.h:235
Belle2::bklm::GeometryPar::instance
static GeometryPar * instance(void)
Static method to get a reference to the singleton GeometryPar instance.
Definition: GeometryPar.cc:30
Belle2::KLM::ScintillatorSimulator::getNPhotoelectrons
double getNPhotoelectrons()
Get number of photoelectrons (fit result).
Definition: ScintillatorSimulator.cc:523
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::KLM::ScintillatorSimulator::generatePhotoelectrons
void generatePhotoelectrons(double stripLen, double distSiPM, int nPhotons, double timeShift, bool isReflected)
Generate photoelectrons.
Definition: ScintillatorSimulator.cc:337
Belle2::StoreObjPtr
Type-safe access to single objects in the data store.
Definition: ParticleList.h:33
Belle2::EKLM::GeometryData::getStripLength
double getStripLength(int strip) const
Get strip length.
Definition: GeometryData.h:73
Belle2::KLMScintillatorDigitizationParameters
Class to store KLM scintillator simulation parameters in the database.
Definition: KLMScintillatorDigitizationParameters.h:33
Belle2::KLM::ScintillatorSimulator::m_amplitudeReflected
float * m_amplitudeReflected
Analog amplitude (reflected).
Definition: ScintillatorSimulator.h:223
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
Belle2::KLM::ScintillatorSimulator::sortPhotoelectrons
int * sortPhotoelectrons(int nPhotoelectrons)
Sort photoelectrons.
Definition: ScintillatorSimulator.cc:289
Belle2::KLM::ScintillatorSimulator::addRandomSiPMNoise
void addRandomSiPMNoise()
Add random noise to the signal (amplitude-dependend).
Definition: ScintillatorSimulator.cc:281
Belle2::KLM::ScintillatorSimulator::m_amplitude
float * m_amplitude
Analog amplitude.
Definition: ScintillatorSimulator.h:226
Belle2::KLM::ScintillatorSimulator::reallocPhotoElectronBuffers
void reallocPhotoElectronBuffers(int size)
Reallocate photoelectron buffers.
Definition: ScintillatorSimulator.cc:47
Belle2::Unit::mm
static const double mm
[millimeters]
Definition: Unit.h:80
Belle2::KLM::ScintillatorSimulator::m_PhotoelectronAmplitude
double m_PhotoelectronAmplitude
Photoelectron amplitude.
Definition: ScintillatorSimulator.h:268
Belle2::KLM::ScintillatorSimulator::fillSiPMOutput
void fillSiPMOutput(float *hist, bool useDirect, bool useReflected)
Fill SiPM output.
Definition: ScintillatorSimulator.cc:446
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::KLM::ScintillatorSimulator::m_PhotoelectronIndex2
int * m_PhotoelectronIndex2
Buffer for photoelectron indices.
Definition: ScintillatorSimulator.h:244
Belle2::BKLMSimHit
Store one simulation hit as a ROOT object.
Definition: BKLMSimHit.h:35
Belle2::KLM::ScintillatorSimulator::m_DigPar
const KLMScintillatorDigitizationParameters * m_DigPar
Parameters.
Definition: ScintillatorSimulator.h:205
Belle2::EKLMChannelData::getPedestal
float getPedestal() const
Get pedestal.
Definition: EKLMChannelData.h:50
Belle2::KLM::ScintillatorSimulator::simulateADC
void simulateADC()
Simulate ADC (create digital signal from analog),.
Definition: ScintillatorSimulator.cc:498
Belle2::KLM::ScintillatorSimulator::prepareSimulation
void prepareSimulation()
Prepare simulation.
Definition: ScintillatorSimulator.cc:154
Belle2::KLMScintillatorFirmwareFitResult
FPGA fit simulation data.
Definition: KLMScintillatorFirmwareFitResult.h:50
Belle2::bklm::GeometryPar::findModule
const Module * findModule(int section, int sector, int layer) const
Get the pointer to the definition of a module.
Definition: GeometryPar.cc:716
Belle2::KLM::ScintillatorSimulator::performSimulation
void performSimulation()
Perform common simulation stage.
Definition: ScintillatorSimulator.cc:259
Belle2::EKLM::GeometryData::Instance
static const GeometryData & Instance(enum DataSource dataSource=c_Database, const GearDir *gearDir=nullptr)
Instantiation.
Definition: GeometryData.cc:35
Belle2::KLMScintillatorDigitizationParameters::getADCThreshold
int getADCThreshold() const
Get ADC readout corresponding to saturation.
Definition: KLMScintillatorDigitizationParameters.h:146
Belle2::KLM::ScintillatorSimulator::m_histRange
double m_histRange
Stands for nDigitizations*ADCSamplingTime.
Definition: ScintillatorSimulator.h:217