Belle II Software  release-05-02-19
BKLMDigitAnalyzerModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2018 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Giacomo De Pietro *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 /* Own header. */
12 #include <klm/bklm/modules/bklmDigitAnalyzer/BKLMDigitAnalyzerModule.h>
13 
14 /* Belle 2 headers. */
15 #include <framework/datastore/StoreObjPtr.h>
16 #include <framework/dataobjects/EventMetaData.h>
17 
18 /* ROOT headers. */
19 #include <TGaxis.h>
20 
21 using namespace Belle2;
22 
23 
24 REG_MODULE(BKLMDigitAnalyzer)
25 
26 
28  m_runNumber(0),
29  m_outputRootFile(nullptr),
30  m_histoList(nullptr)
31 {
32  for (int i = 0; i < 2; ++i) {
33  m_histoLayerVsSector[i] = nullptr;
34  m_histoLayerVsSectorPerPlane[i][0] = nullptr;
35  m_histoLayerVsSectorPerPlane[i][1] = nullptr;
36  for (int s = 0; s < 8; ++s) {
37  m_histoLayer[i][s] = nullptr;
38  for (int p = 0; p < 2; ++p) {
39  m_histoChannel[i][s][p] = nullptr;
40  m_histoStrip[i][s][p] = nullptr;
41  m_histoTdc[i][s][p] = nullptr;
42  m_histoCTimeDiff[i][s][p] = nullptr;
43  }
44  }
45  }
46 
47  // Set module properties
48  setDescription("Module useful to quickly analyze BKLM unpacked data.");
49 
50  // Parameter definitions
51  addParam("outputRootName", m_outputRootName, "Name of output .root file (without .root!)", std::string("bklmHitmap"));
52 
53 }
54 
56 {
57 }
58 
60 {
61  m_digit.isRequired();
62  m_digitRaw.isRequired();
63  m_digitOutOfRange.isRequired("KLMDigitsOutOfRange");
64  m_digitEventInfo.isRequired();
65 }
66 
68 {
69  time_t rawTime;
70  time(&rawTime);
71  struct tm* tm = gmtime(&rawTime);
72 
73  StoreObjPtr<EventMetaData> eventMetaData("EventMetaData", DataStore::c_Event);
74  m_runNumber = eventMetaData->getRun();
75 
76  TString runNumberTString(toString(m_runNumber).c_str());
77 
78  TString outputRootNameTString(m_outputRootName);
79  outputRootNameTString += "_run" + runNumberTString + ".root";
80 
81  TString label[2] = {"BF", "BB"};
82 
83  m_outputRootFile = new TFile(outputRootNameTString, "RECREATE");
84  B2INFO("BKLMDigitAnalyzer:: the output file '" << outputRootNameTString.Data() << "' will be created for run " << m_runNumber);
85 
86  int exp = -1;
87  int run = -1;
88  int year = -1;
89  int month = -1;
90  int day = -1;
91  int hour = -1;
92  int min = -1;
93  int sec = -1;
94  m_extraInfo = new TTree("extraInfo", "Extra informations");
95  m_extraInfo->Branch("exp", &exp, "exp/I");
96  m_extraInfo->Branch("run", &run, "run/I");
97  m_extraInfo->Branch("year", &year, "year/I");
98  m_extraInfo->Branch("month", &month, "month/I");
99  m_extraInfo->Branch("day", &day, "day/I");
100  m_extraInfo->Branch("hour", &hour, "hour/I");
101  m_extraInfo->Branch("min", &min, "min/I");
102  m_extraInfo->Branch("sec", &sec, "sec/I");
103  exp = eventMetaData->getExperiment();
104  run = m_runNumber;
105  year = tm->tm_year + 1900;
106  month = tm->tm_mon + 1;
107  day = tm->tm_mday;
108  hour = tm->tm_hour;
109  min = tm->tm_min;
110  sec = tm->tm_sec;
111  m_extraInfo->Fill();
112 
113  m_histoList = new TList;
114 
115  for (int fb = 0; fb < 2; fb++) {
116 
117  m_histoLayerVsSector[fb] = createTH2("SectLay" + label[fb], label[fb] + " Hitmap -- run" + runNumberTString, 31, -0.5, 15.,
118  "Layer (0-based)", 17, -0.5, 8., "BF Sector", 0, m_histoList);
119 
120  for (int isRPCorPhi = 0; isRPCorPhi < 2; isRPCorPhi++) {
121 
122  if (isRPCorPhi == 0)
123  m_histoLayerVsSectorPerPlane[fb][isRPCorPhi] = createTH2("SectLayPlaneZ" + label[fb],
124  label[fb] + " Hitmap of plane z -- run" + runNumberTString, 31, -0.5, 15., "Layer (0-based)", 17, -0.5, 8., "BF Sector", 0,
125  m_histoList);
126  else
127  m_histoLayerVsSectorPerPlane[fb][isRPCorPhi] = createTH2("SectLayPlanePhi" + label[fb],
128  label[fb] + " Hitmap of plane phi -- run" + runNumberTString, 31, -0.5, 15., "Layer (0-based)", 17, -0.5, 8., "BF Sector", 0,
129  m_histoList);
130 
131  for (int iSector = 0; iSector < 8; iSector++) {
132 
133  TString iSectorTString(toString(iSector).c_str());
134  TString nameSector = label[fb] + iSectorTString;
135 
136  if (isRPCorPhi == 0) {
137 
138  // Create the histogram belowo only one time
139  m_histoLayer[fb][iSector] = createTH1("Layer" + nameSector, nameSector + " Layer -- run" + runNumberTString, 31, -0.5, 15.,
140  "Layer (0-based)", "Counts", 1, m_histoList);
141 
142  m_histoChannel[fb][iSector][isRPCorPhi] = createTH2("PlaneZ" + nameSector,
143  nameSector + " Plane z, electronic channels -- run" + runNumberTString, 31,
144  -0.5, 15., "Layer (0-based)", 130, -0.5, 64.5, "Channel", 1, m_histoList);
145 
146  m_histoStrip[fb][iSector][isRPCorPhi] = createTH2("PlaneZStrip" + nameSector,
147  nameSector + " Plane z, strips -- run" + runNumberTString,
148  31, -0.5, 15., "Layer (0-based)", 130, -0.5, 64.5, "Strip", 1, m_histoList);
149 
150  m_histoTdc[fb][iSector][isRPCorPhi] = createTH1("SciTdc" + nameSector,
151  nameSector + " TDC (Scintillators) -- run" + runNumberTString, 60, 0, 30, "TDC", "Counts", 1, m_histoList);
152 
153  m_histoCTimeDiff[fb][iSector][isRPCorPhi] = createTH1("SciCTimeDiff" + nameSector,
154  nameSector + " CTime diff. (Scintillators) -- run" + runNumberTString, 200, -200, -100, "Sci_CTime - Trg_CTime [ns]", "Counts", 1,
155  m_histoList);
156 
157  } else {
158 
159  m_histoChannel[fb][iSector][isRPCorPhi] = createTH2("PlanePhi" + nameSector,
160  nameSector + " Plane phi, electronic channels -- run" + runNumberTString,
161  31, -0.5, 15., "Layer (0-based)", 130, -0.5, 64.5, "Channel", 1, m_histoList);
162 
163  m_histoStrip[fb][iSector][isRPCorPhi] = createTH2("PlanePhiStrip" + nameSector,
164  nameSector + " Plane phi, strips -- run" + runNumberTString,
165  31, -0.5, 15., "Layer (0-based)", 130, -0.5, 64.5, "Strip", 1, m_histoList);
166 
167  m_histoTdc[fb][iSector][isRPCorPhi] = createTH1("RPCTdc" + nameSector, nameSector + " TDC (RPCs) -- run" + runNumberTString, 230,
168  -100, 2200, "TDC", "Counts", 1, m_histoList);
169 
170  m_histoCTimeDiff[fb][iSector][isRPCorPhi] = createTH1("RPCCTimeDiff" + nameSector,
171  nameSector + " CTime diff. (RPCs) -- run" + runNumberTString, 350, -3000, 500, "RPC_CTime - Trg_CTime [ns]", "Counts", 1,
172  m_histoList);
173 
174  }
175  }
176  }
177  }
178 }
179 
181 {
182  for (int i = 0; i < m_digitEventInfo.getEntries(); i++) {
183 
184  KLMDigitEventInfo* digitEventInfo = m_digitEventInfo[i];
185  // Some warnings (they should never appear, but it's better to be sure)
186  if ((digitEventInfo->getRPCHits() + digitEventInfo->getSciHits()) != (int)digitEventInfo->getRelationsFrom<KLMDigit>().size())
187  B2WARNING("BKLMDigitAnalyzer:: the total number of BKLMDigit differs from the sum of RPC and scintillator hits stored in BKLMEventDigitDebug!");
188  if (digitEventInfo->getOutOfRangeHits() != (int)digitEventInfo->getRelationsTo<KLMDigit>("BKLMDigitsOutOfRange").size())
189  B2WARNING("BKLMDigitAnalyzer:: the total number of BKLMDigit differs from the number of outOfRange-flagged hits stored in BKLMEventDigitDebug!");
190 
191  for (const KLMDigit& digit : digitEventInfo->getRelationsFrom<KLMDigit>()) {
192  if (digit.getSubdetector() != KLMElementNumbers::c_BKLM)
193  continue;
194 
195  KLMDigitRaw* digitRaw = digit.getRelatedTo<KLMDigitRaw>();
196 
197  m_histoLayerVsSector[1 - digit.getSection()]->Fill(digit.getLayer() - 1, digit.getSector() - 1);
198 
199  m_histoLayerVsSectorPerPlane[1 - digit.getSection()][digit.isPhiReadout()]->Fill(digit.getLayer() - 1, digit.getSector() - 1);
200 
201  m_histoLayer[1 - digit.getSection()][digit.getSector() - 1]->Fill(digit.getLayer() - 1);
202 
203  m_histoChannel[1 - digit.getSection()][digit.getSector() - 1][digit.isPhiReadout()]->Fill(digit.getLayer() - 1,
204  digitRaw->getChannel() - 1);
205 
206  m_histoStrip[1 - digit.getSection()][digit.getSector() - 1][digit.isPhiReadout()]->Fill(digit.getLayer() - 1,
207  digit.getStrip() - 1);
208 
209  // getTime() retruns the TDC
210  m_histoTdc[1 - digit.getSection()][digit.getSector() - 1][digit.inRPC()]->Fill(digit.getTime());
211 
212  m_histoCTimeDiff[1 - digit.getSection()][digit.getSector() - 1][digit.inRPC()]->Fill(digit.getCTime() -
213  digitEventInfo->getIntTriggerCTime());
214 
215  }
216  }
217 }
218 
220 {
221  StoreObjPtr<EventMetaData> eventMetaData("EventMetaData", DataStore::c_Event);
222 
223  // Save the .root file
224  if (m_outputRootFile != nullptr) {
225  m_outputRootFile->cd();
226 
227  m_extraInfo->Write();
228 
229  TIter nextHisto(m_histoList);
230  TObject* obj;
231  while ((obj = nextHisto()))
232  obj->Write("", TObject::kWriteDelete);
233  m_outputRootFile->Close();
234  }
235 }
236 
238 {
239 }
240 
241 TH1F* BKLMDigitAnalyzerModule::createTH1(const char* name, const char* title, Int_t nbinsX, Double_t minX, Double_t maxX,
242  const char* titleX, const char* titleY, bool drawStat, TList* histoList)
243 {
244  TH1F* hist = new TH1F(name, title, nbinsX, minX, maxX);
245  hist->GetXaxis()->SetTitle(titleX);
246  hist->GetYaxis()->SetTitle(titleY);
247  TGaxis::SetMaxDigits(3);
248  hist->SetLineWidth(1);
249  hist->SetLineColor(kRed);
250  if (!drawStat)
251  hist->SetStats(0);
252  if (histoList)
253  histoList->Add(hist);
254  return hist;
255 }
256 
257 TH2F* BKLMDigitAnalyzerModule::createTH2(const char* name, const char* title, Int_t nbinsX, Double_t minX, Double_t maxX,
258  const char* titleX, Int_t nbinsY, Double_t minY, Double_t maxY, const char* titleY, bool drawStat, TList* histoList)
259 {
260  TH2F* hist = new TH2F(name, title, nbinsX, minX, maxX, nbinsY, minY, maxY);
261  hist->GetXaxis()->SetTitle(titleX);
262  hist->GetYaxis()->SetTitle(titleY);
263  TGaxis::SetMaxDigits(3);
264  if (!drawStat)
265  hist->SetStats(0);
266  if (histoList)
267  histoList->Add(hist);
268  return hist;
269 }
Belle2::BKLMDigitAnalyzerModule::m_digitOutOfRange
StoreArray< KLMDigit > m_digitOutOfRange
Input BKLMDigits whose time is out of range.
Definition: BKLMDigitAnalyzerModule.h:106
Belle2::KLMDigitEventInfo::getSciHits
int getSciHits() const
Returns the number of scintillator hits in the event.
Definition: KLMDigitEventInfo.h:190
Belle2::BKLMDigitAnalyzerModule::m_outputRootName
std::string m_outputRootName
Output filename.
Definition: BKLMDigitAnalyzerModule.h:112
Belle2::BKLMDigitAnalyzerModule::beginRun
virtual void beginRun() override
Called when entering a new run.
Definition: BKLMDigitAnalyzerModule.cc:67
Belle2::KLMDigitRaw::getChannel
uint16_t getChannel() const
Get the channel number.
Definition: KLMDigitRaw.h:121
Belle2::BKLMDigitAnalyzerModule::m_histoLayerVsSectorPerPlane
TH2F * m_histoLayerVsSectorPerPlane[2][2]
Pointer to occupancy 2D histogram of layer vs sector for each view (=plane)
Definition: BKLMDigitAnalyzerModule.h:127
Belle2::KLMDigitEventInfo::getRPCHits
int getRPCHits() const
Returns the number of RPC hits in the event.
Definition: KLMDigitEventInfo.h:165
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::BKLMDigitAnalyzerModule::m_histoStrip
TH2F * m_histoStrip[2][8][2]
Pointer to occupancy 2D histogram of each strip.
Definition: BKLMDigitAnalyzerModule.h:136
Belle2::KLMDigitEventInfo::getOutOfRangeHits
int getOutOfRangeHits() const
Returns the number of OutOfRange-flagged hits in the event.
Definition: KLMDigitEventInfo.h:215
Belle2::BKLMDigitAnalyzerModule
Module useful to quickly analyze BKLM unpacked data.
Definition: BKLMDigitAnalyzerModule.h:54
Belle2::BKLMDigitAnalyzerModule::createTH2
TH2F * createTH2(const char *name, const char *title, Int_t nBinsX, Double_t minX, Double_t maxX, const char *titleX, Int_t nBinsY, Double_t minY, Double_t maxY, const char *titleY, bool drawStat, TList *histoList=nullptr)
Create a ROOT 2D histogram.
Definition: BKLMDigitAnalyzerModule.cc:257
Belle2::RelationsInterface::getRelatedTo
TO * getRelatedTo(const std::string &name="", const std::string &namedRelation="") const
Get the object to which this object has a relation.
Definition: RelationsObject.h:250
Belle2::BKLMDigitAnalyzerModule::~BKLMDigitAnalyzerModule
virtual ~BKLMDigitAnalyzerModule() override
Destructor.
Definition: BKLMDigitAnalyzerModule.cc:55
Belle2::BKLMDigitAnalyzerModule::initialize
virtual void initialize() override
Initializer.
Definition: BKLMDigitAnalyzerModule.cc:59
Belle2::BKLMDigitAnalyzerModule::m_digit
StoreArray< KLMDigit > m_digit
Input KLMDigits.
Definition: BKLMDigitAnalyzerModule.h:100
Belle2::KLMDigit
KLM digit (class representing a digitized hit in RPCs or scintillators).
Definition: KLMDigit.h:40
Belle2::BKLMDigitAnalyzerModule::event
virtual void event() override
This method is called for each event.
Definition: BKLMDigitAnalyzerModule.cc:180
Belle2::RelationsInterface::getRelationsTo
RelationVector< TO > getRelationsTo(const std::string &name="", const std::string &namedRelation="") const
Get the relations that point from this object to another store array.
Definition: RelationsObject.h:199
Belle2::KLMDigitEventInfo
Class to store debugging informations from the unpacker (event based).
Definition: KLMDigitEventInfo.h:34
Belle2::Module
Base class for Modules.
Definition: Module.h:74
Belle2::BKLMDigitAnalyzerModule::m_histoLayer
TH1F * m_histoLayer[2][8]
Pointer to occupancy 1D histogram of each layer for each sector.
Definition: BKLMDigitAnalyzerModule.h:130
Belle2::BKLMDigitAnalyzerModule::m_runNumber
int m_runNumber
Run number of the current data set.
Definition: BKLMDigitAnalyzerModule.h:97
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::BKLMDigitAnalyzerModule::m_histoChannel
TH2F * m_histoChannel[2][8][2]
Pointer to occupancy 2D histogram of each channel.
Definition: BKLMDigitAnalyzerModule.h:133
Belle2::StoreObjPtr
Type-safe access to single objects in the data store.
Definition: ParticleList.h:33
Belle2::BKLMDigitAnalyzerModule::createTH1
TH1F * createTH1(const char *name, const char *title, Int_t nBinsX, Double_t minX, Double_t maxX, const char *titleX, const char *titleY, bool drawStat, TList *histoList=nullptr)
Create a ROOT 1D histogram.
Definition: BKLMDigitAnalyzerModule.cc:241
Belle2::BKLMDigitAnalyzerModule::m_histoList
TList * m_histoList
Pointer to ROOT list of histograms.
Definition: BKLMDigitAnalyzerModule.h:121
Belle2::KLMDigitEventInfo::getIntTriggerCTime
int getIntTriggerCTime() const
Returns trigger CTIME as int.
Definition: KLMDigitEventInfo.h:85
Belle2::RelationsInterface::getRelationsFrom
RelationVector< FROM > getRelationsFrom(const std::string &name="", const std::string &namedRelation="") const
Get the relations that point from another store array to this object.
Definition: RelationsObject.h:214
Belle2::BKLMDigitAnalyzerModule::m_histoTdc
TH1F * m_histoTdc[2][8][2]
Pointer to TDC 1D histogram of each channel for each sector and view.
Definition: BKLMDigitAnalyzerModule.h:139
Belle2::KLMDigitRaw
Class to store the raw words from the unpacker, digit-by-digit.
Definition: KLMDigitRaw.h:39
Belle2::KLMElementNumbers::c_BKLM
@ c_BKLM
BKLM.
Definition: KLMElementNumbers.h:47
Belle2::BKLMDigitAnalyzerModule::m_outputRootFile
TFile * m_outputRootFile
Pointer to ROOT output file.
Definition: BKLMDigitAnalyzerModule.h:115
Belle2::BKLMDigitAnalyzerModule::m_digitEventInfo
StoreArray< KLMDigitEventInfo > m_digitEventInfo
Output data array of analyzed KLMDigit information.
Definition: BKLMDigitAnalyzerModule.h:109
Belle2::BKLMDigitAnalyzerModule::m_histoLayerVsSector
TH2F * m_histoLayerVsSector[2]
Pointer to occupancy 2D histograms (numerator, denominator) of layer vs sector.
Definition: BKLMDigitAnalyzerModule.h:124
Belle2::BKLMDigitAnalyzerModule::toString
std::string toString(T val)
Convert a number of type T into a string.
Definition: BKLMDigitAnalyzerModule.h:146
Belle2::BKLMDigitAnalyzerModule::endRun
virtual void endRun() override
This method is called if the current run ends.
Definition: BKLMDigitAnalyzerModule.cc:219
Belle2::BKLMDigitAnalyzerModule::m_histoCTimeDiff
TH1F * m_histoCTimeDiff[2][8][2]
Pointer to TDC-difference 1D histogram of each channel for each sector and view.
Definition: BKLMDigitAnalyzerModule.h:142
Belle2::BKLMDigitAnalyzerModule::terminate
virtual void terminate() override
This method is called at the end of the event processing.
Definition: BKLMDigitAnalyzerModule.cc:237
Belle2::BKLMDigitAnalyzerModule::m_digitRaw
StoreArray< KLMDigitRaw > m_digitRaw
Input raw BLMDigits.
Definition: BKLMDigitAnalyzerModule.h:103
Belle2::DataStore::c_Event
@ c_Event
Different object in each event, all objects/arrays are invalidated after event() function has been ca...
Definition: DataStore.h:61
Belle2::BKLMDigitAnalyzerModule::m_extraInfo
TTree * m_extraInfo
Pointer to ROOT tree with extra info.
Definition: BKLMDigitAnalyzerModule.h:118