Belle II Software development
PXDDQMCorrModule.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/PXDDQMCorrModule.h>
10#include <pxd/dataobjects/PXDCluster.h>
11
12#include <framework/core/HistoModule.h>
13
14#include "TDirectory.h"
15
16using namespace std;
17using namespace Belle2;
18
19//-----------------------------------------------------------------
20// Register the Module
21//-----------------------------------------------------------------
22REG_MODULE(PXDDQMCorr);
23
24
25//-----------------------------------------------------------------
26// Implementation
27//-----------------------------------------------------------------
28
30{
31 //Set module properties
32 setDescription("PXD DQM Correlation module");
33 setPropertyFlags(c_ParallelProcessingCertified); // specify this flag if you need parallel processing
34 addParam("histogramDirectoryName", m_histogramDirectoryName, "Name of the directory where histograms will be placed",
35 std::string("pxd"));
36}
37
38
39//------------------------------------------------------------------
40// Function to define histograms
41//-----------------------------------------------------------------
42
44{
45 // Create a separate histogram directory and cd into it.
46 TDirectory* oldDir = gDirectory;
47 if (m_histogramDirectoryName != "") {
48 oldDir->mkdir(m_histogramDirectoryName.c_str());// do not use return value with ->cd(), its ZERO if dir already exists
49 oldDir->cd(m_histogramDirectoryName.c_str());
50 }
51 //----------------------------------------------------------------
52
53 int nPixelsU1 = 1;//getInfo(0).getUCells();// we need an actual sensor ID to use that
54 int nPixelsV1 = 3;//getInfo(0).getVCells();
55 int nPixelsU2 = 1;//getInfo(1).getUCells();
56 int nPixelsV2 = 3;//getInfo(1).getVCells();
57 m_CorrelationU = new TH2F("CorrelationU", "Correlation of U;U1/cm;U2/cm", 25, -nPixelsU1, nPixelsU1, 25, nPixelsU2, nPixelsU2);
58 m_CorrelationV = new TH2F("CorrelationV", "Correlation of V;V1/cm;V2/cm", 50, -nPixelsV1, nPixelsV1, 50, nPixelsV2, nPixelsV2);
59 m_DeltaU = new TH1F("DeltaU", "Correlation of U2-U1;Udiff/cm", 100, -nPixelsU1, nPixelsU2);
60 m_DeltaV = new TH1F("DeltaV", "Correlation of V2-V1;Vdiff/cm", 200, -nPixelsV1, nPixelsV2);
61
62 // cd back to root directory
63 oldDir->cd();
64}
65
66
68{
69 // Register histograms (calls back defineHisto)
70 REG_HISTOGRAM
71
72 //Register collections
74}
75
77{
78 // Just to make sure, reset all the histograms.
79 m_CorrelationU->Reset();
80 m_CorrelationV->Reset();
81 m_DeltaU->Reset();
82 m_DeltaV->Reset();
83}
84
85
87{
88
89 // If there are no clusters, leave
90 if (!m_storeClusters || !m_storeClusters.getEntries()) return;
91
92 for (const PXDCluster& cluster1 : m_storeClusters) {
93 int iPlane1 = cluster1.getSensorID().getLayerNumber();
94 if (iPlane1 == 0) {
95 for (const PXDCluster& cluster2 : m_storeClusters) {
96 int iPlane2 = cluster2.getSensorID().getLayerNumber();
97 if (iPlane2 == 1) {
98 m_CorrelationU->Fill(cluster1.getU(), cluster2.getU());
99 m_CorrelationV->Fill(cluster1.getV(), cluster2.getV());
100 m_DeltaU->Fill(cluster2.getU() - cluster1.getU());
101 m_DeltaV->Fill(cluster2.getV() - cluster1.getV());
102 }
103 }
104 }
105 }
106}
107
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 Cluster class This class stores all information about reconstructed PXD clusters The position...
Definition PXDCluster.h:30
TH2F * m_CorrelationV
Correlation Sensor 1 vs 2.
void initialize() override final
Initialize.
StoreArray< PXDCluster > m_storeClusters
Storearray for clusters.
void defineHisto() override final
Histogram definitions such as TH1(), TH2(), TNtuple(), TTree()....
void event() override final
Event.
std::string m_histogramDirectoryName
Name of the histogram directory in ROOT file.
TH2F * m_CorrelationU
Correlation Sensor 1 vs 2.
TH1F * m_DeltaV
Correlation Sensor 1 vs 2.
std::string m_storeClustersName
PXDClusters StoreArray name.
void beginRun() override final
Begin run.
TH1F * m_DeltaU
Correlation Sensor 1 vs 2.
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
Abstract base class for different kinds of events.
STL namespace.