Belle II Software development
PXDRawDQMModule.cc
1/**************************************************************************
2 * basf2 (Belle II Analysis Software Framework) *
3 * Author: The Belle II Collaboration *
4 * *
5 * See git log for contributors and copyright holders. *
6 * This file is licensed under LGPL-3.0, see LICENSE.md. *
7 **************************************************************************/
8
9#include <pxd/modules/pxdDQM/PXDRawDQMModule.h>
10
11#include <vxd/geometry/GeoCache.h>
12#include <TDirectory.h>
13#include <boost/format.hpp>
14#include <string>
15
16using namespace std;
17using namespace Belle2;
18using namespace Belle2::PXD;
19using namespace Belle2::VXD;
20
21using boost::format;
22
23//-----------------------------------------------------------------
24// Register the Module
25//-----------------------------------------------------------------
26REG_MODULE(PXDRawDQM);
27
28//-----------------------------------------------------------------
29// Implementation
30//-----------------------------------------------------------------
31
33{
34 //Set module properties
35 setDescription("Monitor raw PXD");
37 addParam("histogramDirectoryName", m_histogramDirectoryName, "Name of the directory where histograms will be placed",
38 std::string("pxdraw"));
39
40 addParam("RawPXDsName", m_storeRawPxdrarrayName, "The name of the StoreArray of RawPXDs to be processed", string(""));
41 addParam("PXDRawHitsName", m_storeRawHitsName, "The name of the StoreArray of PXDRawHits to be processed", string(""));
42 addParam("PXDRawAdcsName", m_storeRawAdcsName, "The name of the StoreArray of PXDRawAdcs to be processed", string(""));
43}
44
46{
47 // Create a separate histogram directory and cd into it.
48 TDirectory* oldDir = gDirectory;
49 if (m_histogramDirectoryName != "") {
50 oldDir->mkdir(m_histogramDirectoryName.c_str());// do not use return value with ->cd(), its ZERO if dir already exists
51 oldDir->cd(m_histogramDirectoryName.c_str());
52 }
53
54 hrawPxdPackets = new TH1F("hrawPxdPackets", "Pxd Raw Packet Nr;Nr per Event", 16, 0, 16);
55 hrawPxdPacketSize = new TH1F("hrawPxdPacketSize", "Pxd Raw Packetsize;Words per packet", 1024, 0, 1024);
56 hrawPxdHitMapAll = new TH2F("hrawPxdHitMapAll",
57 "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);
58 // ADC map not supported by DHC anymore ... deactivate filling, later remove
60 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);
61
62 hrawPxdHitsCount = new TH1F("hrawPxdCount", "Pxd Raw Count ;Nr per Event", 8192, 0, 8192);
63
65 std::vector<VxdID> sensors = vxdGeometry.getListOfSensors();
66 for (const VxdID& avxdid : sensors) {
67 const auto& info = vxdGeometry.getSensorInfo(avxdid);
68 if (info.getType() != VXD::SensorInfoBase::PXD) continue;
69 //Only interested in PXD sensors
70 TString buff = (std::string)avxdid;
71 TString buffus = buff;
72 buffus.ReplaceAll(".", "_");
73
74 m_hrawPxdHitMap[avxdid] = new TH2F(("hrawPxdHitMap_" + buffus),
75 ("Pxd Raw Hit Map, " + buff + ";column;row"), 250,
76 0, 250, 768, 0, 768);
77 m_hrawPxdChargeMap[avxdid] = new TH2F(("hrawPxdChargeMap_" + buffus),
78 ("Pxd Raw Charge Map, " + buff + ";column;row"), 250, 0, 250, 768, 0, 768);
79 m_hrawPxdHitsCharge[avxdid] = new TH1F(("hrawPxdHitsCharge_" + buffus),
80 ("Pxd Raw Hit Charge, " + buff + ";Charge"), 256, 0, 256);
81 m_hrawPxdHitTimeWindow[avxdid] = new TH1F(("hrawPxdHitTimeWindow_" + buffus),
82 ("Pxd Raw Hit Time Window (framenr*192-gate_of_hit), " + buff + ";Time [a.u.]"), 2048, -256, 2048 - 256);
83 m_hrawPxdGateTimeWindow[avxdid] = new TH1F(("hrawPxdGateTimeWindow_" + buffus),
84 ("Pxd Raw Gate Time Window (framenr*192-triggergate_of_hit), " + buff + ";Time [a.u.]"), 2048, -256, 2048 - 256);
85 }
86 // cd back to root directory
87 oldDir->cd();
88}
89
91{
92 REG_HISTOGRAM
96 m_storeDAQEvtStats.isRequired();
97}
98
100{
101 // Just to make sure, reset all the histograms.
102 if (hrawPxdPackets) hrawPxdPackets->Reset();
107
108 for (auto& h : m_hrawPxdHitMap) if (h.second) h.second->Reset();
109 for (auto& h : m_hrawPxdChargeMap) if (h.second) h.second->Reset();
110 for (auto& h : m_hrawPxdHitsCharge) if (h.second) h.second->Reset();
111 for (auto& h : m_hrawPxdHitTimeWindow) if (h.second) h.second->Reset();
112 for (auto& h : m_hrawPxdGateTimeWindow) if (h.second) h.second->Reset();
113}
114
116{
117 hrawPxdPackets->Fill(m_storeRawPxdrarray.getEntries());
118
119 for (auto& it : m_storeRawPxdrarray) {
120 if (hrawPxdPacketSize) hrawPxdPacketSize->Fill(it.size());
121 }
122
123 if (hrawPxdHitsCount) hrawPxdHitsCount->Fill(m_storeRawHits.getEntries());
124
125 for (auto& it : m_storeRawHits) {
126
127 VxdID currentVxdId = it.getSensorID();
128 auto layer = currentVxdId.getLayerNumber();
129 auto ladder = currentVxdId.getLadderNumber();
130 auto sensor = currentVxdId.getSensorNumber();
131
132 // Get startrow and DheID from DAQEvtStats
133 const PXDDAQDHEStatus* dhe = (*m_storeDAQEvtStats).findDHE(currentVxdId);
134 if (dhe == nullptr) {
135 B2ERROR("No DHE found for SensorId: " << currentVxdId);
136 continue;
137 }
138 auto startGate = dhe->getTriggerGate();
139
140 if (hrawPxdHitMapAll) hrawPxdHitMapAll->Fill(it.getColumn() + ladder * 300 - 200,
141 100 + it.getRow() + 850 * (layer + layer + sensor - 3));
142
143 if (m_hrawPxdHitMap[currentVxdId]) m_hrawPxdHitMap[currentVxdId]->Fill(it.getColumn(), it.getRow());
144 if (m_hrawPxdChargeMap[currentVxdId]) m_hrawPxdChargeMap[currentVxdId]->Fill(it.getColumn(), it.getRow(), it.getCharge());
145 if (m_hrawPxdHitsCharge[currentVxdId]) m_hrawPxdHitsCharge[currentVxdId]->Fill(it.getCharge());
146 // Is this histogram necessary? we are folding with occupancy of sensor hits here
147 // Think about 1024*framenr-hit_row?
148 if (m_hrawPxdHitTimeWindow[currentVxdId]) m_hrawPxdHitTimeWindow[currentVxdId]->Fill(it.getFrameNr() * 192 - it.getRow() / 4);
149 if (m_hrawPxdGateTimeWindow[currentVxdId]) m_hrawPxdGateTimeWindow[currentVxdId]->Fill(it.getFrameNr() * 192 - startGate);
150 }
151
152 if (hrawPxdAdcMapAll) {
153 for (auto& it : m_storeRawAdcs) {
154 int dhh_id;
155 // calculate DHH id from Vxd Id
156 unsigned int layer, ladder, sensor;//, segment;
157 VxdID currentVxdId;
158 currentVxdId = it.getSensorID();
159 layer = currentVxdId.getLayerNumber();
160 ladder = currentVxdId.getLadderNumber();
161 sensor = currentVxdId.getSensorNumber();
162 // segment = currentVxdId.getSegmentNumber();// Frame nr? ... ignore
163 dhh_id = ((layer - 1) << 5) | ((ladder) << 1) | (sensor - 1);
164 if (dhh_id <= 0 || dhh_id >= 64) {
165 B2ERROR("SensorId (DHH ID) out of range: " << dhh_id);
166 continue;
167 }
168
169 unsigned int chip_offset;
170 chip_offset = it.getChip() * 64;
171 const unsigned char* data = &it.getData()[0];
172 for (int row = 0; row < 768; row++) {
173 for (int col = 0; col < 64; col++) {
174 hrawPxdAdcMapAll->Fill(col + chip_offset + ladder * 300 - 200, 100 + row + 850 * (layer + layer + sensor - 3), *(data++));
175 }
176 }
177 }
178 }
179}
HistoModule()
Constructor.
Definition HistoModule.h:32
void setDescription(const std::string &description)
Sets the description of the module.
Definition Module.cc:214
void setPropertyFlags(unsigned int propertyFlags)
Sets the flags for the module properties.
Definition Module.cc:208
@ c_ParallelProcessingCertified
This module can be run in parallel processing mode safely (All I/O must be done through the data stor...
Definition Module.h:80
The PXD DAQ DHE Status class.
unsigned short getTriggerGate(void) const
get Trigger Gate
TH2F * hrawPxdAdcMapAll
Histogram Adc 2d hitmap (full frames only) (all pxd)
void initialize() override final
Initialize.
std::string m_storeRawAdcsName
RawAdcs StoreArray name.
TH1F * hrawPxdPacketSize
Histogram raw packet size.
StoreObjPtr< PXDDAQStatus > m_storeDAQEvtStats
Input array for DAQ Status.
StoreArray< RawPXD > m_storeRawPxdrarray
Storearray for raw data packets.
std::map< VxdID, TH2F * > m_hrawPxdChargeMap
Histogram 2d chargemap.
void defineHisto() override final
Define histograms.
std::string m_storeRawHitsName
PXDRawHits StoreArray name.
std::map< VxdID, TH1F * > m_hrawPxdGateTimeWindow
Histogram raw pixel trigger gate window.
std::string m_storeRawPxdrarrayName
RawPXD StoreArray name.
void event() override final
Event.
TH1F * hrawPxdHitsCount
Histogram pixelcount/?
TH2F * hrawPxdHitMapAll
Histogram 2d hitmap (all)
std::map< VxdID, TH1F * > m_hrawPxdHitTimeWindow
Histogram raw pixel hit "time" window.
PXDRawDQMModule()
Constructor defining the parameters.
StoreArray< PXDRawAdc > m_storeRawAdcs
Storearray for ADC from full frames.
std::string m_histogramDirectoryName
Name of the histogram directory in ROOT file.
std::map< VxdID, TH2F * > m_hrawPxdHitMap
Histogram 2d hitmap.
TH1F * hrawPxdPackets
Histogram number of raw packets.
void beginRun() override final
Begin run.
std::map< VxdID, TH1F * > m_hrawPxdHitsCharge
Histogram raw pixel charge.
StoreArray< PXDRawHit > m_storeRawHits
Storearray for raw pixels.
Class to facilitate easy access to sensor information of the VXD like coordinate transformations or p...
Definition GeoCache.h:39
const std::vector< VxdID > getListOfSensors() const
Get list of all sensors.
Definition GeoCache.cc:59
const SensorInfoBase & getSensorInfo(Belle2::VxdID id) const
Return a reference to the SensorInfo of a given SensorID.
Definition GeoCache.cc:67
static GeoCache & getInstance()
Return a reference to the singleton instance.
Definition GeoCache.cc:214
Class to uniquely identify a any structure of the PXD and SVD.
Definition VxdID.h:33
baseType getSensorNumber() const
Get the sensor id.
Definition VxdID.h:100
baseType getLadderNumber() const
Get the ladder id.
Definition VxdID.h:98
baseType getLayerNumber() const
Get the layer id.
Definition VxdID.h:96
void addParam(const std::string &name, T &paramVariable, const std::string &description, const T &defaultValue)
Adds a new parameter to the module.
Definition Module.h:559
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition Module.h:649
Namespace to encapsulate code needed for simulation and reconstrucion of the PXD.
Namespace to provide code needed by both Vertex Detectors, PXD and SVD, and also testbeam telescopes.
Definition GeoCache.h:34
Abstract base class for different kinds of events.
STL namespace.