Belle II Software  release-05-01-25
PXDClusterPositionCollectorModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2013 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Benjamin Schwenker *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include <pxd/modules/pxdClusterPositionCollector/PXDClusterPositionCollectorModule.h>
12 #include <pxd/dataobjects/PXDDigit.h>
13 #include <pxd/dataobjects/PXDTrueHit.h>
14 #include <vxd/geometry/GeoCache.h>
15 #include <pxd/reconstruction/PXDClusterPositionEstimator.h>
16 #include <pxd/geometry/SensorInfo.h>
17 
18 #include <TTree.h>
19 #include <TMath.h>
20 #include <TH2F.h>
21 
22 #include <boost/format.hpp>
23 #include <cmath>
24 
25 using namespace std;
26 using boost::format;
27 using namespace Belle2;
28 
29 //-----------------------------------------------------------------
30 // Register the Module
31 //-----------------------------------------------------------------
32 REG_MODULE(PXDClusterPositionCollector)
33 
34 //-----------------------------------------------------------------
35 // Implementation
36 //-----------------------------------------------------------------
37 
39  , m_clusterEta(0), m_positionOffsetU(0), m_positionOffsetV(0), m_sizeV(0), m_pitchV(0)
40 {
41  // Set module properties
42  setDescription("Calibration collector module for PXD cluster position estimation");
43  setPropertyFlags(c_ParallelProcessingCertified);
44 
45  addParam("clustersName", m_storeClustersName, "Name of the collection to use for PXDClusters", string(""));
46  addParam("clusterKind", m_clusterKind, "Collect data for Clusterkind", int(0));
47  addParam("binsU", m_binsU, "Number of bins for thetaU ", int(18));
48  addParam("binsV", m_binsV, "Number of bins for thetaV ", int(18));
49 }
50 
51 void PXDClusterPositionCollectorModule::prepare() // Do your initialise() stuff here
52 {
53  m_pxdCluster.isRequired();
54 
55  // Data object creation --------------------------------------------------
56  string gridname = str(format("GridKind_%1%") % m_clusterKind);
57  auto grid = new TH2F(gridname.c_str(), gridname.c_str(), m_binsU, -90.0, +90.0, m_binsV, -90.0, +90.0);
58  registerObject<TH2F>(gridname, grid);
59 
60  string pitchtreename = "pitchtree";
61  auto ptree = new TTree(pitchtreename.c_str(), pitchtreename.c_str());
62  ptree->Branch<int>("ClusterKind", &m_clusterKind);
63  ptree->Branch<float>("PitchV", &m_pitchV);
64  registerObject<TTree>(pitchtreename, ptree);
65 
66  for (auto uBin = 1; uBin <= grid->GetXaxis()->GetNbins(); uBin++) {
67  for (auto vBin = 1; vBin <= grid->GetYaxis()->GetNbins(); vBin++) {
68  string treename = str(format("tree_%1%_%2%_%3%") % m_clusterKind % uBin % vBin);
69  auto tree = new TTree(treename.c_str(), treename.c_str());
70  tree->Branch<string>("ShapeName", &m_shapeName);
71  tree->Branch<string>("MirroredShapeName", &m_mirroredShapeName);
72  tree->Branch<float>("ClusterEta", &m_clusterEta);
73  tree->Branch<float>("OffsetU", &m_positionOffsetU);
74  tree->Branch<float>("OffsetV", &m_positionOffsetV);
75  tree->Branch<int>("SizeV", &m_sizeV);
76  registerObject<TTree>(treename, tree);
77  }
78  }
79 }
80 
81 void PXDClusterPositionCollectorModule::collect() // Do your event() stuff here
82 {
83  // If no input, nothing to do
84  if (!m_pxdCluster) return;
85 
86  string gridname = str(format("GridKind_%1%") % m_clusterKind);
87  auto grid = getObjectPtr<TH2F>(gridname);
88 
89  for (auto& cluster : m_pxdCluster) {
90 
91  // Ignore clustes with more than one truehit
92  if (cluster.getRelationsTo<PXDTrueHit>().size() > 1) {
93  continue;
94  }
95 
96  // Ignore freak clusters with more than 1000 digits
97  if (cluster.getSize() >= 1000) {
98  continue;
99  }
100 
101  // Ignore clusters with wrong clusterkind
102  if (m_clusterKind != PXD::PXDClusterPositionEstimator::getInstance().getClusterkind(cluster)) {
103  continue;
104  }
105 
106  for (auto& truehit : cluster.getRelationsTo<PXDTrueHit>()) {
107 
108  auto mom = truehit.getMomentum();
109 
110  if (mom[2] == 0 || mom.Mag() < 0.02) continue;
111 
112  double thetaU = TMath::ATan2(mom[0] / mom[2], 1.0) * 180.0 / M_PI;
113  double thetaV = TMath::ATan2(mom[1] / mom[2], 1.0) * 180.0 / M_PI;
114  int uBin = grid->GetXaxis()->FindBin(thetaU);
115  int vBin = grid->GetYaxis()->FindBin(thetaV);
116 
117  if (uBin < 1 || uBin > m_binsU || vBin < 1 || vBin > m_binsV) continue;
118 
119  // Add new entry to sources
120  grid->Fill(thetaU, thetaV);
121 
122  VxdID sensorID = truehit.getSensorID();
123  const PXD::SensorInfo& Info = dynamic_cast<const PXD::SensorInfo&>(VXD::GeoCache::get(sensorID));
124 
125  // Sort all pixels related to the cluster
126  set<PXD::Pixel> pixels;
127  for (int i = 0; i < cluster.getSize(); i++) {
128  const PXDDigit* const storeDigit = cluster.getRelationsTo<PXDDigit>("PXDDigits")[i];
129  pixels.insert(PXD::Pixel(storeDigit, i));
130  }
131 
132  // Fill tree variables
133  m_shapeName = PXD::PXDClusterPositionEstimator::getInstance().getShortName(pixels, cluster.getUStart(), cluster.getVStart(),
134  cluster.getVSize(), thetaU, thetaV);
135  m_mirroredShapeName = PXD::PXDClusterPositionEstimator::getInstance().getMirroredShortName(pixels, cluster.getUStart(),
136  cluster.getVStart(), cluster.getVSize(), thetaU, thetaV);
137  m_clusterEta = PXD::PXDClusterPositionEstimator::getInstance().computeEta(pixels, cluster.getVStart(), cluster.getVSize(), thetaU,
138  thetaV);
139  m_positionOffsetU = truehit.getU() - Info.getUCellPosition(cluster.getUStart());
140  m_positionOffsetV = truehit.getV() - Info.getVCellPosition(cluster.getVStart());
141  m_pitchV = Info.getVPitch(truehit.getV());
142  m_sizeV = cluster.getVSize();
143 
144  // Fill the tree
145  string treename = str(format("tree_%1%_%2%_%3%") % m_clusterKind % uBin % vBin);
146  getObjectPtr<TTree>(treename)->Fill();
147 
148  // Fill the pitch tree (this should happen only once per collector)
149  string pitchtreename = "pitchtree";
150  if (getObjectPtr<TTree>(pitchtreename)->GetEntries() == 0) {
151  getObjectPtr<TTree>(pitchtreename)->Fill();
152  }
153  }
154  }
155 }
Belle2::CalibrationCollectorModule
Calibration collector module base class.
Definition: CalibrationCollectorModule.h:44
Belle2::VxdID
Class to uniquely identify a any structure of the PXD and SVD.
Definition: VxdID.h:43
Belle2::PXDTrueHit
Class PXDTrueHit - Records of tracks that either enter or leave the sensitive volume.
Definition: PXDTrueHit.h:35
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::PXD::Pixel
Class to represent one pixel, used in clustering for fast access.
Definition: Pixel.h:47
Belle2::PXDDigit
The PXD digit class.
Definition: PXDDigit.h:38
Belle2::VXD::SensorInfoBase::getVPitch
double getVPitch(double v=0) const
Return the pitch of the sensor.
Definition: SensorInfoBase.h:150
Belle2::PXDClusterPositionCollectorModule
Calibration collector module for PXD cluster position estimation.
Definition: PXDClusterPositionCollectorModule.h:35
Belle2::PXD::SensorInfo
Specific implementation of SensorInfo for PXD Sensors which provides additional pixel specific inform...
Definition: SensorInfo.h:34
Belle2::VXD::SensorInfoBase::getVCellPosition
double getVCellPosition(int vID) const
Return the position of a specific strip/pixel in v direction.
Definition: SensorInfoBase.h:191
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::VXD::SensorInfoBase::getUCellPosition
double getUCellPosition(int uID, int vID=-1) const
Return the position of a specific strip/pixel in u direction.
Definition: SensorInfoBase.h:179