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