ATLAS Offline Software
PixelNtupleMaker.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
3 */
4 
5 #include "DerivationFrameworkInDet/PixelNtupleMaker.h"
6 
7 #include "xAODCore/ShallowCopy.h"
8 #include "xAODTracking/TrackParticle.h"
9 #include "TrkTrack/TrackStateOnSurface.h"
10 #include "TrkEventPrimitives/TrackStateDefs.h"
11 #include "AthContainers/ConstAccessor.h"
12 
13 #include <vector>
14 #include <string>
15 
16 #include "TLorentzVector.h"
17 
18 DerivationFramework::PixelNtupleMaker::PixelNtupleMaker(const std::string& t, const std::string& n, const IInterface* p) :
19  AthAlgTool(t,n,p)
20 {
21  declareInterface<DerivationFramework::IAugmentationTool>(this);
22 }
23 
24 DerivationFramework::PixelNtupleMaker::~PixelNtupleMaker() = default;
25 
26 StatusCode DerivationFramework::PixelNtupleMaker::initialize() {
27  ATH_CHECK(m_selector.retrieve());
28  ATH_CHECK(m_containerKey.initialize());
29  ATH_CHECK(m_measurementContainerKey.initialize());
30  ATH_CHECK(m_monitoringTracks.initialize());
31  return StatusCode::SUCCESS;
32 }
33 
34 StatusCode DerivationFramework::PixelNtupleMaker::finalize()
35 {
36  return StatusCode::SUCCESS;
37 }
38 
39 StatusCode DerivationFramework::PixelNtupleMaker::addBranches() const {
40  const EventContext& ctx = Gaudi::Hive::currentContext();
41 
42  SG::ReadHandle<xAOD::TrackParticleContainer> tracks(m_containerKey,ctx);
43  if (!tracks.isValid()) { return StatusCode::SUCCESS; }
44 
45  // Check the event contains tracks
46  SG::WriteHandle<xAOD::TrackParticleContainer> PixelMonitoringTrack(m_monitoringTracks,ctx);
47  if (PixelMonitoringTrack.record(std::make_unique<xAOD::TrackParticleContainer>(),
48  std::make_unique<xAOD::TrackParticleAuxContainer>()).isFailure()) {
49  return StatusCode::SUCCESS;
50  }
51 
52  SG::ReadHandle<xAOD::TrackMeasurementValidationContainer> pixClusters(m_measurementContainerKey,ctx);
53  if (!pixClusters.isValid()) { return StatusCode::SUCCESS; }
54 
55  // monitoring variables.
56  const int nbin_charge = 60;
57  const int nbin_tot = 120;
58  const int nbin_eta = 30;
59  const int nbin_dedx = 300;
60  const int nbin_size = 10;
61  const int nbin_reso = 80;
62  const int nbin_occ = 80;
63 
64  std::vector<std::vector<int>> clusCharge(11,std::vector<int>(nbin_charge,0));
65  std::vector<std::vector<int>> clusToT(11,std::vector<int>(nbin_tot,0));
66  std::vector<std::vector<int>> clusEta(11,std::vector<int>(nbin_eta,0));
67  std::vector<std::vector<int>> clusHitEta(11,std::vector<int>(nbin_eta,0));
68  std::vector<std::vector<int>> clusdEdx(11,std::vector<int>(nbin_dedx,0));
69  std::vector<std::vector<int>> clusSizeX(11,std::vector<int>(nbin_size,0));
70  std::vector<std::vector<int>> clusSizeZ(11,std::vector<int>(nbin_size,0));
71  std::vector<std::vector<int>> clusResidualX(11,std::vector<int>(nbin_reso,0));
72  std::vector<std::vector<int>> clusResidualY(11,std::vector<int>(nbin_reso,0));
73  std::vector<std::vector<int>> clusHole(11,std::vector<int>(nbin_eta,0));
74  std::vector<std::vector<int>> clusOcc(11,std::vector<int>(nbin_occ,0));
75 
76  std::vector<int> trkEta(nbin_eta,0);
77  std::vector<int> trkHole(nbin_eta,0);
78  std::vector<int> trkdEdx(nbin_dedx,0);
79 
80  float maxPt = 0.0;
81  const xAOD::TrackParticle* trk_maxpt = nullptr;
82 
83  std::vector<float> tmpCov(15,0.);
84 
85  // StoreGateSvc+InDetTrackParticles.msosLink
86  static const SG::ConstAccessor<MeasurementsOnTrack> acc_MeasurementsOnTrack("msosLink");
87 
88  for (const xAOD::TrackParticle* trk : *tracks) {
89  uint8_t nPixHits = 0; trk->summaryValue(nPixHits,xAOD::numberOfPixelHits);
90  uint8_t nSCTHits = 0; trk->summaryValue(nSCTHits,xAOD::numberOfSCTHits);
91  if (trk->pt()<1000.0 && nPixHits<4) { continue; }
92  if (trk->pt()<1000.0 && nSCTHits<1) { continue; }
93 
94  std::vector<uint64_t> holeIndex;
95  std::vector<int> clusterLayer;
96  std::vector<int> clusterBEC;
97  std::vector<int> clusterModulePhi;
98  std::vector<int> clusterModuleEta;
99  std::vector<float> clusterCharge;
100  std::vector<int> clusterToT;
101  std::vector<int> clusterL1A;
102  std::vector<int> clusterIsSplit;
103  std::vector<int> clusterSize;
104  std::vector<int> clusterSizePhi;
105  std::vector<int> clusterSizeZ;
106  std::vector<bool> isEdge;
107  std::vector<bool> isOverflow;
108  std::vector<float> trackPhi;
109  std::vector<float> trackTheta;
110  std::vector<float> trackX;
111  std::vector<float> trackY;
112  std::vector<float> localX;
113  std::vector<float> localY;
114  std::vector<float> globalX;
115  std::vector<float> globalY;
116  std::vector<float> globalZ;
117  std::vector<float> unbiasedResidualX;
118  std::vector<float> unbiasedResidualY;
119  std::vector<int> clusterIsolation10x2;
120  std::vector<int> clusterIsolation20x4;
121  std::vector<int> numTotalClustersPerModule;
122  std::vector<int> numTotalPixelsPerModule;
123  std::vector<float> moduleLorentzShift;
124 
125  std::vector<std::vector<int>> rdoToT;
126  std::vector<std::vector<float>> rdoCharge;
127  std::vector<std::vector<int>> rdoPhi;
128  std::vector<std::vector<int>> rdoEta;
129 
130  const MeasurementsOnTrack& measurementsOnTrack = acc_MeasurementsOnTrack(*trk);
131  for (const auto & msos_iter : measurementsOnTrack) {
132  if (!msos_iter.isValid()) { continue; }
133  const xAOD::TrackStateValidation* msos = *msos_iter;
134  if (msos->detType()!=Trk::TrackState::Pixel) { continue; } // its a pixel
135  if (msos->type()==Trk::TrackStateOnSurface::Hole) { // hole
136  holeIndex.push_back(msos->detElementId());
137  continue;
138  }
139  if (!msos->trackMeasurementValidationLink().isValid()) { continue; }
140  if (!(*(msos->trackMeasurementValidationLink()))) { continue; }
141 
142  const xAOD::TrackMeasurementValidation* msosClus = *(msos->trackMeasurementValidationLink());
143 
144  static const SG::ConstAccessor<float> chargeAcc("charge");
145  static const SG::ConstAccessor<int> layerAcc("layer");
146  static const SG::ConstAccessor<int> becAcc("bec");
147  static const SG::ConstAccessor<int> phi_moduleAcc("phi_module");
148  static const SG::ConstAccessor<int> eta_moduleAcc("eta_module");
149  static const SG::ConstAccessor<int> ToTAcc("ToT");
150  static const SG::ConstAccessor<int> LVL1AAcc("LVL1A");
151  static const SG::ConstAccessor<char> isSplitAcc("isSplit");
152  static const SG::ConstAccessor<int> nRDOAcc("nRDO");
153  static const SG::ConstAccessor<int> sizePhiAcc("sizePhi");
154  static const SG::ConstAccessor<int> sizeZAcc("sizeZ");
155  static const SG::ConstAccessor<float> LorentzShiftAcc("LorentzShift");
156  static const SG::ConstAccessor<std::vector<int> > rdo_totAcc("rdo_tot");
157  static const SG::ConstAccessor<std::vector<float> > rdo_chargeAcc("rdo_charge");
158 
159  for (const auto *clus_itr : *pixClusters) {
160  if (clus_itr->identifier()!=(msosClus)->identifier()) { continue; }
161  if (chargeAcc(*clus_itr)!=chargeAcc(*msosClus)) { continue; }
162 
163  clusterLayer.push_back(layerAcc(*clus_itr));
164  clusterBEC.push_back(becAcc(*clus_itr));
165  clusterModulePhi.push_back(phi_moduleAcc(*clus_itr));
166  clusterModuleEta.push_back(eta_moduleAcc(*clus_itr));
167  clusterCharge.push_back(chargeAcc(*clus_itr));
168  clusterToT.push_back(ToTAcc(*clus_itr));
169  clusterL1A.push_back(LVL1AAcc(*clus_itr));
170  clusterIsSplit.push_back(isSplitAcc(*clus_itr));
171  clusterSize.push_back(nRDOAcc(*clus_itr));
172  clusterSizePhi.push_back(sizePhiAcc(*clus_itr));
173  clusterSizeZ.push_back(sizeZAcc(*clus_itr));
174 
175  trackPhi.push_back(msos->localPhi());
176  trackTheta.push_back(msos->localTheta());
177  trackX.push_back(msos->localX());
178  trackY.push_back(msos->localY());
179  localX.push_back(clus_itr->localX());
180  localY.push_back(clus_itr->localY());
181  globalX.push_back(clus_itr->globalX());
182  globalY.push_back(clus_itr->globalY());
183  globalZ.push_back(clus_itr->globalZ());
184  unbiasedResidualX.push_back(msos->unbiasedResidualX());
185  unbiasedResidualY.push_back(msos->unbiasedResidualY());
186  moduleLorentzShift.push_back(LorentzShiftAcc(*clus_itr));
187 
188  // cluster isolation IBL:50x250um, PIXEL:50x400um
189  // - isolation region 10x2 = 500x500um for IBL, 500x800um for PIXEL
190  int numNeighborCluster10x2 = 0;
191  int numNeighborCluster20x4 = 0;
192  int nTotalClustersPerModule = 0;
193  int nTotalPixelsPerModule = 0;
194  for (const auto *clus_neighbor : *pixClusters) {
195  if (layerAcc(*clus_neighbor)==layerAcc(*clus_itr)
196  && becAcc(*clus_neighbor)==becAcc(*clus_itr)
197  && phi_moduleAcc(*clus_neighbor)==phi_moduleAcc(*clus_itr)
198  && eta_moduleAcc(*clus_neighbor)==eta_moduleAcc(*clus_itr)) {
199  float deltaX = std::abs(clus_neighbor->localX()-clus_itr->localX());
200  float deltaY = std::abs(clus_neighbor->localY()-clus_itr->localY());
201  nTotalClustersPerModule++;
202  nTotalPixelsPerModule += nRDOAcc(*clus_neighbor);
203  if (deltaX>0.0 && deltaY>0.0) {
204  if (layerAcc(*clus_itr)==0 && becAcc(*clus_itr)==0) { // IBL
205  if (deltaX<0.500 && deltaY<0.500) { numNeighborCluster10x2++; }
206  if (deltaX<1.000 && deltaY<1.000) { numNeighborCluster20x4++; }
207  }
208  else {
209  if (deltaX<0.500 && deltaY<0.800) { numNeighborCluster10x2++; }
210  if (deltaX<1.000 && deltaY<1.600) { numNeighborCluster20x4++; }
211  }
212  }
213  }
214  }
215  clusterIsolation10x2.push_back(numNeighborCluster10x2);
216  clusterIsolation20x4.push_back(numNeighborCluster20x4);
217  numTotalClustersPerModule.push_back(nTotalClustersPerModule);
218  numTotalPixelsPerModule.push_back(nTotalPixelsPerModule);
219 
220  // is edge pixel?
221  // contain overlflow hit?
222  bool checkEdge = false;
223  bool checkOverflow = false;
224 
225  // rdo information
226  std::vector<int> tmpToT;
227  std::vector<float> tmpCharge;
228  std::vector<int> tmpPhi;
229  std::vector<int> tmpEta;
230 
231  static const SG::ConstAccessor<std::vector<int> > rdo_phi_pixel_indexAcc("rdo_phi_pixel_index");
232  int nrdo = rdo_phi_pixel_indexAcc.isAvailable(*clus_itr) ? rdo_phi_pixel_indexAcc(*clus_itr).size() : -1;
233  for (int i=0; i<nrdo; i++) {
234 
235  int phi = rdo_phi_pixel_indexAcc(*clus_itr)[i];
236  if (phi<5) { checkEdge=true; }
237  if (phi>320) { checkEdge=true; }
238 
239  int eta = rdo_phi_pixel_indexAcc(*clus_itr)[i];
240  if (layerAcc(*clus_itr)==0 && becAcc(*clus_itr)==0) { // IBL
241  if (eta_moduleAcc(*clus_itr)>-7 && eta_moduleAcc(*clus_itr)<6) { // IBL Planar
242  if (eta<5) { checkEdge=true; }
243  if (eta>154) { checkEdge=true; }
244  }
245  else { // IBL 3D
246  if (eta<5) { checkEdge=true; }
247  if (eta>74) { checkEdge=true; }
248  }
249  }
250  else {
251  if (eta<5) { checkEdge=true; }
252  if (eta>154) { checkEdge=true; }
253  }
254 
255  int tot = rdo_totAcc(*clus_itr)[i];
256  if (layerAcc(*clus_itr)==0 && becAcc(*clus_itr)==0) { // IBL
257  if (tot==16) { checkOverflow=true; }
258  }
259  else if (layerAcc(*clus_itr)==1 && becAcc(*clus_itr)==0) { // b-layer
260  if (tot==150) { checkOverflow=true; }
261  }
262  else {
263  if (tot==255) { checkOverflow=true; }
264  }
265 
266  float charge = rdo_chargeAcc(*clus_itr)[i];
267  if (trk->pt()>2000.0) {
268  tmpToT.push_back(tot);
269  tmpCharge.push_back(charge);
270  tmpPhi.push_back(phi);
271  tmpEta.push_back(eta);
272  }
273  }
274  isEdge.push_back(checkEdge);
275  isOverflow.push_back(checkOverflow);
276 
277  rdoToT.push_back(tmpToT);
278  rdoCharge.push_back(tmpCharge);
279  rdoPhi.push_back(tmpPhi);
280  rdoEta.push_back(tmpEta);
281  break;
282  }
283  }
284 
285  if (m_storeMode==1) {
286  static const SG::Decorator<float> d0err("d0err");
287  static const SG::Decorator<float> z0err("z0err");
288  static const SG::Decorator<std::vector<uint64_t>> HoleIndex("HoleIndex");
289  static const SG::Decorator<std::vector<int>> ClusterLayer("ClusterLayer");
290  static const SG::Decorator<std::vector<int>> ClusterBEC("ClusterBEC");
291  static const SG::Decorator<std::vector<int>> ClusterModulePhi("ClusterModulePhi");
292  static const SG::Decorator<std::vector<int>> ClusterModuleEta("ClusterModuleEta");
293  static const SG::Decorator<std::vector<float>> ClusterCharge("ClusterCharge");
294  static const SG::Decorator<std::vector<int>> ClusterToT("ClusterToT");
295  static const SG::Decorator<std::vector<int>> ClusterL1A("ClusterL1A");
296  static const SG::Decorator<std::vector<int>> ClusterIsSplit("ClusterIsSplit");
297  static const SG::Decorator<std::vector<int>> ClusterSize("ClusterSize");
298  static const SG::Decorator<std::vector<int>> ClusterSizePhi("ClusterSizePhi");
299  static const SG::Decorator<std::vector<int>> ClusterSizeZ("ClusterSizeZ");
300  static const SG::Decorator<std::vector<bool>> ClusterIsEdge("ClusterIsEdge");
301  static const SG::Decorator<std::vector<bool>> ClusterIsOverflow("ClusterIsOverflow");
302  static const SG::Decorator<std::vector<float>> TrackLocalPhi("TrackLocalPhi");
303  static const SG::Decorator<std::vector<float>> TrackLocalTheta("TrackLocalTheta");
304  static const SG::Decorator<std::vector<float>> TrackLocalX("TrackLocalX");
305  static const SG::Decorator<std::vector<float>> TrackLocalY("TrackLocalY");
306  static const SG::Decorator<std::vector<float>> ClusterLocalX("ClusterLocalX");
307  static const SG::Decorator<std::vector<float>> ClusterLocalY("ClusterLocalY");
308  static const SG::Decorator<std::vector<float>> ClusterGlobalX("ClusterGlobalX");
309  static const SG::Decorator<std::vector<float>> ClusterGlobalY("ClusterGlobalY");
310  static const SG::Decorator<std::vector<float>> ClusterGlobalZ("ClusterGlobalZ");
311  static const SG::Decorator<std::vector<float>> UnbiasedResidualX("UnbiasedResidualX");
312  static const SG::Decorator<std::vector<float>> UnbiasedResidualY("UnbiasedResidualY");
313  static const SG::Decorator<std::vector<int>> ClusterIsolation10x2("ClusterIsolation10x2");
314  static const SG::Decorator<std::vector<int>> ClusterIsolation20x4("ClusterIsolation20x4");
315  static const SG::Decorator<std::vector<int>> NumTotalClustersPerModule("NumTotalClustersPerModule");
316  static const SG::Decorator<std::vector<int>> NumTotalPixelsPerModule("NumTotalPixelsPerModule");
317  static const SG::Decorator<std::vector<float>> ModuleLorentzShift("ModuleLorentzShift");
318  static const SG::Decorator<std::vector<std::vector<int>>> RdoToT("RdoToT");
319  static const SG::Decorator<std::vector<std::vector<float>>> RdoCharge("RdoCharge");
320  static const SG::Decorator<std::vector<std::vector<int>>> RdoPhi("RdoPhi");
321  static const SG::Decorator<std::vector<std::vector<int>>> RdoEta("RdoEta");
322 
323  d0err(*trk) = trk->definingParametersCovMatrixVec().at(0);
324  z0err(*trk) = trk->definingParametersCovMatrixVec().at(2);
325  HoleIndex(*trk) = std::move(holeIndex);
326  ClusterLayer(*trk) = std::move(clusterLayer);
327  ClusterBEC(*trk) = std::move(clusterBEC);
328  ClusterModulePhi(*trk) = std::move(clusterModulePhi);
329  ClusterModuleEta(*trk) = std::move(clusterModuleEta);
330  ClusterCharge(*trk) = std::move(clusterCharge);
331  ClusterToT(*trk) = std::move(clusterToT);
332  ClusterL1A(*trk) = std::move(clusterL1A);
333  ClusterIsSplit(*trk) = std::move(clusterIsSplit);
334  ClusterSize(*trk) = std::move(clusterSize);
335  ClusterSizePhi(*trk) = std::move(clusterSizePhi);
336  ClusterSizeZ(*trk) = std::move(clusterSizeZ);
337  ClusterIsEdge(*trk) = std::move(isEdge);
338  ClusterIsOverflow(*trk) = std::move(isOverflow);
339  TrackLocalPhi(*trk) = std::move(trackPhi);
340  TrackLocalTheta(*trk) = std::move(trackTheta);
341  TrackLocalX(*trk) = std::move(trackX);
342  TrackLocalY(*trk) = std::move(trackY);
343  ClusterLocalX(*trk) = std::move(localX);
344  ClusterLocalY(*trk) = std::move(localY);
345  ClusterGlobalX(*trk) = std::move(globalX);
346  ClusterGlobalY(*trk) = std::move(globalY);
347  ClusterGlobalZ(*trk) = std::move(globalZ);
348  UnbiasedResidualX(*trk) = std::move(unbiasedResidualX);
349  UnbiasedResidualY(*trk) = std::move(unbiasedResidualY);
350  ClusterIsolation10x2(*trk) = std::move(clusterIsolation10x2);
351  ClusterIsolation20x4(*trk) = std::move(clusterIsolation20x4);
352  NumTotalClustersPerModule(*trk) = std::move(numTotalClustersPerModule);
353  NumTotalPixelsPerModule(*trk) = std::move(numTotalPixelsPerModule);
354  ModuleLorentzShift(*trk) = std::move(moduleLorentzShift);
355  RdoToT(*trk) = std::move(rdoToT);
356  RdoCharge(*trk) = std::move(rdoCharge);
357  RdoPhi(*trk) = std::move(rdoPhi);
358  RdoEta(*trk) = std::move(rdoEta);
359  }
360 
361  // further track selection for slimmed track variables
362  uint8_t nPixelDeadSensor = 0; trk->summaryValue(nPixelDeadSensor,xAOD::numberOfPixelDeadSensors);
363  uint8_t numberOfPixelHits= 0; trk->summaryValue(numberOfPixelHits,xAOD::numberOfPixelHits);
364  uint8_t numberOfPixelHoles= 0; trk->summaryValue(numberOfPixelHoles,xAOD::numberOfPixelHoles);
365  uint8_t numberOfPixelOutliers= 0; trk->summaryValue(numberOfPixelOutliers,xAOD::numberOfPixelOutliers);
366  uint8_t nSCTHoles = 0; trk->summaryValue(nSCTHoles,xAOD::numberOfSCTHoles);
367  uint8_t nPixSharedHits = 0; trk->summaryValue(nPixSharedHits,xAOD::numberOfPixelSharedHits);
368  uint8_t nSCTSharedHits = 0; trk->summaryValue(nSCTSharedHits,xAOD::numberOfSCTSharedHits);
369 
370  // loose track selection
371  bool passLooseCut = static_cast<bool>(m_selector->accept(trk));
372  if (!passLooseCut) { continue; }
373 
374  //==================
375  // Efficiency check
376  //==================
377  int checkIBL(0),check3D(0),checkBLY(0),checkLY1(0),checkLY2(0),checkEA1(0),checkEA2(0),checkEA3(0),checkEC1(0),checkEC2(0),checkEC3(0);
378  for (int i=0; i<(int)clusterLayer.size(); i++) {
379  if (clusterBEC.at(i)== 0 && clusterLayer.at(i)==0) {
380  if (clusterModuleEta.at(i)>-7 && clusterModuleEta.at(i)<6) { checkIBL++; }
381  else { check3D++; }
382  }
383  else if (clusterBEC.at(i)== 0 && clusterLayer.at(i)==1) { checkBLY++; }
384  else if (clusterBEC.at(i)== 0 && clusterLayer.at(i)==2) { checkLY1++; }
385  else if (clusterBEC.at(i)== 0 && clusterLayer.at(i)==3) { checkLY2++; }
386  else if (clusterBEC.at(i)== 2 && clusterLayer.at(i)==0) { checkEA1++; }
387  else if (clusterBEC.at(i)== 2 && clusterLayer.at(i)==1) { checkEA2++; }
388  else if (clusterBEC.at(i)== 2 && clusterLayer.at(i)==2) { checkEA3++; }
389  else if (clusterBEC.at(i)==-2 && clusterLayer.at(i)==0) { checkEC1++; }
390  else if (clusterBEC.at(i)==-2 && clusterLayer.at(i)==1) { checkEC2++; }
391  else if (clusterBEC.at(i)==-2 && clusterLayer.at(i)==2) { checkEC3++; }
392  }
393 
394  if (trk->pt()>2000.0 && nPixelDeadSensor==0) {
395  int ietabin = trunc((trk->eta()+3.0)/0.2);
396 
397  std::for_each((trkEta.begin()+ietabin),(trkEta.begin()+ietabin+1),[](int &n){ n++; });
398  if (numberOfPixelHoles>0) { std::for_each((trkHole.begin()+ietabin),(trkHole.begin()+ietabin+1),[](int &n){ n++; }); }
399 
400  if (checkBLY>0 && checkLY1>0 && checkLY2>0) { std::for_each((clusEta[0].begin()+ietabin),(clusEta[0].begin()+ietabin+1),[](int &n){ n++; }); if (checkIBL>0) { std::for_each((clusHitEta[0].begin()+ietabin),(clusHitEta[0].begin()+ietabin+1),[](int &n){ n++; }); }}
401  if (checkBLY>0 && checkEC2>0 && checkEC3>0) { std::for_each((clusEta[1].begin()+ietabin),(clusEta[1].begin()+ietabin+1),[](int &n){ n++; }); if (check3D>0) { std::for_each((clusHitEta[1].begin()+ietabin),(clusHitEta[1].begin()+ietabin+1),[](int &n){ n++; }); }}
402  if (checkIBL>0 && checkLY1>0 && checkLY2>0) { std::for_each((clusEta[2].begin()+ietabin),(clusEta[2].begin()+ietabin+1),[](int &n){ n++; }); if (checkBLY>0) { std::for_each((clusHitEta[2].begin()+ietabin),(clusHitEta[2].begin()+ietabin+1),[](int &n){ n++; }); }}
403  if (checkIBL>0 && checkBLY>0 && checkLY2>0) { std::for_each((clusEta[3].begin()+ietabin),(clusEta[3].begin()+ietabin+1),[](int &n){ n++; }); if (checkLY1>0) { std::for_each((clusHitEta[3].begin()+ietabin),(clusHitEta[3].begin()+ietabin+1),[](int &n){ n++; }); }}
404  if (checkIBL>0 && checkBLY>0 && checkLY1>0) { std::for_each((clusEta[4].begin()+ietabin),(clusEta[4].begin()+ietabin+1),[](int &n){ n++; }); if (checkLY2>0) { std::for_each((clusHitEta[4].begin()+ietabin),(clusHitEta[4].begin()+ietabin+1),[](int &n){ n++; }); }}
405  if (checkIBL>0 && checkEA2>0 && checkEA3>0) { std::for_each((clusEta[5].begin()+ietabin),(clusEta[5].begin()+ietabin+1),[](int &n){ n++; }); if (checkEA1>0) { std::for_each((clusHitEta[5].begin()+ietabin),(clusHitEta[5].begin()+ietabin+1),[](int &n){ n++; }); }}
406  if (checkIBL>0 && checkEA1>0 && checkEA3>0) { std::for_each((clusEta[6].begin()+ietabin),(clusEta[6].begin()+ietabin+1),[](int &n){ n++; }); if (checkEA2>0) { std::for_each((clusHitEta[6].begin()+ietabin),(clusHitEta[6].begin()+ietabin+1),[](int &n){ n++; }); }}
407  if (checkIBL>0 && checkEA1>0 && checkEA2>0) { std::for_each((clusEta[7].begin()+ietabin),(clusEta[7].begin()+ietabin+1),[](int &n){ n++; }); if (checkEA3>0) { std::for_each((clusHitEta[7].begin()+ietabin),(clusHitEta[7].begin()+ietabin+1),[](int &n){ n++; }); }}
408  if (checkIBL>0 && checkEC2>0 && checkEC3>0) { std::for_each((clusEta[8].begin()+ietabin),(clusEta[8].begin()+ietabin+1),[](int &n){ n++; }); if (checkEC1>0) { std::for_each((clusHitEta[8].begin()+ietabin),(clusHitEta[8].begin()+ietabin+1),[](int &n){ n++; }); }}
409  if (checkIBL>0 && checkEC1>0 && checkEC3>0) { std::for_each((clusEta[9].begin()+ietabin),(clusEta[9].begin()+ietabin+1),[](int &n){ n++; }); if (checkEC2>0) { std::for_each((clusHitEta[9].begin()+ietabin),(clusHitEta[9].begin()+ietabin+1),[](int &n){ n++; }); }}
410  if (checkIBL>0 && checkEC1>0 && checkEC2>0) { std::for_each((clusEta[10].begin()+ietabin),(clusEta[10].begin()+ietabin+1),[](int &n){ n++; }); if (checkEC3>0) { std::for_each((clusHitEta[10].begin()+ietabin),(clusHitEta[10].begin()+ietabin+1),[](int &n){ n++; }); }}
411 
412  for (int i=0; i<(int)holeIndex.size(); i++) {
413  int becH, layerH, etaH, phiH;
414  GetLayerEtaPhiFromId(holeIndex.at(i),&becH,&layerH,&etaH,&phiH);
415  if (becH== 0 && layerH==0) {
416  if (etaH>-7 && etaH<6) { std::for_each((clusHole[0].begin()+ietabin),(clusHole[0].begin()+ietabin+1),[](int &n){ n++; }); }
417  else { std::for_each((clusHole[1].begin()+ietabin),(clusHole[1].begin()+ietabin+1),[](int &n){ n++; }); }
418  }
419  if (becH== 0 && layerH==1) { std::for_each((clusHole[2].begin()+ietabin),(clusHole[2].begin()+ietabin+1),[](int &n){ n++; }); }
420  if (becH== 0 && layerH==2) { std::for_each((clusHole[3].begin()+ietabin),(clusHole[3].begin()+ietabin+1),[](int &n){ n++; }); }
421  if (becH== 0 && layerH==3) { std::for_each((clusHole[4].begin()+ietabin),(clusHole[4].begin()+ietabin+1),[](int &n){ n++; }); }
422  if (becH== 2 && layerH==0) { std::for_each((clusHole[5].begin()+ietabin),(clusHole[5].begin()+ietabin+1),[](int &n){ n++; }); }
423  if (becH== 2 && layerH==1) { std::for_each((clusHole[6].begin()+ietabin),(clusHole[6].begin()+ietabin+1),[](int &n){ n++; }); }
424  if (becH== 2 && layerH==2) { std::for_each((clusHole[7].begin()+ietabin),(clusHole[7].begin()+ietabin+1),[](int &n){ n++; }); }
425  if (becH==-2 && layerH==0) { std::for_each((clusHole[8].begin()+ietabin),(clusHole[8].begin()+ietabin+1),[](int &n){ n++; }); }
426  if (becH==-2 && layerH==1) { std::for_each((clusHole[9].begin()+ietabin),(clusHole[9].begin()+ietabin+1),[](int &n){ n++; }); }
427  if (becH==-2 && layerH==2) { std::for_each((clusHole[10].begin()+ietabin),(clusHole[10].begin()+ietabin+1),[](int &n){ n++; }); }
428  }
429  }
430 
431  // Cluster should not contain edge pixels
432  int clusterEdge = 0;
433  for (int i=0; i<(int)clusterLayer.size(); i++) { clusterEdge+=isEdge.at(i); }
434 
435  // Cluster should not contain overflow pixels
436  int clusterOverflow = 0;
437  for (int i=0; i<(int)clusterLayer.size(); i++) { clusterOverflow+=isOverflow.at(i); }
438 
439  // Cluster should not contain split state
440  int isSplit = 0;
441  for (int i=0; i<(int)clusterLayer.size(); i++) { isSplit+=clusterIsSplit.at(i); }
442 
443  // Strong isolation
444  int iso20x4 = 0;
445  for (int i=0; i<(int)clusterLayer.size(); i++) { iso20x4+=clusterIsolation20x4.at(i); }
446 
447  // Good tracks must be required for the dE/dx and Lorentz angle measurements
448  bool passCut = false;
449  if (numberOfPixelHits>3 && numberOfPixelHoles==0 && nPixelDeadSensor==0 && numberOfPixelOutliers==0 && clusterEdge==0 && clusterOverflow==0 && isSplit==0 && iso20x4==0) {
450  if (checkIBL>0 && checkBLY>0 && checkLY1>0 && checkLY2>0) { // Barrel
451  passCut = true;
452  }
453  else if (checkIBL>0 && checkEA1>0 && checkEA2>0 && checkEA3>0) { // Endcap
454  passCut = true;
455  }
456  else if (checkIBL>0 && checkEC1>0 && checkEC2>0 && checkEC3>0) { // Endcap
457  passCut = true;
458  }
459  else if (check3D>0 && checkEC2>0 && checkEC3>0) { // 3D
460  passCut = true;
461  }
462  }
463  if (!passCut) { continue; }
464 
465  // Cut on angle alpha
466  bool isAlphaCut = false;
467  for (int i=0; i<(int)clusterLayer.size(); i++) {
468  float alpha = std::atan(std::hypot(std::tan(trackTheta[i]),std::tan(trackPhi[i])));
469  if (std::cos(alpha)<0.16) { isAlphaCut=true; break; }
470  }
471  if (isAlphaCut) { continue; }
472 
473  if (trk->pt()>maxPt) {
474  maxPt = trk->pt();
475  trk_maxpt = trk;
476  }
477 
478  // Cluster study
479  float energyPair = 3.62e-6; // Electron-hole pair creation energy e-h/[MeV]
480  float siDensity = 2.329; // Silicon density [g/cm^3]
481  float thicknessIBL = 0.0200; // thickness for IBL planar [cm]
482  float thickness3D = 0.0230; // thickness for IBL 3D [cm]
483  float thicknessPIX = 0.0250; // thickness for PIXEL [cm]
484  float sumE = 0.0;
485  float sumX = 0.0;
486  for (int i=0; i<(int)clusterLayer.size(); i++) {
487 
488  float alpha = std::atan(std::hypot(std::tan(trackTheta[i]),std::tan(trackPhi[i])));
489  float thickness = thicknessPIX;
490  if (clusterBEC.at(i)==0 && clusterLayer.at(i)==0) {
491  thickness = (clusterModuleEta[i]<6 && clusterModuleEta[i]>-7) ? thicknessIBL : thickness3D;
492  }
493  float recodEdx = clusterCharge[i]*energyPair/thickness/siDensity*std::cos(alpha);
494  sumE += clusterCharge[i]*energyPair;
495  sumX += thickness*siDensity/std::cos(alpha);
496 
497  int ibc = (clusterCharge[i]>0.0 && clusterCharge[i]<120000.0) ? trunc(clusterCharge[i]/2000.0) : nbin_charge-1;
498  int ibt = (clusterToT[i]>0 && clusterToT[i]<120) ? trunc(clusterToT[i]) : nbin_tot-1;
499  int ibe = (recodEdx>0.0 && recodEdx<3.0) ? trunc(recodEdx*100) : nbin_dedx-1;
500  int ibx = (clusterSizePhi[i]>-0.5 && clusterSizePhi[i]<9.5) ? trunc(clusterSizePhi[i]) : nbin_size-1;
501  int ibz = (clusterSizeZ[i]>-0.5 && clusterSizeZ[i]<9.5) ? trunc(clusterSizeZ[i]) : nbin_size-1;
502  int irx = (unbiasedResidualX[i]>-0.4 && unbiasedResidualX[i]<0.4) ? trunc((unbiasedResidualX[i]+0.4)/0.01) : nbin_reso-1;
503  int iry = (unbiasedResidualY[i]>-0.4 && unbiasedResidualY[i]<0.4) ? trunc((unbiasedResidualY[i]+0.4)/0.01) : nbin_reso-1;
504  int ibo = (numTotalPixelsPerModule[i]>-0.1 && numTotalPixelsPerModule[i]<160) ? trunc(numTotalPixelsPerModule[i]/2.0) : nbin_occ-1;
505 
506  int idx = 0;
507  if (clusterBEC.at(i)==0 && clusterLayer.at(i)==0) {
508  if (clusterModuleEta[i]<6 && clusterModuleEta[i]>-7) { idx=0; } // IBL-planar
509  else { idx=1; } // 3D sensor
510  }
511  else if (clusterBEC.at(i)==0 && clusterLayer.at(i)==1) { idx=2; } // B-layer
512  else if (clusterBEC.at(i)==0 && clusterLayer.at(i)==2) { idx=3; } // Layer-1
513  else if (clusterBEC.at(i)==0 && clusterLayer.at(i)==3) { idx=4; } // Layer-2
514  else if (clusterBEC.at(i)==2 && clusterLayer.at(i)==0) { idx=5; } // Endcap-A1
515  else if (clusterBEC.at(i)==2 && clusterLayer.at(i)==1) { idx=6; } // Endcap-A2
516  else if (clusterBEC.at(i)==2 && clusterLayer.at(i)==2) { idx=7; } // Endcap-A3
517  else if (clusterBEC.at(i)==-2 && clusterLayer.at(i)==0) { idx=8; } // Endcap-C1
518  else if (clusterBEC.at(i)==-2 && clusterLayer.at(i)==1) { idx=9; } // Endcap-C2
519  else if (clusterBEC.at(i)==-2 && clusterLayer.at(i)==2) { idx=10; } // Endcap-C3
520 
521  std::for_each((clusCharge[idx].begin()+ibc),(clusCharge[idx].begin()+ibc+1),[](int &n){ n++; });
522  std::for_each((clusToT[idx].begin()+ibt),(clusToT[idx].begin()+ibt+1),[](int &n){ n++; });
523  std::for_each((clusdEdx[idx].begin()+ibe),(clusdEdx[idx].begin()+ibe+1),[](int &n){ n++; });
524  std::for_each((clusSizeX[idx].begin()+ibx),(clusSizeX[idx].begin()+ibx+1),[](int &n){ n++; });
525  std::for_each((clusSizeZ[idx].begin()+ibz),(clusSizeZ[idx].begin()+ibz+1),[](int &n){ n++; });
526  std::for_each((clusResidualX[idx].begin()+irx),(clusResidualX[idx].begin()+irx+1),[](int &n){ n++; });
527  std::for_each((clusResidualY[idx].begin()+iry),(clusResidualY[idx].begin()+iry+1),[](int &n){ n++; });
528  std::for_each((clusOcc[idx].begin()+ibo),(clusOcc[idx].begin()+ibo+1),[](int &n){ n++; });
529  }
530 
531  if (sumX>0.0) {
532  int iby = (sumE/sumX>0.0 && sumE/sumX<3.0) ? trunc(sumE/sumX*100) : 299;
533  std::for_each((trkdEdx.begin()+iby),(trkdEdx.begin()+iby+1),[](int &n){ n++; });
534  }
535  }
536 
537  if (m_storeMode==2) {
538  if (maxPt>0.0) {
539  xAOD::TrackParticle* tp = new xAOD::TrackParticle();
540  tp->makePrivateStore(*trk_maxpt);
541  tp->setDefiningParametersCovMatrixVec(tmpCov);
542 
543  static const SG::Decorator<std::vector<int>> TrackEtaIBL("TrackEtaIBL");
544  static const SG::Decorator<std::vector<int>> TrackEta3D("TrackEta3D");
545  static const SG::Decorator<std::vector<int>> TrackEtaBL("TrackEtaBL");
546  static const SG::Decorator<std::vector<int>> TrackEtaL1("TrackEtaL1");
547  static const SG::Decorator<std::vector<int>> TrackEtaL2("TrackEtaL2");
548  static const SG::Decorator<std::vector<int>> TrackEtaEA1("TrackEtaEA1");
549  static const SG::Decorator<std::vector<int>> TrackEtaEA2("TrackEtaEA2");
550  static const SG::Decorator<std::vector<int>> TrackEtaEA3("TrackEtaEA3");
551  static const SG::Decorator<std::vector<int>> TrackEtaEC1("TrackEtaEC1");
552  static const SG::Decorator<std::vector<int>> TrackEtaEC2("TrackEtaEC2");
553  static const SG::Decorator<std::vector<int>> TrackEtaEC3("TrackEtaEC3");
554  static const SG::Decorator<std::vector<int>> TrackHitEtaIBL("TrackHitEtaIBL");
555  static const SG::Decorator<std::vector<int>> TrackHitEta3D("TrackHitEta3D");
556  static const SG::Decorator<std::vector<int>> TrackHitEtaBL("TrackHitEtaBL");
557  static const SG::Decorator<std::vector<int>> TrackHitEtaL1("TrackHitEtaL1");
558  static const SG::Decorator<std::vector<int>> TrackHitEtaL2("TrackHitEtaL2");
559  static const SG::Decorator<std::vector<int>> TrackHitEtaEA1("TrackHitEtaEA1");
560  static const SG::Decorator<std::vector<int>> TrackHitEtaEA2("TrackHitEtaEA2");
561  static const SG::Decorator<std::vector<int>> TrackHitEtaEA3("TrackHitEtaEA3");
562  static const SG::Decorator<std::vector<int>> TrackHitEtaEC1("TrackHitEtaEC1");
563  static const SG::Decorator<std::vector<int>> TrackHitEtaEC2("TrackHitEtaEC2");
564  static const SG::Decorator<std::vector<int>> TrackHitEtaEC3("TrackHitEtaEC3");
565  static const SG::Decorator<std::vector<int>> ClusterChargeIBL("ClusterChargeIBL");
566  static const SG::Decorator<std::vector<int>> ClusterCharge3D("ClusterCharge3D");
567  static const SG::Decorator<std::vector<int>> ClusterChargeBL("ClusterChargeBL");
568  static const SG::Decorator<std::vector<int>> ClusterChargeL1("ClusterChargeL1");
569  static const SG::Decorator<std::vector<int>> ClusterChargeL2("ClusterChargeL2");
570  static const SG::Decorator<std::vector<int>> ClusterChargeEA1("ClusterChargeEA1");
571  static const SG::Decorator<std::vector<int>> ClusterChargeEA2("ClusterChargeEA2");
572  static const SG::Decorator<std::vector<int>> ClusterChargeEA3("ClusterChargeEA3");
573  static const SG::Decorator<std::vector<int>> ClusterChargeEC1("ClusterChargeEC1");
574  static const SG::Decorator<std::vector<int>> ClusterChargeEC2("ClusterChargeEC2");
575  static const SG::Decorator<std::vector<int>> ClusterChargeEC3("ClusterChargeEC3");
576  static const SG::Decorator<std::vector<int>> ClusterToTIBL("ClusterToTIBL");
577  static const SG::Decorator<std::vector<int>> ClusterToT3D("ClusterToT3D");
578  static const SG::Decorator<std::vector<int>> ClusterToTBL("ClusterToTBL");
579  static const SG::Decorator<std::vector<int>> ClusterToTL1("ClusterToTL1");
580  static const SG::Decorator<std::vector<int>> ClusterToTL2("ClusterToTL2");
581  static const SG::Decorator<std::vector<int>> ClusterToTEA1("ClusterToTEA1");
582  static const SG::Decorator<std::vector<int>> ClusterToTEA2("ClusterToTEA2");
583  static const SG::Decorator<std::vector<int>> ClusterToTEA3("ClusterToTEA3");
584  static const SG::Decorator<std::vector<int>> ClusterToTEC1("ClusterToTEC1");
585  static const SG::Decorator<std::vector<int>> ClusterToTEC2("ClusterToTEC2");
586  static const SG::Decorator<std::vector<int>> ClusterToTEC3("ClusterToTEC3");
587  static const SG::Decorator<std::vector<int>> ClusterdEdxIBL("ClusterdEdxIBL");
588  static const SG::Decorator<std::vector<int>> ClusterdEdx3D("ClusterdEdx3D");
589  static const SG::Decorator<std::vector<int>> ClusterdEdxBL("ClusterdEdxBL");
590  static const SG::Decorator<std::vector<int>> ClusterdEdxL1("ClusterdEdxL1");
591  static const SG::Decorator<std::vector<int>> ClusterdEdxL2("ClusterdEdxL2");
592  static const SG::Decorator<std::vector<int>> ClusterdEdxEA1("ClusterdEdxEA1");
593  static const SG::Decorator<std::vector<int>> ClusterdEdxEA2("ClusterdEdxEA2");
594  static const SG::Decorator<std::vector<int>> ClusterdEdxEA3("ClusterdEdxEA3");
595  static const SG::Decorator<std::vector<int>> ClusterdEdxEC1("ClusterdEdxEC1");
596  static const SG::Decorator<std::vector<int>> ClusterdEdxEC2("ClusterdEdxEC2");
597  static const SG::Decorator<std::vector<int>> ClusterdEdxEC3("ClusterdEdxEC3");
598  static const SG::Decorator<std::vector<int>> ClusterSizeXIBL("ClusterSizeXIBL");
599  static const SG::Decorator<std::vector<int>> ClusterSizeX3D("ClusterSizeX3D");
600  static const SG::Decorator<std::vector<int>> ClusterSizeXBL("ClusterSizeXBL");
601  static const SG::Decorator<std::vector<int>> ClusterSizeXL1("ClusterSizeXL1");
602  static const SG::Decorator<std::vector<int>> ClusterSizeXL2("ClusterSizeXL2");
603  static const SG::Decorator<std::vector<int>> ClusterSizeXEA1("ClusterSizeXEA1");
604  static const SG::Decorator<std::vector<int>> ClusterSizeXEA2("ClusterSizeXEA2");
605  static const SG::Decorator<std::vector<int>> ClusterSizeXEA3("ClusterSizeXEA3");
606  static const SG::Decorator<std::vector<int>> ClusterSizeXEC1("ClusterSizeXEC1");
607  static const SG::Decorator<std::vector<int>> ClusterSizeXEC2("ClusterSizeXEC2");
608  static const SG::Decorator<std::vector<int>> ClusterSizeXEC3("ClusterSizeXEC3");
609  static const SG::Decorator<std::vector<int>> ClusterSizeZIBL("ClusterSizeZIBL");
610  static const SG::Decorator<std::vector<int>> ClusterSizeZ3D("ClusterSizeZ3D");
611  static const SG::Decorator<std::vector<int>> ClusterSizeZBL("ClusterSizeZBL");
612  static const SG::Decorator<std::vector<int>> ClusterSizeZL1("ClusterSizeZL1");
613  static const SG::Decorator<std::vector<int>> ClusterSizeZL2("ClusterSizeZL2");
614  static const SG::Decorator<std::vector<int>> ClusterSizeZEA1("ClusterSizeZEA1");
615  static const SG::Decorator<std::vector<int>> ClusterSizeZEA2("ClusterSizeZEA2");
616  static const SG::Decorator<std::vector<int>> ClusterSizeZEA3("ClusterSizeZEA3");
617  static const SG::Decorator<std::vector<int>> ClusterSizeZEC1("ClusterSizeZEC1");
618  static const SG::Decorator<std::vector<int>> ClusterSizeZEC2("ClusterSizeZEC2");
619  static const SG::Decorator<std::vector<int>> ClusterSizeZEC3("ClusterSizeZEC3");
620  static const SG::Decorator<std::vector<int>> ClusterResidualXIBL("ClusterResidualXIBL");
621  static const SG::Decorator<std::vector<int>> ClusterResidualX3D("ClusterResidualX3D");
622  static const SG::Decorator<std::vector<int>> ClusterResidualXBL("ClusterResidualXBL");
623  static const SG::Decorator<std::vector<int>> ClusterResidualXL1("ClusterResidualXL1");
624  static const SG::Decorator<std::vector<int>> ClusterResidualXL2("ClusterResidualXL2");
625  static const SG::Decorator<std::vector<int>> ClusterResidualXEA1("ClusterResidualXEA1");
626  static const SG::Decorator<std::vector<int>> ClusterResidualXEA2("ClusterResidualXEA2");
627  static const SG::Decorator<std::vector<int>> ClusterResidualXEA3("ClusterResidualXEA3");
628  static const SG::Decorator<std::vector<int>> ClusterResidualXEC1("ClusterResidualXEC1");
629  static const SG::Decorator<std::vector<int>> ClusterResidualXEC2("ClusterResidualXEC2");
630  static const SG::Decorator<std::vector<int>> ClusterResidualXEC3("ClusterResidualXEC3");
631  static const SG::Decorator<std::vector<int>> ClusterResidualYIBL("ClusterResidualYIBL");
632  static const SG::Decorator<std::vector<int>> ClusterResidualY3D("ClusterResidualY3D");
633  static const SG::Decorator<std::vector<int>> ClusterResidualYBL("ClusterResidualYBL");
634  static const SG::Decorator<std::vector<int>> ClusterResidualYL1("ClusterResidualYL1");
635  static const SG::Decorator<std::vector<int>> ClusterResidualYL2("ClusterResidualYL2");
636  static const SG::Decorator<std::vector<int>> ClusterResidualYEA1("ClusterResidualYEA1");
637  static const SG::Decorator<std::vector<int>> ClusterResidualYEA2("ClusterResidualYEA2");
638  static const SG::Decorator<std::vector<int>> ClusterResidualYEA3("ClusterResidualYEA3");
639  static const SG::Decorator<std::vector<int>> ClusterResidualYEC1("ClusterResidualYEC1");
640  static const SG::Decorator<std::vector<int>> ClusterResidualYEC2("ClusterResidualYEC2");
641  static const SG::Decorator<std::vector<int>> ClusterResidualYEC3("ClusterResidualYEC3");
642  static const SG::Decorator<std::vector<int>> ClusterHoleIBL("ClusterHoleIBL");
643  static const SG::Decorator<std::vector<int>> ClusterHole3D("ClusterHole3D");
644  static const SG::Decorator<std::vector<int>> ClusterHoleBL("ClusterHoleBL");
645  static const SG::Decorator<std::vector<int>> ClusterHoleL1("ClusterHoleL1");
646  static const SG::Decorator<std::vector<int>> ClusterHoleL2("ClusterHoleL2");
647  static const SG::Decorator<std::vector<int>> ClusterHoleEA1("ClusterHoleEA1");
648  static const SG::Decorator<std::vector<int>> ClusterHoleEA2("ClusterHoleEA2");
649  static const SG::Decorator<std::vector<int>> ClusterHoleEA3("ClusterHoleEA3");
650  static const SG::Decorator<std::vector<int>> ClusterHoleEC1("ClusterHoleEC1");
651  static const SG::Decorator<std::vector<int>> ClusterHoleEC2("ClusterHoleEC2");
652  static const SG::Decorator<std::vector<int>> ClusterHoleEC3("ClusterHoleEC3");
653  static const SG::Decorator<std::vector<int>> ClusterOccIBL("ClusterOccIBL");
654  static const SG::Decorator<std::vector<int>> ClusterOcc3D("ClusterOcc3D");
655  static const SG::Decorator<std::vector<int>> ClusterOccBL("ClusterOccBL");
656  static const SG::Decorator<std::vector<int>> ClusterOccL1("ClusterOccL1");
657  static const SG::Decorator<std::vector<int>> ClusterOccL2("ClusterOccL2");
658  static const SG::Decorator<std::vector<int>> ClusterOccEA1("ClusterOccEA1");
659  static const SG::Decorator<std::vector<int>> ClusterOccEA2("ClusterOccEA2");
660  static const SG::Decorator<std::vector<int>> ClusterOccEA3("ClusterOccEA3");
661  static const SG::Decorator<std::vector<int>> ClusterOccEC1("ClusterOccEC1");
662  static const SG::Decorator<std::vector<int>> ClusterOccEC2("ClusterOccEC2");
663  static const SG::Decorator<std::vector<int>> ClusterOccEC3("ClusterOccEC3");
664  static const SG::Decorator<std::vector<int>> TrackALL("TrackALL");
665  static const SG::Decorator<std::vector<int>> TrackHOLE("TrackHOLE");
666  static const SG::Decorator<std::vector<int>> TrackdEdx("TrackdEdx");
667 
668  TrackEtaIBL(*tp) = std::move(clusEta[0]);
669  TrackEta3D(*tp) = std::move(clusEta[1]);
670  TrackEtaBL(*tp) = std::move(clusEta[2]);
671  TrackEtaL1(*tp) = std::move(clusEta[3]);
672  TrackEtaL2(*tp) = std::move(clusEta[4]);
673  TrackEtaEA1(*tp) = std::move(clusEta[5]);
674  TrackEtaEA2(*tp) = std::move(clusEta[6]);
675  TrackEtaEA3(*tp) = std::move(clusEta[7]);
676  TrackEtaEC1(*tp) = std::move(clusEta[8]);
677  TrackEtaEC2(*tp) = std::move(clusEta[9]);
678  TrackEtaEC3(*tp) = std::move(clusEta[10]);
679  TrackHitEtaIBL(*tp) = std::move(clusHitEta[0]);
680  TrackHitEta3D(*tp) = std::move(clusHitEta[1]);
681  TrackHitEtaBL(*tp) = std::move(clusHitEta[2]);
682  TrackHitEtaL1(*tp) = std::move(clusHitEta[3]);
683  TrackHitEtaL2(*tp) = std::move(clusHitEta[4]);
684  TrackHitEtaEA1(*tp) = std::move(clusHitEta[5]);
685  TrackHitEtaEA2(*tp) = std::move(clusHitEta[6]);
686  TrackHitEtaEA3(*tp) = std::move(clusHitEta[7]);
687  TrackHitEtaEC1(*tp) = std::move(clusHitEta[8]);
688  TrackHitEtaEC2(*tp) = std::move(clusHitEta[9]);
689  TrackHitEtaEC3(*tp) = std::move(clusHitEta[10]);
690  ClusterChargeIBL(*tp) = std::move(clusCharge[0]);
691  ClusterCharge3D(*tp) = std::move(clusCharge[1]);
692  ClusterChargeBL(*tp) = std::move(clusCharge[2]);
693  ClusterChargeL1(*tp) = std::move(clusCharge[3]);
694  ClusterChargeL2(*tp) = std::move(clusCharge[4]);
695  ClusterChargeEA1(*tp) = std::move(clusCharge[5]);
696  ClusterChargeEA2(*tp) = std::move(clusCharge[6]);
697  ClusterChargeEA3(*tp) = std::move(clusCharge[7]);
698  ClusterChargeEC1(*tp) = std::move(clusCharge[8]);
699  ClusterChargeEC2(*tp) = std::move(clusCharge[9]);
700  ClusterChargeEC3(*tp) = std::move(clusCharge[10]);
701  ClusterToTIBL(*tp) = std::move(clusToT[0]);
702  ClusterToT3D(*tp) = std::move(clusToT[1]);
703  ClusterToTBL(*tp) = std::move(clusToT[2]);
704  ClusterToTL1(*tp) = std::move(clusToT[3]);
705  ClusterToTL2(*tp) = std::move(clusToT[4]);
706  ClusterToTEA1(*tp) = std::move(clusToT[5]);
707  ClusterToTEA2(*tp) = std::move(clusToT[6]);
708  ClusterToTEA3(*tp) = std::move(clusToT[7]);
709  ClusterToTEC1(*tp) = std::move(clusToT[8]);
710  ClusterToTEC2(*tp) = std::move(clusToT[9]);
711  ClusterToTEC3(*tp) = std::move(clusToT[10]);
712  ClusterdEdxIBL(*tp) = std::move(clusdEdx[0]);
713  ClusterdEdx3D(*tp) = std::move(clusdEdx[1]);
714  ClusterdEdxBL(*tp) = std::move(clusdEdx[2]);
715  ClusterdEdxL1(*tp) = std::move(clusdEdx[3]);
716  ClusterdEdxL2(*tp) = std::move(clusdEdx[4]);
717  ClusterdEdxEA1(*tp) = std::move(clusdEdx[5]);
718  ClusterdEdxEA2(*tp) = std::move(clusdEdx[6]);
719  ClusterdEdxEA3(*tp) = std::move(clusdEdx[7]);
720  ClusterdEdxEC1(*tp) = std::move(clusdEdx[8]);
721  ClusterdEdxEC2(*tp) = std::move(clusdEdx[9]);
722  ClusterdEdxEC3(*tp) = std::move(clusdEdx[10]);
723  ClusterSizeXIBL(*tp) = std::move(clusSizeX[0]);
724  ClusterSizeX3D(*tp) = std::move(clusSizeX[1]);
725  ClusterSizeXBL(*tp) = std::move(clusSizeX[2]);
726  ClusterSizeXL1(*tp) = std::move(clusSizeX[3]);
727  ClusterSizeXL2(*tp) = std::move(clusSizeX[4]);
728  ClusterSizeXEA1(*tp) = std::move(clusSizeX[5]);
729  ClusterSizeXEA2(*tp) = std::move(clusSizeX[6]);
730  ClusterSizeXEA3(*tp) = std::move(clusSizeX[7]);
731  ClusterSizeXEC1(*tp) = std::move(clusSizeX[8]);
732  ClusterSizeXEC2(*tp) = std::move(clusSizeX[9]);
733  ClusterSizeXEC3(*tp) = std::move(clusSizeX[10]);
734  ClusterSizeZIBL(*tp) = std::move(clusSizeZ[0]);
735  ClusterSizeZ3D(*tp) = std::move(clusSizeZ[1]);
736  ClusterSizeZBL(*tp) = std::move(clusSizeZ[2]);
737  ClusterSizeZL1(*tp) = std::move(clusSizeZ[3]);
738  ClusterSizeZL2(*tp) = std::move(clusSizeZ[4]);
739  ClusterSizeZEA1(*tp) = std::move(clusSizeZ[5]);
740  ClusterSizeZEA2(*tp) = std::move(clusSizeZ[6]);
741  ClusterSizeZEA3(*tp) = std::move(clusSizeZ[7]);
742  ClusterSizeZEC1(*tp) = std::move(clusSizeZ[8]);
743  ClusterSizeZEC2(*tp) = std::move(clusSizeZ[9]);
744  ClusterSizeZEC3(*tp) = std::move(clusSizeZ[10]);
745  ClusterResidualXIBL(*tp) = std::move(clusResidualX[0]);
746  ClusterResidualX3D(*tp) = std::move(clusResidualX[1]);
747  ClusterResidualXBL(*tp) = std::move(clusResidualX[2]);
748  ClusterResidualXL1(*tp) = std::move(clusResidualX[3]);
749  ClusterResidualXL2(*tp) = std::move(clusResidualX[4]);
750  ClusterResidualXEA1(*tp) = std::move(clusResidualX[5]);
751  ClusterResidualXEA2(*tp) = std::move(clusResidualX[6]);
752  ClusterResidualXEA3(*tp) = std::move(clusResidualX[7]);
753  ClusterResidualXEC1(*tp) = std::move(clusResidualX[8]);
754  ClusterResidualXEC2(*tp) = std::move(clusResidualX[9]);
755  ClusterResidualXEC3(*tp) = std::move(clusResidualX[10]);
756  ClusterResidualYIBL(*tp) = std::move(clusResidualY[0]);
757  ClusterResidualY3D(*tp) = std::move(clusResidualY[1]);
758  ClusterResidualYBL(*tp) = std::move(clusResidualY[2]);
759  ClusterResidualYL1(*tp) = std::move(clusResidualY[3]);
760  ClusterResidualYL2(*tp) = std::move(clusResidualY[4]);
761  ClusterResidualYEA1(*tp) = std::move(clusResidualY[5]);
762  ClusterResidualYEA2(*tp) = std::move(clusResidualY[6]);
763  ClusterResidualYEA3(*tp) = std::move(clusResidualY[7]);
764  ClusterResidualYEC1(*tp) = std::move(clusResidualY[8]);
765  ClusterResidualYEC2(*tp) = std::move(clusResidualY[9]);
766  ClusterResidualYEC3(*tp) = std::move(clusResidualY[10]);
767  ClusterHoleIBL(*tp) = std::move(clusHole[0]);
768  ClusterHole3D(*tp) = std::move(clusHole[1]);
769  ClusterHoleBL(*tp) = std::move(clusHole[2]);
770  ClusterHoleL1(*tp) = std::move(clusHole[3]);
771  ClusterHoleL2(*tp) = std::move(clusHole[4]);
772  ClusterHoleEA1(*tp) = std::move(clusHole[5]);
773  ClusterHoleEA2(*tp) = std::move(clusHole[6]);
774  ClusterHoleEA3(*tp) = std::move(clusHole[7]);
775  ClusterHoleEC1(*tp) = std::move(clusHole[8]);
776  ClusterHoleEC2(*tp) = std::move(clusHole[9]);
777  ClusterHoleEC3(*tp) = std::move(clusHole[10]);
778  ClusterOccIBL(*tp) = std::move(clusOcc[0]);
779  ClusterOcc3D(*tp) = std::move(clusOcc[1]);
780  ClusterOccBL(*tp) = std::move(clusOcc[2]);
781  ClusterOccL1(*tp) = std::move(clusOcc[3]);
782  ClusterOccL2(*tp) = std::move(clusOcc[4]);
783  ClusterOccEA1(*tp) = std::move(clusOcc[5]);
784  ClusterOccEA2(*tp) = std::move(clusOcc[6]);
785  ClusterOccEA3(*tp) = std::move(clusOcc[7]);
786  ClusterOccEC1(*tp) = std::move(clusOcc[8]);
787  ClusterOccEC2(*tp) = std::move(clusOcc[9]);
788  ClusterOccEC3(*tp) = std::move(clusOcc[10]);
789  TrackALL(*tp) = std::move(trkEta);
790  TrackHOLE(*tp) = std::move(trkHole);
791  TrackdEdx(*tp) = std::move(trkdEdx);
792 
793  PixelMonitoringTrack->push_back(tp);
794  }
795  }
796  return StatusCode::SUCCESS;
797 }
798 
799 void DerivationFramework::PixelNtupleMaker::GetLayerEtaPhiFromId(uint64_t id,int *barrelEC, int *layer, int *eta, int *phi) {
800  *eta =(id>>43) % 0x20 - 10;
801  *phi =(id>>43) % 0x800 / 0x20;
802  int layer_index = ((id>>43) / 0x800) & 0xF;
803 
804  //A possibility is to use a bit by bit AND: layer_index & 0xF
805  switch (layer_index) {
806  case 4:
807  *barrelEC=-2;
808  *layer = 0;
809  break;
810  case 5:
811  *barrelEC=-2;
812  *layer = 1;
813  break;
814  case 6:
815  *barrelEC=-2;
816  *layer = 2;
817  break;
818  case 8:
819  *barrelEC=0;
820  *layer=0;
821  break;
822  case 9:
823  *layer=1;
824  *barrelEC=0;
825  break;
826  case 10:
827  *layer=2;
828  *barrelEC=0;
829  break;
830  case 11:
831  *barrelEC=0;
832  *layer =3;
833  break;
834  case 12:
835  *barrelEC=2;
836  *layer=0;
837  break;
838  case 13:
839  *layer =1;
840  *barrelEC=2;
841  break;
842  case 14:
843  *layer=2;
844  *barrelEC=2;
845  break;
846  }
847 }
848 
ShallowCopy.h
xAOD::numberOfPixelHoles
@ numberOfPixelHoles
number of pixel layers on track with absence of hits [unit8_t].
Definition: TrackingPrimitives.h:261
Trk::TrackState::Pixel
@ Pixel
Definition: TrackStateDefs.h:28
python.PerfMonSerializer.p
def p
Definition: PerfMonSerializer.py:743
xAOD::numberOfSCTSharedHits
@ numberOfSCTSharedHits
number of SCT hits shared by several tracks [unit8_t].
Definition: TrackingPrimitives.h:272
phi
Scalar phi() const
phi method
Definition: AmgMatrixBasePlugin.h:64
xAOD::uint8_t
uint8_t
Definition: Muon_v1.cxx:575
xAOD::TrackStateValidation_v1::detType
char detType() const
Returns the detector type.
CaloCellPos2Ntuple.int
int
Definition: CaloCellPos2Ntuple.py:24
makeComparison.deltaY
int deltaY
Definition: makeComparison.py:44
eta
Scalar eta() const
pseudorapidity method
Definition: AmgMatrixBasePlugin.h:79
TrackStateDefs.h
SG::ReadHandle
Definition: StoreGate/StoreGate/ReadHandle.h:70
PlotCalibFromCool.begin
begin
Definition: PlotCalibFromCool.py:94
TRT::Hit::globalZ
@ globalZ
Definition: HitInfo.h:39
ParticleTest.tp
tp
Definition: ParticleTest.py:25
xAOD::TrackStateValidation_v1
Class describing a TrackStateValidation.
Definition: TrackStateValidation_v1.h:28
xAOD::TrackStateValidation_v1::detElementId
uint64_t detElementId() const
Returns the detector element identifier.
xAOD::numberOfPixelHits
@ numberOfPixelHits
these are the pixel hits, including the b-layer [unit8_t].
Definition: TrackingPrimitives.h:259
xAOD::TrackStateValidation_v1::trackMeasurementValidationLink
ElementLink< xAOD::TrackMeasurementValidationContainer > trackMeasurementValidationLink() const
SG::ConstAccessor
Helper class to provide constant type-safe access to aux data.
Definition: ConstAccessor.h:54
InDetSecVtxTruthMatchUtils::isSplit
bool isSplit(int matchInfo)
Definition: InDetSecVtxTruthMatchTool.h:56
read_hist_ntuple.t
t
Definition: read_hist_ntuple.py:5
drawFromPickle.cos
cos
Definition: drawFromPickle.py:36
xAOD::TrackStateValidation_v1::localX
float localX() const
DerivationFramework::PixelNtupleMaker::~PixelNtupleMaker
~PixelNtupleMaker()
drawFromPickle.atan
atan
Definition: drawFromPickle.py:36
DerivationFramework::PixelNtupleMaker::PixelNtupleMaker
PixelNtupleMaker(const std::string &t, const std::string &n, const IInterface *p)
Definition: PixelNtupleMaker.cxx:18
xAOD::TrackParticle
TrackParticle_v1 TrackParticle
Reference the current persistent version:
Definition: Event/xAOD/xAODTracking/xAODTracking/TrackParticle.h:13
xAOD::TrackStateValidation_v1::unbiasedResidualY
float unbiasedResidualY() const
Returns the unbiased x residual.
xAOD::TrackMeasurementValidation_v1
Class describing a TrackMeasurementValidation.
Definition: TrackMeasurementValidation_v1.h:27
xAOD::numberOfPixelSharedHits
@ numberOfPixelSharedHits
number of Pixel all-layer hits shared by several tracks [unit8_t].
Definition: TrackingPrimitives.h:262
SG::Decorator< float >
Trk::TrackStateOnSurface::Hole
@ Hole
A hole on the track - this is defined in the following way.
Definition: TrackStateOnSurface.h:128
lumiFormat.i
int i
Definition: lumiFormat.py:92
beamspotman.n
n
Definition: beamspotman.py:731
EL::StatusCode
::StatusCode StatusCode
StatusCode definition for legacy code.
Definition: PhysicsAnalysis/D3PDTools/EventLoop/EventLoop/StatusCode.h:22
TRT::Hit::layer
@ layer
Definition: HitInfo.h:79
DerivationFramework::PixelNtupleMaker::MeasurementsOnTrack
std::vector< ElementLink< xAOD::TrackStateValidationContainer > > MeasurementsOnTrack
Definition: PixelNtupleMaker.h:55
xAOD::uint64_t
uint64_t
Definition: EventInfo_v1.cxx:123
xAOD::numberOfPixelOutliers
@ numberOfPixelOutliers
these are the pixel outliers, including the b-layer [unit8_t].
Definition: TrackingPrimitives.h:260
ATH_CHECK
#define ATH_CHECK
Definition: AthCheckMacros.h:40
drawFromPickle.tan
tan
Definition: drawFromPickle.py:36
xAOD::numberOfSCTHoles
@ numberOfSCTHoles
number of SCT holes [unit8_t].
Definition: TrackingPrimitives.h:270
compareGeometries.deltaX
float deltaX
Definition: compareGeometries.py:32
xAOD::TrackStateValidation_v1::localY
float localY() const
Returns the y position.
LB_AnalMapSplitter.tot
tot
Definition: LB_AnalMapSplitter.py:46
SG::ReadHandle::isValid
virtual bool isValid() override final
Can the handle be successfully dereferenced?
DerivationFramework::PixelNtupleMaker::finalize
StatusCode finalize()
Definition: PixelNtupleMaker.cxx:34
charge
double charge(const T &p)
Definition: AtlasPID.h:494
TauGNNUtils::Variables::Track::nSCTHits
bool nSCTHits(const xAOD::TauJet &, const xAOD::TauTrack &track, double &out)
Definition: TauGNNUtils.cxx:549
DataVector::push_back
value_type push_back(value_type pElem)
Add an element to the end of the collection.
xAOD::TrackParticle_v1::e
virtual double e() const override final
The total energy of the particle.
Definition: TrackParticle_v1.cxx:109
TrackParticle.h
SG::WriteHandle
Definition: StoreGate/StoreGate/WriteHandle.h:76
xAOD::TrackStateValidation_v1::localPhi
float localPhi() const
Returns the phi position.
PixelNtupleMaker.h
SG::WriteHandle::record
StatusCode record(std::unique_ptr< T > data)
Record a const object to the store.
xAOD::TrackStateValidation_v1::unbiasedResidualX
float unbiasedResidualX() const
Returns the unbiased x residual.
DerivationFramework::PixelNtupleMaker::addBranches
virtual StatusCode addBranches() const
Pass the thinning service
Definition: PixelNtupleMaker.cxx:39
LArNewCalib_DelayDump_OFC_Cali.idx
idx
Definition: LArNewCalib_DelayDump_OFC_Cali.py:69
xAOD::numberOfSCTHits
@ numberOfSCTHits
number of hits in SCT [unit8_t].
Definition: TrackingPrimitives.h:268
SG::ConstAccessor::isAvailable
bool isAvailable(const ELT &e) const
Test to see if this variable exists in the store.
TauGNNUtils::Variables::Track::trackPhi
bool trackPhi(const xAOD::TauJet &, const xAOD::TauTrack &track, double &out)
Definition: TauGNNUtils.cxx:480
xAOD::numberOfPixelDeadSensors
@ numberOfPixelDeadSensors
number of dead pixel sensors crossed [unit8_t].
Definition: TrackingPrimitives.h:266
xAOD::TrackParticle_v1
Class describing a TrackParticle.
Definition: TrackParticle_v1.h:43
DerivationFramework::PixelNtupleMaker::GetLayerEtaPhiFromId
static void GetLayerEtaPhiFromId(uint64_t id, int *barrelEC, int *layer, int *eta, int *phi)
Definition: PixelNtupleMaker.cxx:799
ConstAccessor.h
Helper class to provide constant type-safe access to aux data.
DerivationFramework::PixelNtupleMaker::initialize
StatusCode initialize()
Definition: PixelNtupleMaker.cxx:26
AthAlgTool
Definition: AthAlgTool.h:26
xAOD::TrackStateValidation_v1::type
int type() const
xAOD::TrackStateValidation_v1::localTheta
float localTheta() const
Returns the theta position.
TrackStateOnSurface.h
Generated on Sun Jun 30 2024 21:24:32 for ATLAS Offline Software by   doxygen 1.8.18