Belle II Software  release-05-01-25
PXDTrackClusterDQMModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2019 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Björn Spruck *
7  * Use: PXD-Tracking-Cluster DQM *
8  * *
9  * This software is provided "as is" without any warranty. *
10  **************************************************************************/
11 
12 #include <pxd/modules/pxdDQM/PXDTrackClusterDQMModule.h>
13 #include <pxd/unpacking/PXDMappingLookup.h>
14 #include <pxd/geometry/SensorInfo.h>
15 #include <TDirectory.h>
16 
17 using namespace Belle2;
18 using namespace Belle2::PXD;
19 
20 //-----------------------------------------------------------------
21 // Register the Module
22 //-----------------------------------------------------------------
23 REG_MODULE(PXDTrackClusterDQM)
24 
25 
26 //-----------------------------------------------------------------
27 // Implementation
28 //-----------------------------------------------------------------
29 
30 PXDTrackClusterDQMModule::PXDTrackClusterDQMModule() : HistoModule(), m_vxdGeometry(VXD::GeoCache::getInstance())
31 {
32  //Set module properties
33  setDescription("DQM for PXD Cluster matched to a Track");
34  setPropertyFlags(c_ParallelProcessingCertified);
35 
36  addParam("histogramDirectoryName", m_histogramDirectoryName, "Name of the directory where histograms will be placed",
37  std::string("PXDER"));
38  addParam("moreHistos", m_moreHistos, "Fill additional histograms (not for ereco)", false);
39  addParam("ASICHistos", m_ASICHistos, "Fill additional histograms ASIC combination", true);
40 }
41 
42 
43 //------------------------------------------------------------------
44 // Function to define histograms
45 //-----------------------------------------------------------------
46 
48 {
49  m_recoTracks.isOptional(m_RecoTracksStoreArrayName);
50  m_tracks.isOptional(m_TracksStoreArrayName);
51 
52  // Register histograms (calls back defineHisto)
53  REG_HISTOGRAM
54 
55 }
56 
58 {
59  // Create a separate histogram directories and cd into it.
60  TDirectory* oldDir = gDirectory;
61  if (m_histogramDirectoryName != "") {
62  oldDir->mkdir(m_histogramDirectoryName.c_str());// do not use return value with ->cd(), its ZERO if dir already exists
63  oldDir->cd(m_histogramDirectoryName.c_str());
64  }
65 
66  std::vector<VxdID> sensors = m_vxdGeometry.getListOfSensors();
67  for (VxdID& avxdid : sensors) {
68  VXD::SensorInfoBase info = m_vxdGeometry.getSensorInfo(avxdid);
69  if (info.getType() != VXD::SensorInfoBase::PXD) continue;
70  //Only interested in PXD sensors
71 
72  TString buff = (std::string)avxdid;
73  buff.ReplaceAll(".", "_");
74 
75  m_trackClusterCharge[avxdid] = new TH1F("PXD_Track_Cluster_Charge_" + buff, "PXD Track Cluster Charge " + buff + ";Charge/ADU;",
76  100, 0, 100);
77  if (m_moreHistos) {
78  m_trackClusterChargeUC[avxdid] = new TH1F("PXD_Track_Cluster_Charge_UC_" + buff,
79  "PXD Track Cluster Charge (uncorrected)" + buff + ";Charge/ADU;", 100, 0, 100);
80  }
81  if (m_ASICHistos) {
82  // for now, we only want to have this module in
83  if (avxdid == VxdID("1.5.1")) {
84  for (int s = 0; s < 6; s++) {
85  for (int d = 0; d < 4; d++) {
86  m_trackASICClusterCharge[avxdid][s][d] = new TH1F("PXD_Track_Cluster_Charge_" + buff + Form("_sw%d_dcd%d", s + 1, d + 1),
87  "PXD Track Cluster Charge " + buff + Form(" sw%d dcd%d ", s + 1, d + 1) + ";Charge/ADU;", 100, 0, 100);
88  }
89  }
90  }
91  }
92  }
93 
94  m_trackedClusters = new TH1F("PXD_Tracked_Clusters", "PXD_Tracked_Clusters", 64, 0, 64);
95 
96  for (auto i = 0; i < 64; i++) {
97  auto layer = (((i >> 5) & 0x1) + 1);
98  auto ladder = ((i >> 1) & 0xF);
99  auto sensor = ((i & 0x1) + 1);
100 
101  auto id = Belle2::VxdID(layer, ladder, sensor);
102  // Check if sensor exist
103  if (Belle2::VXD::GeoCache::getInstance().validSensorID(id)) {
104  m_vxd_to_dhe[id] = i;
105  }
106  }
107 
108  oldDir->cd();
109 
110 }
111 
113 {
114  for (auto& it : m_trackClusterCharge) if (it.second) it.second->Reset();
115  for (auto& it : m_trackClusterChargeUC) if (it.second) it.second->Reset();
116  if (m_trackedClusters) m_trackedClusters->Reset();
117  for (const auto& it1 : m_trackASICClusterCharge) {
118  for (const auto& it2 : it1.second) {
119  for (const auto& it3 : it2) {
120  if (it3) it3->Reset();
121  }
122  }
123  }
124 }
125 
126 
128 {
129  if (m_trackedClusters) m_trackedClusters->Fill(-1); // Underflow as event counter
130  for (const Track& track : m_tracks) {
134  RelationVector<RecoTrack> recoTrack = track.getRelationsTo<RecoTrack>(m_RecoTracksStoreArrayName);
135  if (!recoTrack.size()) continue;
136  RelationVector<PXDCluster> pxdClustersTrack = DataStore::getRelationsWithObj<PXDCluster>(recoTrack[0]);
137 
138  const TrackFitResult* tfr = track.getTrackFitResultWithClosestMass(Const::pion);
139  double correction = 1.0;
140  if (tfr) correction = sin(tfr->getMomentum().Theta());
141  for (auto& cluster : pxdClustersTrack) {
142  m_trackedClusters->Fill(m_vxd_to_dhe[cluster.getSensorID()]);
143  if (m_trackClusterChargeUC[cluster.getSensorID()]) m_trackClusterChargeUC[cluster.getSensorID()]->Fill(cluster.getCharge());
144  if (tfr && m_trackClusterCharge[cluster.getSensorID()]) m_trackClusterCharge[cluster.getSensorID()]->Fill(
145  cluster.getCharge()*correction);
146  if (m_ASICHistos && tfr) {
147  // for now, we only want to have this module in
148  if (cluster.getSensorID() == VxdID("1.5.1")) {
149  auto SensorInfo = dynamic_cast<const PXD::SensorInfo&>(VXD::GeoCache::get(cluster.getSensorID()));
150  auto d = PXDMappingLookup::getDCDID(SensorInfo.getUCellID(cluster.getU()), SensorInfo.getVCellID(cluster.getV()),
151  cluster.getSensorID());
152  auto s = PXDMappingLookup::getSWBID(SensorInfo.getVCellID(cluster.getV()));
153 
154  TH1F* h = nullptr;
155  try {
156  h = m_trackASICClusterCharge[cluster.getSensorID()].at(s - 1).at(d - 1);
157  } catch (...) {
158  }
159  if (h) h->Fill(cluster.getCharge()*correction);
160  }
161  }
162  }
163  }
164 }
Belle2::RelationVector::size
size_t size() const
Get number of relations.
Definition: RelationVector.h:98
Belle2::VxdID
Class to uniquely identify a any structure of the PXD and SVD.
Definition: VxdID.h:43
Belle2::VXD::SensorInfoBase::getUCellID
int getUCellID(double u, double v=0, bool clamp=false) const
Return the corresponding pixel/strip ID of a given u coordinate.
Definition: SensorInfoBase.h:203
Belle2::PXDTrackClusterDQMModule::initialize
void initialize() override final
Module functions.
Definition: PXDTrackClusterDQMModule.cc:47
Belle2::TrackFitResult::getMomentum
TVector3 getMomentum() const
Getter for vector of momentum at closest approach of track in r/phi projection.
Definition: TrackFitResult.h:116
Belle2::VXD::GeoCache::get
static const SensorInfoBase & get(Belle2::VxdID id)
Return a reference to the SensorInfo of a given SensorID.
Definition: GeoCache.h:141
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::PXDTrackClusterDQMModule::event
void event() override final
Function to process event record.
Definition: PXDTrackClusterDQMModule.cc:127
Belle2::PXD::PXDMappingLookup::getSWBID
static int getSWBID(const int v)
get ID of SWB for giving pixel, range: 1..6.
Definition: PXDMappingLookup.cc:48
Belle2::VXD::SensorInfoBase
Base class to provide Sensor Information for PXD and SVD.
Definition: SensorInfoBase.h:40
Belle2::PXDTrackClusterDQMModule
DQM of cluster in matched tracks.
Definition: PXDTrackClusterDQMModule.h:46
Belle2::TrackFitResult
Values of the result of a track fit with a given particle hypothesis.
Definition: TrackFitResult.h:59
Belle2::RecoTrack
This is the Reconstruction Event-Data Model Track.
Definition: RecoTrack.h:78
Belle2::Const::pion
static const ChargedStable pion
charged pion particle
Definition: Const.h:535
Belle2::PXDTrackClusterDQMModule::defineHisto
void defineHisto() override final
Function to define histograms.
Definition: PXDTrackClusterDQMModule.cc:57
Belle2::PXD::SensorInfo
Specific implementation of SensorInfo for PXD Sensors which provides additional pixel specific inform...
Definition: SensorInfo.h:34
Belle2::RelationVector
Class for type safe access to objects that are referred to in relations.
Definition: DataStore.h:38
Belle2::PXD
Namespace to encapsulate code needed for simulation and reconstrucion of the PXD.
Definition: PXDCalibrationUtilities.h:28
Belle2::VXD::SensorInfoBase::getVCellID
int getVCellID(double v, bool clamp=false) const
Return the corresponding pixel/strip ID of a given v coordinate.
Definition: SensorInfoBase.h:213
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::VXD::SensorInfoBase::PXD
@ PXD
PXD Sensor.
Definition: SensorInfoBase.h:44
Belle2::Track
Class that bundles various TrackFitResults.
Definition: Track.h:35
Belle2::HistoModule
HistoModule.h is supposed to be used instead of Module.h for the modules with histogram definitions t...
Definition: HistoModule.h:29
Belle2::PXDTrackClusterDQMModule::beginRun
void beginRun() override final
Function to process begin_run record.
Definition: PXDTrackClusterDQMModule.cc:112
Belle2::PXD::PXDMappingLookup::getDCDID
static int getDCDID(const int u, const int v, const VxdID sensorID)
get ID of DCD for giving pixel, range: 1..4.
Definition: PXDMappingLookup.cc:22