Belle II Software development
PXDROIPlotModule.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/pxdHelper/PXDROIPlotModule.h>
10#include <framework/dataobjects/EventMetaData.h>
11#include <tracking/dataobjects/ROIid.h>
12#include <pxd/dataobjects/PXDRawHit.h>
13
14#include <TBox.h>
15#include <TCanvas.h>
16#include <TStyle.h>
17#include <TH2.h>
18#include <boost/format.hpp>
19
20using namespace std;
21using namespace Belle2;
22using namespace Belle2::PXD;
23
24//-----------------------------------------------------------------
25// Register the Module
26//-----------------------------------------------------------------
27REG_MODULE(PXDROIPlot);
28
29//-----------------------------------------------------------------
30// Implementation
31//-----------------------------------------------------------------
32
34{
35 //Set module properties
36 setDescription("Plot ROIs on PXD Hit/Charge Maps and write pictures");
38 addParam("PXDRawHitsName", m_storeRawHitsName, "The name of the StoreArray of PXDRawHits to be processed", string(""));
39 // we do not want to have a default name, as the user should be very sure what he is plotting!
40 addParam("ROIsName", m_ROIsName, "name of the list of ROIs (plotted in red)", string("__unusedROIs__"));
41 addParam("DCROIsName", m_DCROIsName, "name of the list of DatCon ROIs (optional, plotted in blue)", string("__unusedROIs__"));
42 addParam("HLTROIsName", m_HLTROIsName, "name of the list of HLT ROIs (optional, plotted in green)", string("__unusedROIs__"));
43 addParam("Prefix", m_prefix, "prefix for plots", string(""));
44}
45
47{
48 m_eventMetaData.isRequired();
50 // all ROIs are optional, allowing to use always the same color for comparison even if one category is not there
51 m_storeROIs.isOptional(m_ROIsName);
52 m_storeDCROIs.isOptional(m_DCROIsName);
53 m_storeHLTROIs.isOptional(m_HLTROIsName);
54}
55
57{
58 unsigned int evtNr = m_eventMetaData->getEvent();
59 unsigned int evtRun = m_eventMetaData->getRun();
60 map <VxdID, bool> flag; // flag sensors with ROIs
61 map <VxdID, vector <ROIid>> list;
62 map <VxdID, vector <ROIid>> listDC;
63 map <VxdID, vector <ROIid>> listHLT;
64 string histoname;
65 TCanvas* c1;
66 gStyle->SetPalette(55);
67 gStyle->SetOptStat(0);
68 c1 = new TCanvas("c1", "c1", 1200, 600);
69 TH2F* h = new TH2F("ChargemapWithROIs", "PXD Module Chargemap;VCell;UCell", 768, 0, 768, 250, 0, 250);
70 h->SetContour(100);
71 for (auto& it : m_storeROIs) {
72 list[it.getSensorID()].push_back(it);
73 flag[it.getSensorID()] = true;
74 }
75 for (auto& it : m_storeDCROIs) {
76 listDC[it.getSensorID()].push_back(it);
77 flag[it.getSensorID()] = true;
78 }
79 for (auto& it : m_storeHLTROIs) {
80 listHLT[it.getSensorID()].push_back(it);
81 flag[it.getSensorID()] = true;
82 }
83
84 for (auto& f : flag) {
85 c1->Clear();
86 c1->cd();
87 h->Reset();
88
89 VxdID currentVxdId = f.first;
90 histoname = m_prefix + boost::str(boost::format("Run_%d_Evt_%d_") % evtRun % evtNr) + string(currentVxdId);
91 h->SetTitle(histoname.data());
92
93 for (auto& pix : m_storeRawHits) {
94
95 if (currentVxdId != pix.getSensorID()) continue;
96 h->Fill(pix.getRow(), pix.getColumn(), pix.getCharge());
97
98 }
99
100 h->Draw("colz");
101
102 // Alpha Blending only works in ROOT's Batch Mode, but LineStyle won't.
103 // LineStyle works in normal mode, but Alpha won't
104 for (auto& it : list[currentVxdId]) {
105 TBox* b;
106 b = new TBox(it.getMinVid(), it.getMinUid(), it.getMaxVid(), it.getMaxUid());
107 b->SetLineColorAlpha(kRed, 0.3);
108 b->SetLineWidth(3);
109 b->SetFillStyle(0);// Hollow
110 b->Draw();
111 }
112 // we move the other box by more than half a bin. this is needed as alpha seems not to work
113 for (auto& it : listDC[currentVxdId]) {
114 TBox* b;
115 b = new TBox(it.getMinVid() + 0.7, it.getMinUid() + 0.7, it.getMaxVid() + 0.7, it.getMaxUid() + 0.7);
116 b->SetLineColorAlpha(kBlue, 0.3);
117 b->SetLineWidth(2);
118 b->SetLineStyle(2);
119 b->SetFillStyle(0);// Hollow
120 b->Draw();
121 }
122 // we move the other box by half a bin. this is needed as alpha seems not to work, in addition we use a dashed style
123 // dashed style doesnt work with png export, thus if all ROIs are identical, lines might overlap completely
124 for (auto& it : listHLT[currentVxdId]) {
125 TBox* b;
126 b = new TBox(it.getMinVid() - 0.7, it.getMinUid() - 0.7, it.getMaxVid() - 0.7, it.getMaxUid() - 0.7);
127 b->SetLineColorAlpha(kGreen, 0.3);
128 b->SetLineWidth(2);
129 b->SetLineStyle(7);
130 b->SetFillStyle(0);// Hollow
131 b->Draw();
132 }
133 c1->Print((histoname + ".png").data());
134 c1->Print((histoname + ".root").data());
135 }
136 delete h;
137 delete c1;
138}
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
Module()
Constructor.
Definition Module.cc:30
@ 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
void initialize() override final
Initialize.
StoreArray< ROIid > m_storeROIs
Storearray for ROIs.
std::string m_prefix
prefix for plot names
std::string m_DCROIsName
DC ROIs StoreArray name.
std::string m_storeRawHitsName
PXDRawHits StoreArray name.
PXDROIPlotModule()
Constructor defining the parameters.
StoreObjPtr< EventMetaData > m_eventMetaData
Event Metadata StorePtr.
void event() override final
Event.
std::string m_ROIsName
ROIs StoreArray name.
StoreArray< ROIid > m_storeDCROIs
Storearray for ROIs (DC)
StoreArray< PXDRawHit > m_storeRawHits
Storearray for raw pixels.
std::string m_HLTROIsName
HLT ROIs StoreArray name.
StoreArray< ROIid > m_storeHLTROIs
Storearray for ROIs (HLT)
Class to uniquely identify a any structure of the PXD and SVD.
Definition VxdID.h:32
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.
Abstract base class for different kinds of events.
STL namespace.