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>
22 #include <boost/format.hpp>
39 , m_clusterEta(0), m_positionOffsetU(0), m_positionOffsetV(0), m_sizeV(0), m_pitchV(0)
42 setDescription(
"Calibration collector module for PXD cluster position estimation");
43 setPropertyFlags(c_ParallelProcessingCertified);
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));
51 void PXDClusterPositionCollectorModule::prepare()
53 m_pxdCluster.isRequired();
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);
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);
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);
81 void PXDClusterPositionCollectorModule::collect()
84 if (!m_pxdCluster)
return;
86 string gridname = str(format(
"GridKind_%1%") % m_clusterKind);
87 auto grid = getObjectPtr<TH2F>(gridname);
89 for (
auto& cluster : m_pxdCluster) {
92 if (cluster.getRelationsTo<
PXDTrueHit>().size() > 1) {
97 if (cluster.getSize() >= 1000) {
102 if (m_clusterKind != PXD::PXDClusterPositionEstimator::getInstance().getClusterkind(cluster)) {
106 for (
auto& truehit : cluster.getRelationsTo<
PXDTrueHit>()) {
108 auto mom = truehit.getMomentum();
110 if (mom[2] == 0 || mom.Mag() < 0.02)
continue;
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);
117 if (uBin < 1 || uBin > m_binsU || vBin < 1 || vBin > m_binsV)
continue;
120 grid->Fill(thetaU, thetaV);
122 VxdID sensorID = truehit.getSensorID();
126 set<PXD::Pixel> pixels;
127 for (
int i = 0; i < cluster.getSize(); i++) {
128 const PXDDigit*
const storeDigit = cluster.getRelationsTo<
PXDDigit>(
"PXDDigits")[i];
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,
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();
145 string treename = str(format(
"tree_%1%_%2%_%3%") % m_clusterKind % uBin % vBin);
146 getObjectPtr<TTree>(treename)->Fill();
149 string pitchtreename =
"pitchtree";
150 if (getObjectPtr<TTree>(pitchtreename)->GetEntries() == 0) {
151 getObjectPtr<TTree>(pitchtreename)->Fill();