Belle II Software  release-05-01-25
PXDRawDQMModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2013 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Klemens Lautenbach, Bjoern Spruck *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include <pxd/modules/pxdDQM/PXDRawDQMModule.h>
12 #include <vxd/geometry/GeoCache.h>
13 
14 #include <TDirectory.h>
15 #include <boost/format.hpp>
16 #include <string>
17 
18 using namespace std;
19 using namespace Belle2;
20 using namespace Belle2::PXD;
21 using namespace Belle2::VXD;
22 
23 using boost::format;
24 
25 //-----------------------------------------------------------------
26 // Register the Module
27 //-----------------------------------------------------------------
28 REG_MODULE(PXDRawDQM)
29 
30 //-----------------------------------------------------------------
31 // Implementation
32 //-----------------------------------------------------------------
33 
34 PXDRawDQMModule::PXDRawDQMModule() : HistoModule() , m_storeRawPxdrarray() , m_storeRawHits(), m_storeRawAdcs()
35 {
36  //Set module properties
37  setDescription("Monitor raw PXD");
38  setPropertyFlags(c_ParallelProcessingCertified);
39  addParam("histogramDirectoryName", m_histogramDirectoryName, "Name of the directory where histograms will be placed",
40  std::string("pxdraw"));
41 
42  addParam("RawPXDsName", m_storeRawPxdrarrayName, "The name of the StoreArray of RawPXDs to be processed", string(""));
43  addParam("PXDRawHitsName", m_storeRawHitsName, "The name of the StoreArray of PXDRawHits to be processed", string(""));
44  addParam("PXDRawAdcsName", m_storeRawAdcsName, "The name of the StoreArray of PXDRawAdcs to be processed", string(""));
45 }
46 
47 void PXDRawDQMModule::defineHisto()
48 {
49  // Create a separate histogram directory and cd into it.
50  TDirectory* oldDir = gDirectory;
51  if (m_histogramDirectoryName != "") {
52  oldDir->mkdir(m_histogramDirectoryName.c_str());// do not use return value with ->cd(), its ZERO if dir already exists
53  oldDir->cd(m_histogramDirectoryName.c_str());
54  }
55 
56  hrawPxdPackets = new TH1F("hrawPxdPackets", "Pxd Raw Packet Nr;Nr per Event", 16, 0, 16);
57  hrawPxdPacketSize = new TH1F("hrawPxdPacketSize", "Pxd Raw Packetsize;Words per packet", 1024, 0, 1024);
58  hrawPxdHitMapAll = new TH2F("hrawPxdHitMapAll",
59  "Pxd Raw Hit Map Overview;column+(ladder-1)*300+100;row+850*((layer-1)*2+(sensor-1))", 3700 / 50, 0, 3700, 3500 / 50, 0, 3500);
60  // ADC map not supported by DHC anymore ... deactive filling, later remove
61  hrawPxdAdcMapAll =
62  NULL;// new TH2F("hrawPxdAdcMapAll", "Pxd Raw Adc Map Overview;column+(ladder-1)*300+100;row+850*((layer-1)*2+(sensor-1))", 370/*0*/, 0, 3700, 350/*0*/, 0, 3500);
63 
64  hrawPxdHitsCount = new TH1F("hrawPxdCount", "Pxd Raw Count ;Nr per Event", 8192, 0, 8192);
65  for (auto i = 0; i < 64; i++) {
66  auto layer = (((i >> 5) & 0x1) + 1);
67  auto ladder = ((i >> 1) & 0xF);
68  auto sensor = ((i & 0x1) + 1);
69 
70  // Check if sensor exist
71  if (Belle2::VXD::GeoCache::getInstance().validSensorID(Belle2::VxdID(layer, ladder, sensor))) {
72  string s = str(format("Sensor %d:%d:%d (DHH ID %02Xh)") % layer % ladder % sensor % i);
73  string s2 = str(format("_%d.%d.%d") % layer % ladder % sensor);
74 
75  hrawPxdHitMap[i] = new TH2F(("hrawPxdHitMap" + s2).c_str(),
76  ("Pxd Raw Hit Map, " + s + ";column;row").c_str(), 250,
77  0, 250, 768, 0, 768);
78  hrawPxdChargeMap[i] = new TH2F(("hrawPxdChargeMap" + s2).c_str(),
79  ("Pxd Raw Charge Map, " + s + ";column;row").c_str(), 250, 0, 250, 768, 0, 768);
80  hrawPxdHitsCharge[i] = new TH1F(("hrawPxdHitsCharge" + s2).c_str(),
81  ("Pxd Raw Hit Charge, " + s + ";Charge").c_str(), 256, 0, 256);
82  hrawPxdHitTimeWindow[i] = new TH1F(("hrawPxdHitTimeWindow" + s2).c_str(),
83  ("Pxd Raw Hit Time Window (framenr*192-gate_of_hit), " + s + ";Time [a.u.]").c_str(), 2048, -256, 2048 - 256);
84  hrawPxdGateTimeWindow[i] = new TH1F(("hrawPxdGateTimeWindow" + s2).c_str(),
85  ("Pxd Raw Gate Time Window (framenr*192-triggergate_of_hit), " + s + ";Time [a.u.]").c_str(), 2048, -256, 2048 - 256);
86  } else {
87  hrawPxdHitMap[i] = NULL;
88  hrawPxdChargeMap[i] = NULL;
89  hrawPxdHitsCharge[i] = NULL;
90  hrawPxdHitTimeWindow[i] = NULL;
91  hrawPxdGateTimeWindow[i] = NULL;
92  }
93  }
94 
95  // cd back to root directory
96  oldDir->cd();
97 }
98 
99 void PXDRawDQMModule::initialize()
100 {
101  REG_HISTOGRAM
102  m_storeRawPxdrarray.isOptional(m_storeRawPxdrarrayName);
103  m_storeRawHits.isRequired(m_storeRawHitsName);
104  m_storeRawAdcs.isRequired(m_storeRawAdcsName);
105  m_storeDAQEvtStats.isRequired();
106 }
107 
108 void PXDRawDQMModule::beginRun()
109 {
110  // Just to make sure, reset all the histograms.
111  if (hrawPxdPackets) hrawPxdPackets->Reset();
112  if (hrawPxdPacketSize) hrawPxdPacketSize->Reset();
113  if (hrawPxdHitsCount) hrawPxdHitsCount->Reset();
114  if (hrawPxdHitMapAll) hrawPxdHitMapAll->Reset();
115  if (hrawPxdAdcMapAll) hrawPxdAdcMapAll->Reset();
116  for (int i = 0; i < 64; i++) {
117  if (hrawPxdHitMap[i]) hrawPxdHitMap[i]->Reset();
118  if (hrawPxdChargeMap[i]) hrawPxdChargeMap[i]->Reset();
119  if (hrawPxdHitsCharge[i]) hrawPxdHitsCharge[i]->Reset();
120  if (hrawPxdHitTimeWindow[i]) hrawPxdHitTimeWindow[i]->Reset();
121  if (hrawPxdGateTimeWindow[i]) hrawPxdGateTimeWindow[i]->Reset();
122  }
123 }
124 
125 void PXDRawDQMModule::event()
126 {
127  hrawPxdPackets->Fill(m_storeRawPxdrarray.getEntries());
128 
129  for (auto& it : m_storeRawPxdrarray) {
130  if (hrawPxdPacketSize) hrawPxdPacketSize->Fill(it.size());
131  }
132 
133  if (hrawPxdHitsCount) hrawPxdHitsCount->Fill(m_storeRawHits.getEntries());
134 
135  for (auto& it : m_storeRawHits) {
136 
137  VxdID currentVxdId = it.getSensorID();
138  auto layer = currentVxdId.getLayerNumber();
139  auto ladder = currentVxdId.getLadderNumber();
140  auto sensor = currentVxdId.getSensorNumber();
141 
142  // Get startrow and DheID from DAQEvtStats
143  const PXDDAQDHEStatus* dhe = (*m_storeDAQEvtStats).findDHE(currentVxdId);
144  if (dhe == nullptr) {
145  B2ERROR("No DHE found for SensorId: " << currentVxdId);
146  continue;
147  }
148 
149  auto dhh_id = dhe->getDHEID();
150  auto startGate = dhe->getTriggerGate();
151 
152  if (dhh_id == 0 || dhh_id >= 64) {
153  B2ERROR("SensorId (DHH ID) out of range: " << dhh_id);
154  continue;
155  }
156  if (hrawPxdHitMap[dhh_id]) hrawPxdHitMap[dhh_id]->Fill(it.getColumn(), it.getRow());
157  if (hrawPxdHitMapAll) hrawPxdHitMapAll->Fill(it.getColumn() + ladder * 300 - 200,
158  100 + it.getRow() + 850 * (layer + layer + sensor - 3));
159  if (hrawPxdChargeMap[dhh_id]) hrawPxdChargeMap[dhh_id]->Fill(it.getColumn(), it.getRow(), it.getCharge());
160  if (hrawPxdHitsCharge[dhh_id]) hrawPxdHitsCharge[dhh_id]->Fill(it.getCharge());
161  // Is this histogram necessary? we are folding with occupancy of sensor hits here
162  // Think about 1024*framenr-hit_row?
163  if (hrawPxdHitTimeWindow[dhh_id]) hrawPxdHitTimeWindow[dhh_id]->Fill(it.getFrameNr() * 192 - it.getRow() / 4);
164  if (hrawPxdGateTimeWindow[dhh_id]) hrawPxdGateTimeWindow[dhh_id]->Fill(it.getFrameNr() * 192 - startGate);
165  }
166 
167  if (hrawPxdAdcMapAll) {
168  for (auto& it : m_storeRawAdcs) {
169  int dhh_id;
170  // calculate DHH id from Vxd Id
171  unsigned int layer, ladder, sensor;//, segment;
172  VxdID currentVxdId;
173  currentVxdId = it.getSensorID();
174  layer = currentVxdId.getLayerNumber();
175  ladder = currentVxdId.getLadderNumber();
176  sensor = currentVxdId.getSensorNumber();
177  // segment = currentVxdId.getSegmentNumber();// Frame nr? ... ignore
178  dhh_id = ((layer - 1) << 5) | ((ladder) << 1) | (sensor - 1);
179  if (dhh_id <= 0 || dhh_id >= 64) {
180  B2ERROR("SensorId (DHH ID) out of range: " << dhh_id);
181  continue;
182  }
183 
184  unsigned int chip_offset;
185  chip_offset = it.getChip() * 64;
186  const unsigned char* data = &it.getData()[0];
187  for (int row = 0; row < 768; row++) {
188  for (int col = 0; col < 64; col++) {
189  hrawPxdAdcMapAll->Fill(col + chip_offset + ladder * 300 - 200, 100 + row + 850 * (layer + layer + sensor - 3), *(data++));
190  }
191  }
192  }
193  }
194 }
Belle2::VxdID
Class to uniquely identify a any structure of the PXD and SVD.
Definition: VxdID.h:43
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::VxdID::getLadderNumber
baseType getLadderNumber() const
Get the ladder id.
Definition: VxdID.h:108
Belle2::VXD::GeoCache::validSensorID
bool validSensorID(Belle2::VxdID id) const
Check that id is a valid sensor number.
Definition: GeoCache.cc:53
Belle2::PXD::PXDRawDQMModule
The raw PXD DQM module.
Definition: PXDRawDQMModule.h:45
Belle2::VXD
Namespace to provide code needed by both Vertex Detectors, PXD and SVD, and also testbeam telescopes.
Definition: GeoCache.h:36
Belle2::PXD
Namespace to encapsulate code needed for simulation and reconstrucion of the PXD.
Definition: PXDCalibrationUtilities.h:28
Belle2::VXD::GeoCache::getInstance
static GeoCache & getInstance()
Return a reference to the singleton instance.
Definition: GeoCache.cc:215
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::VxdID::getSensorNumber
baseType getSensorNumber() const
Get the sensor id.
Definition: VxdID.h:110
Belle2::PXDDAQDHEStatus
The PXD DAQ DHE Status class.
Definition: PXDDAQDHEStatus.h:46
Belle2::PXDDAQDHEStatus::getDHEID
unsigned short getDHEID(void) const
Get DHE ID of sensor.
Definition: PXDDAQDHEStatus.h:105
Belle2::PXDDAQDHEStatus::getTriggerGate
unsigned short getTriggerGate(void) const
get Trigger Gate
Definition: PXDDAQDHEStatus.h:120
Belle2::VxdID::getLayerNumber
baseType getLayerNumber() const
Get the layer id.
Definition: VxdID.h:106
Belle2::HistoModule
HistoModule.h is supposed to be used instead of Module.h for the modules with histogram definitions t...
Definition: HistoModule.h:29