Belle II Software  release-05-02-19
EKLMADCModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2013 Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Kirill Chilikin *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 /* Own header. */
12 #include <klm/eklm/modules/EKLMADC/EKLMADCModule.h>
13 
14 /* KLM headers. */
15 #include <klm/eklm/geometry/GeometryData.h>
16 #include <klm/simulation/ScintillatorSimulator.h>
17 
18 /* Belle 2 headers. */
19 #include <framework/gearbox/Unit.h>
20 
21 /* ROOT headers. */
22 #include <TH1F.h>
23 
24 using namespace Belle2;
25 
26 REG_MODULE(EKLMADC)
27 
28 static const char MemErr[] = "Memory allocation error.";
29 
31  Module(),
32  m_fout(nullptr),
33  m_SciSimPar(nullptr),
34  m_hDir(nullptr),
35  m_hRef(nullptr)
36 {
37  setDescription("Standalone generation and studies of ADC output.");
39  addParam("OutputFile", m_out, "Output file.", std::string("EKLMADC.root"));
40  addParam("Mode", m_mode, "Mode (\"Shape\" or \"Strips\").",
41  std::string("Strips"));
42 }
43 
45 {
46 }
47 
48 void EKLMADCModule::generateHistogram(const char* name, double l, double d,
49  int npe)
50 {
51  int j;
52  double t, s;
53  KLM::ScintillatorSimulator fe(m_SciSimPar, nullptr, 0, false);
54  TH1F* h = nullptr;
56  try {
57  h = new TH1F(name, "", m_SciSimPar->getNDigitizations(), 0, t);
58  } catch (std::bad_alloc& ba) {
59  B2FATAL(MemErr);
60  }
61  for (j = 0; j < m_SciSimPar->getNDigitizations(); j++) {
62  m_hDir[j] = 0;
63  m_hRef[j] = 0;
64  }
65  fe.generatePhotoelectrons(l, d, npe, 0, false);
66  fe.generatePhotoelectrons(l, d, npe, 0, true);
67  fe.fillSiPMOutput(m_hDir, true, false);
68  fe.fillSiPMOutput(m_hRef, false, true);
69  s = 0;
70  for (j = 0; j < m_SciSimPar->getNDigitizations(); j++)
71  s = s + m_hDir[j] + m_hRef[j];
72  for (j = 1; j <= m_SciSimPar->getNDigitizations(); j++)
73  h->SetBinContent(j, (m_hDir[j - 1] + m_hRef[j - 1]) / s);
74  h->Write();
75  delete h;
76 }
77 
79 {
80  /* cppcheck-suppress variableScope */
81  char str[32];
82  /* cppcheck-suppress variableScope */
83  int i;
84  /* cppcheck-suppress variableScope */
85  double l;
86  if (!m_SciSimParDatabase.isValid())
87  B2FATAL("EKLM digitization parameters are not available.");
90  try {
91  m_fout = new TFile(m_out.c_str(), "recreate");
92  } catch (std::bad_alloc& ba) {
93  B2FATAL(MemErr);
94  }
95  m_hDir = (float*)malloc(m_SciSimPar->getNDigitizations() * sizeof(float));
96  if (m_hDir == nullptr)
97  B2FATAL(MemErr);
98  m_hRef = (float*)malloc(m_SciSimPar->getNDigitizations() * sizeof(float));
99  if (m_hRef == nullptr)
100  B2FATAL(MemErr);
101  if (m_mode.compare("Strips") == 0) {
102  for (i = 1; i <= geoDat->getNStrips(); i++) {
103  l = geoDat->getStripLength(i) / CLHEP::mm * Unit::mm;
104  snprintf(str, 32, "h%d_near", i);
105  generateHistogram(str, l, 0, 10000);
106  snprintf(str, 32, "h%d_far", i);
107  generateHistogram(str, l, l, 10000);
108  }
109  } else if (m_mode.compare("Shape") == 0) {
111  generateHistogram("FitShape", 0, 0, 1000000);
112  } else
113  B2FATAL("Unknown operation mode.");
114  free(m_hDir);
115  free(m_hRef);
116  m_fout->Close();
117  delete m_fout;
118 }
119 
121 {
122 }
123 
125 {
126 }
127 
129 {
130 }
131 
133 {
134  delete m_SciSimPar;
135 }
136 
Belle2::EKLMADCModule::m_hDir
float * m_hDir
Direct histogram.
Definition: EKLMADCModule.h:108
Belle2::KLM::ScintillatorSimulator
Digitize EKLMSim2Hits to get EKLM StripHits.
Definition: ScintillatorSimulator.h:39
Belle2::Module::setDescription
void setDescription(const std::string &description)
Sets the description of the module.
Definition: Module.cc:216
Belle2::EKLMADCModule::m_mode
std::string m_mode
Operation mode.
Definition: EKLMADCModule.h:93
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::Module::c_ParallelProcessingCertified
@ 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:82
Belle2::EKLMADCModule::m_fout
TFile * m_fout
Output file.
Definition: EKLMADCModule.h:99
Belle2::EKLMADCModule::m_out
std::string m_out
Name of output file.
Definition: EKLMADCModule.h:96
Belle2::KLMScintillatorDigitizationParameters::getNDigitizations
int getNDigitizations() const
Get number of digitizations (points) in one sample.
Definition: KLMScintillatorDigitizationParameters.h:98
Belle2::EKLMADCModule::beginRun
void beginRun() override
Called when entering a new run.
Definition: EKLMADCModule.cc:120
Belle2::Module
Base class for Modules.
Definition: Module.h:74
Belle2::KLMScintillatorDigitizationParameters::getADCSamplingTime
float getADCSamplingTime() const
Get ADC sampling time in ns.
Definition: KLMScintillatorDigitizationParameters.h:82
Belle2::Module::setPropertyFlags
void setPropertyFlags(unsigned int propertyFlags)
Sets the flags for the module properties.
Definition: Module.cc:210
Belle2::EKLMADCModule::event
void event() override
This method is called for each event.
Definition: EKLMADCModule.cc:124
Belle2::EKLMGeometry::getNStrips
int getNStrips() const
Get number of strips.
Definition: EKLMGeometry.h:1741
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
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::EKLMADCModule::EKLMADCModule
EKLMADCModule()
Constructor.
Definition: EKLMADCModule.cc:30
Belle2::EKLMADCModule::m_SciSimPar
KLMScintillatorDigitizationParameters * m_SciSimPar
Scintillator simulation parameters.
Definition: EKLMADCModule.h:105
Belle2::EKLMADCModule::m_SciSimParDatabase
DBObjPtr< KLMScintillatorDigitizationParameters > m_SciSimParDatabase
Scintillator simulation parameters.
Definition: EKLMADCModule.h:102
Belle2::EKLM::GeometryData
EKLM geometry data.
Definition: GeometryData.h:40
Belle2::EKLMADCModule::generateHistogram
void generateHistogram(const char *name, double l, double d, int npe)
Generate output histogram.
Definition: EKLMADCModule.cc:48
Belle2::Module::addParam
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:562
Belle2::Unit::mm
static const double mm
[millimeters]
Definition: Unit.h:80
Belle2::EKLMADCModule::terminate
void terminate() override
This method is called at the end of the event processing.
Definition: EKLMADCModule.cc:132
Belle2::KLMScintillatorDigitizationParameters::setMirrorReflectiveIndex
void setMirrorReflectiveIndex(float reflectiveIndex)
Set mirror reflective index.
Definition: KLMScintillatorDigitizationParameters.h:218
Belle2::EKLMADCModule::~EKLMADCModule
~EKLMADCModule()
Destructor.
Definition: EKLMADCModule.cc:44
Belle2::EKLMADCModule::endRun
void endRun() override
This method is called if the current run ends.
Definition: EKLMADCModule.cc:128
Belle2::EKLM::GeometryData::Instance
static const GeometryData & Instance(enum DataSource dataSource=c_Database, const GearDir *gearDir=nullptr)
Instantiation.
Definition: GeometryData.cc:35
Belle2::EKLMADCModule::initialize
void initialize() override
Initializer.
Definition: EKLMADCModule.cc:78
Belle2::EKLMADCModule::m_hRef
float * m_hRef
Reflected histogram.
Definition: EKLMADCModule.h:111