Belle II Software development
PXDClusterPositionCollectorModule.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/pxdClusterPositionCollector/PXDClusterPositionCollectorModule.h>
10#include <pxd/dataobjects/PXDCluster.h>
11#include <pxd/dataobjects/PXDDigit.h>
12#include <pxd/dataobjects/PXDTrueHit.h>
13#include <vxd/geometry/GeoCache.h>
14#include <pxd/reconstruction/PXDClusterPositionEstimator.h>
15#include <pxd/reconstruction/Pixel.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
25using namespace std;
26using boost::format;
27using namespace Belle2;
28
29//-----------------------------------------------------------------
30// Register the Module
31//-----------------------------------------------------------------
32REG_MODULE(PXDClusterPositionCollector);
33
34//-----------------------------------------------------------------
35// Implementation
36//-----------------------------------------------------------------
37
40{
41 // Set module properties
42 setDescription("Calibration collector module for PXD cluster position estimation");
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
51void 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
81void 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.Z() == 0 || mom.R() < 0.02) continue;
111
112 double thetaU = TMath::ATan2(mom.X() / mom.Z(), 1.0) * TMath::RadToDeg();
113 double thetaV = TMath::ATan2(mom.Y() / mom.Z(), 1.0) * TMath::RadToDeg();
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::getInstance().getSensorInfo(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);
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}
void registerObject(std::string name, T *obj)
Register object with a name, takes ownership, do not access the pointer beyond prepare()
CalibrationCollectorModule()
Constructor. Sets the default prefix for calibration dataobjects.
T * getObjectPtr(std::string name)
Calls the CalibObjManager to get the requested stored collector data.
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
std::string m_mirroredShapeName
Name of mirrored cluster shape.
StoreArray< PXDCluster > m_pxdCluster
< Required input for PXDCluster
std::string m_storeClustersName
Name of the collection to use for PXDClusters.
PXDClusterPositionCollectorModule()
Constructor: Sets the description, the properties and the parameters of the module.
The PXD digit class.
Definition PXDDigit.h:27
Class PXDTrueHit - Records of tracks that either enter or leave the sensitive volume.
Definition PXDTrueHit.h:31
static PXDClusterPositionEstimator & getInstance()
Main (and only) way to access the PXDClusterPositionEstimator.
const std::string getMirroredShortName(const std::set< Pixel > &pixels, int uStart, int vStart, int vSize, double thetaU, double thetaV) const
Return the mirrored name for the pixel set.
float computeEta(const std::set< Pixel > &pixels, int vStart, int vSize, double thetaU, double thetaV) const
Return the normed charge ratio between head and tail pixels (size>=2) or the charge of the seed (size...
const std::string getShortName(const std::set< Pixel > &pixels, int uStart, int vStart, int vSize, double thetaU, double thetaV) const
Return the name for the pixel set.
Class to represent one pixel, used in clustering for fast access.
Definition Pixel.h:36
Specific implementation of SensorInfo for PXD Sensors which provides additional pixel specific inform...
Definition SensorInfo.h:23
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: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
Abstract base class for different kinds of events.
STL namespace.