Belle II Software  release-05-01-25
ARICHHitRateCounter.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2019 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Marko Staric, Luka Santelj *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 // Own include
12 #include <background/modules/BeamBkgHitRateMonitor/ARICHHitRateCounter.h>
13 
14 // framework aux
15 #include <framework/logging/Logger.h>
16 
17 using namespace std;
18 
19 namespace Belle2 {
24  namespace Background {
25 
26  void ARICHHitRateCounter::initialize(TTree* tree)
27  {
28  // register collection(s) as optional, your detector might be excluded in DAQ
29  m_hits.isOptional();
30 
31  // set branch address
32  tree->Branch("arich", &m_rates, "segmentRates[18]/F:averageRate/F:numEvents/I:valid/O");
33 
34  // make map of modules to 18 segments
35  // (first hapd ring is 6 segments (by sector), for the rest, each segment merges 3 hapd rings (again by sector))
36  int nModInRing[7] = {0, 42, 90, 144, 204, 270, 342};
37  int iRing = 0;
38  for (int i = 0; i < 420; i++) {
39  if (i == nModInRing[iRing + 1]) iRing++;
40  int segment = (i - nModInRing[iRing]) / (7 + iRing);
41  if (iRing > 0) segment += 6;
42  if (iRing > 3) segment += 6;
43  m_segmentMap[i] = segment;
44  }
45 
46  // set fractions of active channels
47  setActiveHapds();
48  }
49 
50  void ARICHHitRateCounter::clear()
51  {
52  m_buffer.clear();
53  }
54 
55  void ARICHHitRateCounter::accumulate(unsigned timeStamp)
56  {
57  // check if data are available
58  if (not m_hits.isValid()) return;
59 
60  // get buffer element
61  auto& rates = m_buffer[timeStamp];
62 
63  // increment event counter
64  rates.numEvents++;
65 
66  // count and weight hits accoring to channel efficiecny
67  for (const auto& hit : m_hits) {
68  if (hit.getModule() < 1 || hit.getModule() > 420) continue;
69  auto effi = m_modulesInfo->getChannelQE(hit.getModule(), hit.getChannel());
70  float wt = std::min(1.0 / effi, 100.);
71  rates.segmentRates[m_segmentMap[hit.getModule() - 1]] += wt;
72  rates.averageRate += wt;
73  }
74 
75  // set flag to true to indicate the rates are valid
76  rates.valid = true;
77 
78  }
79 
80  void ARICHHitRateCounter::normalize(unsigned timeStamp)
81  {
82  // copy buffer element
83  m_rates = m_buffer[timeStamp];
84 
85  if (not m_rates.valid) return;
86 
87  // normalize
88  m_rates.normalize();
89 
90  // correct rates for masked-out channels
91  if (m_channelMask.hasChanged()) setActiveHapds();
92 
93  for (int iSegment = 0; iSegment < 18; iSegment++) {
94  double nHapds = m_activeHapds[iSegment];
95  if (nHapds > 0) m_rates.segmentRates[iSegment] /= nHapds;
96  else m_rates.segmentRates[iSegment] = 0;
97  }
98  m_rates.averageRate /= m_activeTotal;
99  }
100 
101  void ARICHHitRateCounter::setActiveHapds()
102  {
103  for (auto& nactive : m_activeHapds) nactive = 0;
104 
105  if (not m_channelMask.isValid()) {
106  for (unsigned imod = 1; imod < 421; imod++) m_activeHapds[m_segmentMap[imod - 1]] += 1.;
107  m_activeTotal = 420;
108  B2WARNING("ARICHHitRateCounter: no valid channel mask - all HAPDs set to active");
109  return;
110  }
111 
112  int nactiveTotal = 0;
113  for (unsigned imod = 1; imod < 421; imod++) {
114  int nactive = 0;
115  for (unsigned ichn = 0; ichn < 144; ichn++) {
116  if (m_channelMask->isActive(imod, ichn)) nactive++;
117  }
118  nactiveTotal += nactive;
119  m_activeHapds[m_segmentMap[imod - 1]] += (float)nactive / 144.;
120  }
121  m_activeTotal = (float)nactiveTotal / 144.;
122  }
123 
124  } // Background namespace
126 } // Belle2 namespace
127 
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19