Belle II Software  release-05-01-25
eclDQMInjection.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2019 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Bjoern Spruck, Dmitry Matvienko *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 //This module
12 #include <ecl/modules/eclDQMInjection/eclDQMInjection.h>
13 
14 //Boost
15 #include <boost/format.hpp>
16 #include <boost/range/combine.hpp>
17 
18 //ECL
19 #include <ecl/dataobjects/ECLDsp.h>
20 #include <ecl/utility/ECLDspUtilities.h>
21 
22 //ROOT
23 #include "TDirectory.h"
24 
25 using namespace std;
26 using namespace Belle2;
27 using namespace Belle2::ECL;
28 
29 //-----------------------------------------------------------------
30 // Register the Module
31 //-----------------------------------------------------------------
32 REG_MODULE(ECLDQMInjection)
33 
34 //-----------------------------------------------------------------
35 // Implementation
36 //-----------------------------------------------------------------
37 
39  : HistoModule(),
40  m_calibrationThrApsd("ECL_FPGA_StoreWaveform")
41 {
42  //Set module properties
43  setDescription("Monitor Occupancy after Injection");
44  setPropertyFlags(c_ParallelProcessingCertified);
45  addParam("histogramDirectoryName", m_histogramDirectoryName, "Name of the directory where histograms will be placed",
46  std::string("ECLINJ"));
47  addParam("ECLDigitsName", m_ECLDigitsName, "Name of ECL hits", std::string(""));
48  // BeamRevolutionCycle is set based on 'Timing distribution for the Belle II
49  // data acquistion system'. RF clock of 508 MHz is synchronized to
50  // beam-revolution cycle (5120 RF bunches in one cycle).
51  addParam("BeamRevolutionCycle", m_revolutionTime, "Beam revolution cycle in musec", 5120 / 508.);
52  addParam("ECLThresholdforVetoTuning", m_ECLThresholdforVetoTuning, "ECL Threshold for injection veto tuning, ADC channels", 400.);
53  addParam("DPHYTTYP", m_DPHYTTYP,
54  "Flag to control trigger of delayed bhabha events; 0 - select events by 'bha_delay' trigger bit, 1 - select by TTYP_DPHY", false);
55 
56 }
57 
58 void ECLDQMInjectionModule::defineHisto()
59 {
60  TDirectory* oldDir = gDirectory;
61  oldDir->mkdir(m_histogramDirectoryName.c_str());// do not rely on return value, might be ZERO
62  oldDir->cd(m_histogramDirectoryName.c_str());//changing to the right directory
63 
64  hHitsAfterInjLER = new TH1F("ECLHitsInjLER", "ECLHitsInjLER/Time;Time in #mus;Count/Time (5 #mus bins)", 4000, 0, 20000);
65  hHitsAfterInjHER = new TH1F("ECLHitsInjHER", "ECLHitsInjHER/Time;Time in #mus;Count/Time (5 #mus bins)", 4000, 0, 20000);
66  hEHitsAfterInjLER = new TH1F("ECLEHitsInjLER", "ECLEHitsInjLER/Time;Time in #mus;Triggers/Time (5 #mus bins)", 4000, 0, 20000);
67  hEHitsAfterInjHER = new TH1F("ECLEHitsInjHER", "ECLEHitsInjHER/Time;Time in #mus;Triggers/Time (5 #mus bins)", 4000, 0, 20000);
68  hBurstsAfterInjLER = new TH1F("ECLBurstsInjLER", "ECLBurstsInjLER/Time;Time in #mus;Count/Time (1 #mus bins)", 20000, 0, 20000);
69  hBurstsAfterInjHER = new TH1F("ECLBurstsInjHER", "ECLBurstsInjHER/Time;Time in #mus;Count/Time (1 #mus bins)", 20000, 0, 20000);
70  hEBurstsAfterInjLER = new TH1F("ECLEBurstsInjLER", "ECLEBurstsInjLER/Time;Time in #mus;Triggers/Time (1 #mus bins)", 20000, 0,
71  20000);
72  hEBurstsAfterInjHER = new TH1F("ECLEBurstsInjHER", "ECLEBurstsInjHER/Time;Time in #mus;Triggers/Time (1 #mus bins)", 20000, 0,
73  20000);
74  hVetoAfterInjLER = new TH2F("ECLVetoAfterInjLER",
75  "ECL Hits for LER veto tuning (E > 20 MeV);"
76  "Time since last injection in #mus;Time within beam cycle in #mus", 500, 0, 30000, 100, 0,
77  m_revolutionTime);
78  hVetoAfterInjHER = new TH2F("ECLVetoAfterInjHER",
79  "ECL Hits for HER veto tuning (E > 20 MeV);"
80  "Time since last injection in #mus;Time within beam cycle in #mus", 500, 0, 30000, 100, 0,
81  m_revolutionTime);
82  hOccAfterInjLER = new TH2F("ECLOccAfterInjLER",
83  "ECL Occupancy after LER injection (E > 1 MeV); Time since last injection in #mus;Occupancy (Nhits/8736) [%]",
84  100, 0, 20000, 98, 2, 100);
85  hOccAfterInjHER = new TH2F("ECLOccAfterInjHER",
86  "ECL Occupancy after HER injection (E > 1 MeV); Time since last injection in #mus;Occupancy (Nhits/8736) [%]",
87  100, 0, 20000, 98, 2, 100);
88 
89 
90  //== Fill h_ped_peak vector
91 
92  m_ped_peak_range = {
93  1.0, 1.5, 2.0, 4, 6, 8, 10
94  };
95  int ped_peak_range_count = m_ped_peak_range.size() - 1;
96 
97  static const std::string part_names[] = {"fwd", "bar", "bwd"};
98  static const std::string title_suffix[] = {
99  "in fwd endcap", "in barrel", "in bwd endcap"
100  };
101 
102  for (int ler_her = 0; ler_her < 2; ler_her++) {
103  std::string ring_name = (ler_her == 0) ? "LER" : "HER";
104  for (int part = 0; part < 3; part++) {
105  std::string suffix = title_suffix[part];
106  for (int i = 0; i < ped_peak_range_count; i++) {
107  float min_time = m_ped_peak_range[i];
108  float max_time = m_ped_peak_range[i + 1];
109  std::string name, title;
110  name = str(boost::format("ped_peak_%s_%s_%d") %
111  ring_name % part_names[part] % i);
112  title = str(boost::format("Peak height %.1f-%.1f ms after %s inj %s") %
113  min_time % max_time % ring_name % suffix);
114 
115  auto h = new TH1F(name.c_str(), title.c_str(), 300, 0.0, 0.3);
116  h->GetXaxis()->SetTitle("Peak height in first 16 points [GeV]");
117 
118  h_ped_peak.push_back(h);
119  }
120  }
121  }
122 
123  // Initialize coefficients used by pedestalFit function
124  ECLDspUtilities::initPedestalFit();
125 
126  // cd back to root directory
127  oldDir->cd();
128 }
129 
130 void ECLDQMInjectionModule::initialize()
131 {
132  REG_HISTOGRAM
133  m_rawTTD.isOptional();
134  m_storeHits.isRequired(m_ECLDigitsName);
135  m_ECLTrigs.isOptional();
136  m_ECLDsps.isOptional();
137  m_l1Trigger.isOptional();
138 
139  if (!mapper.initFromDB()) B2FATAL("ECL Display:: Can't initialize eclChannelMapper");
140 
141  v_totalthrApsd.resize((m_calibrationThrApsd->getCalibVector()).size());
142  for (size_t i = 0; i < v_totalthrApsd.size(); i++) v_totalthrApsd[i] = (int)(m_calibrationThrApsd->getCalibVector())[i];
143 }
144 
145 void ECLDQMInjectionModule::beginRun()
146 {
147  // Assume that everthing is non-yero ;-)
148  hHitsAfterInjLER->Reset();
149  hHitsAfterInjHER->Reset();
150  hEHitsAfterInjLER->Reset();
151  hEHitsAfterInjHER->Reset();
152  hBurstsAfterInjLER->Reset();
153  hBurstsAfterInjHER->Reset();
154  hEBurstsAfterInjLER->Reset();
155  hEBurstsAfterInjHER->Reset();
156  hVetoAfterInjLER->Reset();
157  hVetoAfterInjHER->Reset();
158  hOccAfterInjHER->Reset();
159  hOccAfterInjLER->Reset();
160 }
161 
162 void ECLDQMInjectionModule::event()
163 {
164  bool bhatrig = false;
165 
166  if (m_l1Trigger.isValid() && m_DPHYTTYP) bhatrig = m_l1Trigger->getTimType() == TRGSummary::ETimingType::TTYP_DPHY;
167  else if (m_l1Trigger.isValid() && !m_DPHYTTYP) bhatrig = m_l1Trigger->testInput("bha_delay");
168 
169  if (m_eventmetadata.isValid() && m_eventmetadata->getErrorFlag() != 0x10) {
170  m_iEvent = m_eventmetadata->getEvent();
171  } else m_iEvent = -1;
172  int discarded_wfs = 0;
173  for (auto& aECLTrig : m_ECLTrigs) {
174  int crate = aECLTrig.getTrigId();
175  int suppress = aECLTrig.getBurstSuppressionMask();
176  int shaper_pos = 0;
177  while (suppress) {
178  shaper_pos ++;
179  bool shaper_bit = suppress & 1;
180  if (shaper_bit) {
181  if (m_iEvent % 1000 == 999 || (m_l1Trigger.isValid() && m_l1Trigger->getTimType() == TRGSummary::ETimingType::TTYP_RAND) ||
182  (m_l1Trigger.isValid() && bhatrig)) {
183  for (int channel_pos = 0; channel_pos < 16; channel_pos ++) {
184  if (mapper.getCellId(crate, shaper_pos, channel_pos) > 0) discarded_wfs += 1;
185  }
186  } else {
187  for (auto& aECLDigit : m_storeHits) {
188  if (crate == mapper.getCrateID(aECLDigit.getCellId()) && shaper_pos == mapper.getShaperPosition(aECLDigit.getCellId()) &&
189  aECLDigit.getAmp() >= (v_totalthrApsd[aECLDigit.getCellId() - 1] / 4 * 4)) discarded_wfs += 1;
190  }
191  }
192  }
193  suppress >>= 1;
194  }
195  }
196 
197  unsigned int ECLDigitsAboveThr = 0; // Threshold is set to 20 MeV
198  unsigned int ECLDigitsAboveThr1MeV = 0;
199  for (auto& aECLDigit : m_storeHits) {
200  if (aECLDigit.getAmp() > m_ECLThresholdforVetoTuning) ECLDigitsAboveThr += 1;
201  if (aECLDigit.getAmp() > 20) ECLDigitsAboveThr1MeV += 1;
202  }
203 
204  for (auto& it : m_rawTTD) {
205  B2DEBUG(29, "TTD FTSW : " << hex << it.GetTTUtime(0) << " " << it.GetTTCtime(0) << " EvtNr " << it.GetEveNo(0) << " Type " <<
206  (it.GetTTCtimeTRGType(0) & 0xF) << " TimeSincePrev " << it.GetTimeSincePrevTrigger(0) << " TimeSinceInj " <<
207  it.GetTimeSinceLastInjection(0) << " IsHER " << it.GetIsHER(0) << " Bunch " << it.GetBunchNumber(0));
208 
209  // get last injection time
210  auto difference = it.GetTimeSinceLastInjection(0);
211  // check time overflow, too long ago
212  if (difference != 0x7FFFFFFF) {
213  unsigned int all = m_storeHits.getEntries();
214  float diff2 = difference / 127.; // 127MHz clock ticks to us, inexact rounding
215  int is_her = it.GetIsHER(0);
216  if (is_her) {
217  hHitsAfterInjHER->Fill(diff2, all);
218  hEHitsAfterInjHER->Fill(diff2);
219  hBurstsAfterInjHER->Fill(diff2, discarded_wfs);
220  hEBurstsAfterInjHER->Fill(diff2);
221  hVetoAfterInjHER->Fill(diff2, diff2 - int(diff2 / m_revolutionTime)*m_revolutionTime, ECLDigitsAboveThr);
222  if (all > 0) hOccAfterInjHER->Fill(diff2, ECLDigitsAboveThr1MeV / 8736.*100.);
223  } else {
224  hHitsAfterInjLER->Fill(diff2, all);
225  hEHitsAfterInjLER->Fill(diff2);
226  hBurstsAfterInjLER->Fill(diff2, discarded_wfs);
227  hEBurstsAfterInjLER->Fill(diff2);
228  hVetoAfterInjLER->Fill(diff2, diff2 - int(diff2 / m_revolutionTime)*m_revolutionTime, ECLDigitsAboveThr);
229  if (all > 0) hOccAfterInjLER->Fill(diff2, ECLDigitsAboveThr1MeV / 8736.*100.);
230  }
231 
232  //== Filling h_ped_peak histograms
233  int range_count = m_ped_peak_range.size() - 1;
234  if (diff2 < m_ped_peak_range[range_count] * 1000) {
235  //== Identify which histogram to fill (according to inj time range)
236  int range_id;
237  for (range_id = 0; range_id < range_count; range_id++) {
238  // Converting from ms to us
239  float min_time = m_ped_peak_range[range_id ] * 1000;
240  float max_time = m_ped_peak_range[range_id + 1] * 1000;
241  if (diff2 > min_time && diff2 < max_time) break;
242  }
243  //== Find pedestal peaks in all available waveforms
244  if (range_id < range_count) {
245  for (auto& aECLDsp : m_ECLDsps) {
246  auto result = ECLDspUtilities::pedestalFit(aECLDsp.getDspA());
247 
248  //== Identify which histogram to fill (HER/LER,FWD/BAR/BWD)
249  int cid = aECLDsp.getCellId();
250  int part_id = 0; // forward endcap
251  if (cid >= 1153) part_id = 1; // barrel
252  if (cid >= 7777) part_id = 2; // backward endcap
253 
254  int hist_id = is_her * 3 * range_count + part_id * range_count + range_id;
255  // NOTE: We are using the approximate conversion to energy here.
256  // (20'000 ADC counts ~= 1 GeV)
257  h_ped_peak[hist_id]->Fill(result.amp / 2e4);
258  }
259  }
260  }
261  }
262 
263  break;
264  }
265 }
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::ECL::ECLDQMInjectionModule
The ECL Occ after Injection DQM module.
Definition: eclDQMInjection.h:53
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::HistoModule
HistoModule.h is supposed to be used instead of Module.h for the modules with histogram definitions t...
Definition: HistoModule.h:29