Belle II Software  release-05-01-25
cdcDQM7.cc
1 /* Nanae Taniguchi 2017.07.12 */
2 /* Nanae Taniguchi 2018.02.06 */
3 /* add occupancy plot 2019.03 */
4 /* add occupancy hist for PV 2019.04 */
5 /* add occupancy hist for PV, updated 2019.04.28 */
6 /* update shifter's plots */
7 
8 #include "cdc/modules/cdcDQM/cdcDQM7.h"
9 
10 #include <framework/datastore/StoreArray.h>
11 #include <framework/core/HistoModule.h>
12 
13 // framework aux
14 #include <framework/gearbox/Unit.h>
15 #include <framework/gearbox/Const.h>
16 #include <framework/logging/Logger.h>
17 
18 #include <cdc/dataobjects/CDCHit.h>
19 #include <cdc/dataobjects/CDCRawHit.h>
20 
21 #include <TDirectory.h>
22 
23 using namespace std;
24 using namespace Belle2;
25 
26 //-----------------------------------------------------------------
27 // Register the Module
28 //-----------------------------------------------------------------
29 REG_MODULE(cdcDQM7)
30 
31 //-----------------------------------------------------------------
32 // Implementation
33 //-----------------------------------------------------------------
34 
36 {
37  //Set module properties
38  setDescription("CDC DQM module");
39  setPropertyFlags(c_ParallelProcessingCertified);
40 }
41 
42 cdcDQM7Module::~cdcDQM7Module()
43 {
44 }
45 
46 //------------------------------------------------------------------
47 // Function to define histograms
48 //-----------------------------------------------------------------
49 
50 TH1D* h_nhits_L[56] = {0}; // hit in each layer
51 
52 TH1D* h_tdc_sL[9] = {0}; // tdc each super layer
53 TH1D* h_adc_sL[9] = {0}; // adc each super layer
54 
55 TH1D* h_fast_tdc; // fastest TDC in each event
56 TH1D* h_board_out_tdc; // out of range TDC
57 TH2D* bmap_2; // board(copper-finess) status map 2D
58 
59 // add 20190205
60 TH1D* h_occ; // occupancy
61 TH1D* h_occ_L[56]; // occupancy layer dependent
62 TH1D* h_hit_cell; // hit vs cell id
63 
64 TH1D* h_EoccAfterInjLER; /*<nhits after LER injection>*/
65 TH1D* h_EoccAfterInjHER; /*<nhits after HER injection>*/
66 TH1D* h_occAfterInjLER; /*<occupancy after LER injection>*/
67 TH1D* h_occAfterInjHER; /*<occupancy after HER injection>*/
68 
69 
70 void cdcDQM7Module::defineHisto()
71 {
72 
73  TDirectory* oldDir = gDirectory;
74  TDirectory* dirDAQ = oldDir->mkdir("CDC");
75 
76  dirDAQ->cd();
77 
78  int bmin = 0;
79  int bmax = 0;
80  int ndiv[9] = {160, 160, 192, 224, 256, 288, 320, 352, 384};
81 
82  for (int b = 0; b < 56; b++) {
83 
84  if (b < 8) {
85  bmax = ndiv[0];
86  } else if (b >= 8 && b < 14) {
87  bmax = ndiv[1];
88  } else if (b >= 14 && b < 20) {
89  bmax = ndiv[2];
90  } else if (b >= 20 && b < 26) {
91  bmax = ndiv[3];
92  } else if (b >= 26 && b < 32) {
93  bmax = ndiv[4];
94  } else if (b >= 32 && b < 38) {
95  bmax = ndiv[5];
96  } else if (b >= 38 && b < 44) {
97  bmax = ndiv[6];
98  } else if (b >= 44 && b < 50) {
99  bmax = ndiv[7];
100  } else if (b >= 50 && b < 56) {
101  bmax = ndiv[8];
102  }
103 
104  // nhits_L
105  h_nhits_L[b] = new TH1D(Form("nhits_L%d", b), Form("nhits Layer %d", b), bmax - bmin, bmin, bmax);
106  h_nhits_L[b]->SetMinimum(0);
107  h_nhits_L[b]->SetFillColor(7);
108 
109  h_occ_L[b] = new TH1D(Form("h_occ_L%d", b), "occ. each layer", 100, 0, 1); //
110  h_occ_L[b]->SetFillColor(96);
111  h_occ_L[b]->SetMinimum(0);
112  h_occ_L[b]->SetStats(0);
113 
114  }
115 
116  for (int s = 0; s < 9; s++) {
117  h_tdc_sL[s] = new TH1D(Form("tdc_sL%d", s), Form("tdc sLayer %d", s), 250, 4200, 5200);
118  h_tdc_sL[s]->SetMinimum(0);
119  h_tdc_sL[s]->SetFillColor(6);
120 
121  h_adc_sL[s] = new TH1D(Form("adc_sL%d", s), Form("adc sLayer %d", s), 100, 0, 500);
122  h_adc_sL[s]->SetMinimum(0);
123  h_adc_sL[s]->SetFillColor(8);
124  }
125 
126  h_fast_tdc = new TH1D("fast_tdc", "fastest TDC", 50, 4800, 5000);
127  h_fast_tdc->SetFillColor(6);
128 
129  h_board_out_tdc = new TH1D("h_board_out_tdc", "out of range TDC", 300, 0, 300);
130  h_board_out_tdc->SetFillColor(95);
131 
132  // 20190205
133  h_occ = new TH1D("occ", "occ. total", 100, 0, 1.);
134  h_occ->SetFillColor(95);
135 
136  // 20191108
137  h_hit_cell = new TH1D("h_hit_cell", "Hit of each cell", 14336, 0, 14335);
138  h_hit_cell->SetFillColor(20);
139  //
140  bmap_2 = new TH2D("bmap_2", "", 75, 0, 75, 4, 0, 4);
141 
142  // LIVE
143  h_tdc_sL[6]->SetOption("LIVE");
144  h_tdc_sL[6]->SetOption("hist");
145 
146  h_board_out_tdc->SetOption("LIVE");
147  h_board_out_tdc->SetOption("hist");
148 
149  // set
150  h_fast_tdc->SetStats(1);
151  bmap_2->SetOption("zcol"); //
152  bmap_2->SetStats(0);
153 
154  h_EoccAfterInjLER = new TH1D("EoccAfterInjLer", "Eocc after LER injection", 4000, 0, 20000);
155  h_EoccAfterInjLER->SetMinimum(0);
156  h_EoccAfterInjLER->SetFillColor(7);
157 
158  h_EoccAfterInjHER = new TH1D("EoccAfterInjHer", "Eocc after HER injection", 4000, 0, 20000);
159  h_EoccAfterInjHER->SetMinimum(0);
160  h_EoccAfterInjHER->SetFillColor(7);
161 
162  h_occAfterInjLER = new TH1D("occAfterInjLer", "occupancy after LER injection", 4000, 0, 20000);
163  h_occAfterInjLER->SetMinimum(0);
164  h_occAfterInjLER->SetFillColor(7);
165 
166  h_occAfterInjHER = new TH1D("occAfterInjHer", "occupancy after HER injection", 4000, 0, 20000);
167  h_occAfterInjHER->SetMinimum(0);
168  h_occAfterInjHER->SetFillColor(7);
169 
170  oldDir->cd();//
171 
172 }
173 
174 void cdcDQM7Module::initialize()
175 {
176  REG_HISTOGRAM // required to register histograms to HistoManager
177  // register dataobjects
178  m_rawFTSW.isOptional();
179 
180 }
181 
182 void cdcDQM7Module::beginRun()
183 {
184  for (int i = 0; i < 56; i++) {
185  h_nhits_L[i]->Reset();
186  h_occ_L[i]->Reset();
187 
188  }
189 
190  for (int j = 0; j < 9; j++) {
191  h_tdc_sL[j]->Reset();
192  h_adc_sL[j]->Reset();
193  }
194 
195  h_fast_tdc->Reset();
196  h_board_out_tdc->Reset();
197  bmap_2->Reset();
198  h_occ->Reset();
199  h_hit_cell->Reset();
200 
201  h_EoccAfterInjLER->Reset();
202  h_EoccAfterInjHER->Reset();
203  h_occAfterInjLER->Reset();
204  h_occAfterInjHER->Reset();
205 
206 }
207 
208 void cdcDQM7Module::event()
209 {
210 
211  StoreArray<CDCHit> cdcHits;
212  int nent = cdcHits.getEntries();
213  int ftdc = 0;
214 
215  // occ total
216  double occ_total = nent / 14336.;
217  h_occ->Fill(occ_total);
218 
219  // for layer dependent occupancy
220  int whits_L[56] = {}; // wire hits
221  double occ_L[56] = {}; // occupancy
222 
223  int ndiv[9] = {160, 160, 192, 224, 256, 288, 320, 352, 384};
224 
225  for (int i = 0; i < nent; i++) {
226  CDCHit* cdchit = static_cast<CDCHit*>(cdcHits[i]);
227 
228  int sL = cdchit->getISuperLayer();
229  int iL = cdchit->getILayer();
230  int wid = cdchit->getIWire();
231  int adcsum = cdchit->getADCCount();
232  int vtdc = cdchit->getTDCCount();
233 
234  if (sL > 8) continue; // error
235  if (iL > 8) continue; // error
236 
237  int num = sL * 6 + iL + 2;
238  if (num > 55) continue; // error
239 
240  // wire hits
241  if (sL == 0) {
242  whits_L[iL]++;
243  } else {
244  whits_L[num]++;
245  }
246 
247  if (adcsum > -1) {
248 
249  if (sL == 0) {
250  h_nhits_L[iL]->Fill(wid);
251  } else {
252  h_nhits_L[num]->Fill(wid);
253  }
254 
255  //
256  if (vtdc > ftdc && adcsum > 20) {
257  ftdc = vtdc;
258  }// fastest
259 
260  }// adc
261 
262  // add by J.H. Yin
263  if (adcsum > 25) {
264  int cid(0);
265  if (sL == 0) {
266  cid = iL * ndiv[sL] + wid;
267  } else {
268  for (int isl = 0; isl < sL; isl ++) {
269  cid += 6 * ndiv[isl];
270  }
271  cid += 2 * ndiv[0] + iL * ndiv[sL] + wid;
272  }
273  h_hit_cell -> Fill(cid);
274  }
275  }// cdchit
276 
277  h_fast_tdc->Fill(ftdc);
278 
279  // each layer
280  for (int b = 0; b < 56; b++) {
281  int n_wire = 0;
282  if (b < 8) {
283  n_wire = ndiv[0];
284  } else if (b >= 8 && b < 14) {
285  n_wire = ndiv[1];
286  } else if (b >= 14 && b < 20) {
287  n_wire = ndiv[2];
288  } else if (b >= 20 && b < 26) {
289  n_wire = ndiv[3];
290  } else if (b >= 26 && b < 32) {
291  n_wire = ndiv[4];
292  } else if (b >= 32 && b < 38) {
293  n_wire = ndiv[5];
294  } else if (b >= 38 && b < 44) {
295  n_wire = ndiv[6];
296  } else if (b >= 44 && b < 50) {
297  n_wire = ndiv[7];
298  } else if (b >= 50 && b < 56) {
299  n_wire = ndiv[8];
300  }
301 
302  // cal. occupancy
303  occ_L[b] = (double)whits_L[b] / n_wire;
304  h_occ_L[b]->Fill(occ_L[b]);
305  }
306 
307 
308  //
309  StoreArray<CDCRawHit> cdcRawHits;
310  int r_nent = cdcRawHits.getEntries();
311 
312  // new
313  for (int j = 0; j < r_nent; j++) {
314  CDCRawHit* cdcrawhit = static_cast<CDCRawHit*>(cdcRawHits[j]);
315 
316  int brd = cdcrawhit->getBoardId();
317  int v_adc = cdcrawhit->getFADC();
318  int v_tdc = cdcrawhit->getTDC();
319  int n_tot = cdcrawhit->getTOT();
320  int n_node = cdcrawhit->getNode();
321  int n_fns = cdcrawhit->getFiness();
322 
323  bmap_2->Fill(n_node, n_fns);
324 
325  if (v_tdc > 5200 || v_tdc < 4200) {
326  h_board_out_tdc->Fill(brd);
327  }
328 
329  // printf("%d, %d:\n", j, board);
330  if (brd > 299) continue;
331  if (n_tot < 4) continue;
332  if (v_adc < 35) continue;
333 
334  // each sL
335  if (brd < 28) {
336  h_tdc_sL[0]->Fill(v_tdc);
337  h_adc_sL[0]->Fill(v_adc);
338  } else if (brd > 27 && brd < 48) {
339  h_tdc_sL[1]->Fill(v_tdc);
340  h_adc_sL[1]->Fill(v_adc);
341  } else if (brd > 47 && brd < 72) {
342  h_tdc_sL[2]->Fill(v_tdc);
343  h_adc_sL[2]->Fill(v_adc);
344  } else if (brd > 71 && brd < 100) {
345  h_tdc_sL[3]->Fill(v_tdc);
346  h_adc_sL[3]->Fill(v_adc);
347  } else if (brd > 99 && brd < 132) {
348  h_tdc_sL[4]->Fill(v_tdc);
349  h_adc_sL[4]->Fill(v_adc);
350  } else if (brd > 131 && brd < 168) {
351  h_tdc_sL[5]->Fill(v_tdc);
352  h_adc_sL[5]->Fill(v_adc);
353  } else if (brd > 167 && brd < 208) {
354  h_tdc_sL[6]->Fill(v_tdc);
355  h_adc_sL[6]->Fill(v_adc);
356  } else if (brd > 207 && brd < 252) {
357  h_tdc_sL[7]->Fill(v_tdc);
358  h_adc_sL[7]->Fill(v_adc);
359  } else if (brd > 251) {
360  h_tdc_sL[8]->Fill(v_tdc);
361  h_adc_sL[8]->Fill(v_adc);
362  }
363 
364 
365  }// cdcrawhits
366 
367  for (auto& it : m_rawFTSW) {
368  B2DEBUG(29, "TTD FTSW : " << hex << it.GetTTUtime(0) << " " << it.GetTTCtime(0) << " EvtNr " << it.GetEveNo(0) << " Type " <<
369  (it.GetTTCtimeTRGType(0) & 0xF) << " TimeSincePrev " << it.GetTimeSincePrevTrigger(0) << " TimeSinceInj " <<
370  it.GetTimeSinceLastInjection(0) << " IsHER " << it.GetIsHER(0) << " Bunch " << it.GetBunchNumber(0));
371  auto difference = it.GetTimeSinceLastInjection(0);
372  if (difference != 0x7FFFFFFF) {
373  //unsigned int nentries = m_digits.getEntries();
374  float diff2 = difference / 127.; // 127MHz clock ticks to us, inexact rounding
375  if (it.GetIsHER(0)) {
376  h_occAfterInjHER->Fill(diff2, occ_total);
377  h_EoccAfterInjHER->Fill(diff2);
378  } else {
379  h_occAfterInjLER->Fill(diff2, occ_total);
380  h_EoccAfterInjLER->Fill(diff2);
381  }
382  }
383  }
384 
385 
386 
387 
388 }
389 
390 
391 void cdcDQM7Module::endRun()
392 {
393  //
394 }
395 
396 
397 void cdcDQM7Module::terminate()
398 {
399  //
400 }
Belle2::CDCHit::getISuperLayer
unsigned short getISuperLayer() const
Getter for iSuperLayer.
Definition: CDCHit.h:195
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::CDCHit
Class containing the result of the unpacker in raw data and the result of the digitizer in simulation...
Definition: CDCHit.h:51
Belle2::CDCHit::getTDCCount
short getTDCCount() const
Getter for TDC count.
Definition: CDCHit.h:230
Belle2::CDCRawHit::getFADC
unsigned short getFADC(void) const
Getter for FADC value.
Definition: CDCRawHit.h:111
Belle2::CDCHit::getADCCount
unsigned short getADCCount() const
Getter for integrated charge.
Definition: CDCHit.h:241
Belle2::CDCRawHit::getBoardId
unsigned short getBoardId(void) const
Getter for boar ID.
Definition: CDCRawHit.h:103
Belle2::cdcDQM7Module
The module for Data Quality Monitor.
Definition: cdcDQM7.h:28
Belle2::CDCHit::getILayer
unsigned short getILayer() const
Getter for iLayer.
Definition: CDCHit.h:183
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::CDCHit::getIWire
unsigned short getIWire() const
Getter for iWire.
Definition: CDCHit.h:177
Belle2::CDCRawHit::getTOT
unsigned short getTOT(void) const
Getter for TOT value.
Definition: CDCRawHit.h:127
Belle2::CDCRawHit::getNode
unsigned short getNode(void) const
Getter for Node ID.
Definition: CDCRawHit.h:70
Belle2::StoreArray
Accessor to arrays stored in the data store.
Definition: ECLMatchingPerformanceExpertModule.h:33
Belle2::CDCRawHit
The CDCRawHit (suppressed mode) class.
Definition: CDCRawHit.h:36
Belle2::CDCRawHit::getFiness
unsigned short getFiness(void) const
Getter fot Finess ID.
Definition: CDCRawHit.h:79
Belle2::CDCRawHit::getTDC
unsigned short getTDC(void) const
Getter for TDC value.
Definition: CDCRawHit.h:119
Belle2::StoreArray::getEntries
int getEntries() const
Get the number of objects in the array.
Definition: StoreArray.h:226
Belle2::HistoModule
HistoModule.h is supposed to be used instead of Module.h for the modules with histogram definitions t...
Definition: HistoModule.h:29