ATLAS Offline Software
TrackClusterAssValidation.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration
3 */
4 
5 #include "GaudiKernel/IPartPropSvc.h"
6 #include "GaudiKernel/ServiceHandle.h"
12 #include "StoreGate/ReadHandle.h"
13 #include "HepPDT/ParticleDataTable.hh"
14 #include "AtlasHepMC/GenVertex.h"
15 #include "AtlasHepMC/GenParticle.h"
16 
17 #include <cmath>
18 
20 // Constructor
22 
24 (const std::string& name,ISvcLocator* pSvcLocator)
25  : AthReentrantAlgorithm(name,pSvcLocator),
26  m_spacepointsSCTname("SCT_SpacePoints"),
27  m_spacepointsPixelname("PixelSpacePoints"),
28  m_spacepointsOverlapname("OverlapSpacePoints"),
29  m_clustersSCTname("SCT_Clusters"),
30  m_clustersPixelname("PixelClusters"),
31  m_clustersTRTname("TRT_DriftCircles"),
32  m_truth_locationPixel( "PRD_MultiTruthPixel" ),
33  m_truth_locationSCT( "PRD_MultiTruthSCT" ),
34  m_truth_locationTRT( "PRD_MultiTruthTRT" )
35 {
36 
37  // TrackClusterAssValidation steering parameters
38  //
39  m_ptcut = 1000. ;
40  m_ptcutmax = 1.e20 ;
41  m_rapcut = 2.6 ;
42  m_rmin = 0. ;
43  m_rmax = 30. ;
44  m_clcut = 6 ;
45  m_clcutTRT = 0 ;
46  m_spcut = 3 ;
47  m_usePIX = true ;
48  m_useSCT = true ;
49  m_useTRT = true ;
50  m_useOutliers = false ;
51  m_pdg = 0 ;
52  m_tcut = 0. ;
53  m_particleDataTable = nullptr ;
54 
55  declareProperty("TracksLocation" ,m_tracklocation );
56  declareProperty("SpacePointsSCTName" ,m_spacepointsSCTname );
57  declareProperty("SpacePointsPixelName" ,m_spacepointsPixelname );
58  declareProperty("SpacePointsOverlapName",m_spacepointsOverlapname);
59  declareProperty("PixelClusterContainer" ,m_clustersPixelname );
60  declareProperty("SCT_ClusterContainer" ,m_clustersSCTname );
61  declareProperty("TRT_ClusterContainer" ,m_clustersTRTname );
62  declareProperty("TruthLocationSCT" ,m_truth_locationSCT );
63  declareProperty("TruthLocationPixel" ,m_truth_locationPixel );
64  declareProperty("TruthLocationTRT" ,m_truth_locationTRT );
65  declareProperty("MomentumCut" ,m_ptcut );
66  declareProperty("MomentumMaxCut" ,m_ptcutmax );
67  declareProperty("RapidityCut" ,m_rapcut );
68  declareProperty("RadiusMin" ,m_rmin );
69  declareProperty("RadiusMax" ,m_rmax );
70  declareProperty("MinNumberClusters" ,m_clcut );
71  declareProperty("MinNumberClustersTRT" ,m_clcutTRT );
72  declareProperty("MinNumberSpacePoints" ,m_spcut );
73  declareProperty("usePixel" ,m_usePIX );
74  declareProperty("useSCT" ,m_useSCT );
75  declareProperty("useTRT" ,m_useTRT );
76  declareProperty("useOutliers" ,m_useOutliers );
77  declareProperty("pdgParticle" ,m_pdg );
78 }
79 
81 
82 
84 // Initialisation
86 
88 {
89 
90  StatusCode sc;
91 
92  if (m_rapcut == 0) {
93  m_tcut = 0;
94  }
95  else {
96  double den = tan(2.*atan(exp(-m_rapcut)));
97  if (den > 0) {
98  m_tcut = 1./den;
99  }
100  else {
102  }
103  }
104 
105  // get the Particle Properties Service
106  //
107  IPartPropSvc* partPropSvc = nullptr;
108  sc = service("PartPropSvc", partPropSvc, true);
109  if (sc.isFailure()) {
110  msg(MSG::FATAL) << " Could not initialize Particle Properties Service" << endmsg;
111  return StatusCode::FAILURE;
112  }
113 
114  // Particle Data Table
115  //
116  m_particleDataTable = partPropSvc->PDT();
117  if(!m_particleDataTable) {
118  msg(MSG::FATAL) << " Could not initialize Particle Properties Service" << endmsg;
119  return StatusCode::FAILURE;
120  }
121 
122  // Erase statistic information
123  //
124  m_pdg = std::abs(m_pdg) ;
125 
126  m_trackCollectionStat.resize(m_tracklocation.size());
127  m_eventStat = EventStat_t();
128 
129  if(!m_useTRT) m_clcutTRT = 0;
130  if(!m_clcutTRT) m_useTRT = false;
131 
132  // Read Handle Key
136 
140 
142 
145 
146  ATH_CHECK( m_tracklocation.initialize());
147 
148  // Read Cond Handle Key
151 
152  if (msgLvl(MSG::DEBUG)) {
154  msg() << endmsg;
155  }
156 
157  return sc;
158 }
159 
161 // Execute
163 
165 {
166 
167  if(!m_usePIX && !m_useSCT) return StatusCode::SUCCESS;
168  EventData_t event_data(m_tracklocation.size() );
169 
170  std::vector<SG::ReadHandle<PRD_MultiTruthCollection> > read_handle;
171  read_handle.reserve(3);
172  if(m_usePIX) {
173  read_handle.emplace_back(m_truth_locationPixel,ctx);
174  if (not read_handle.back().isValid()) {
175  ATH_MSG_FATAL( "Could not find TruthPIX" );
176  return StatusCode::FAILURE;
177  }
178  event_data.m_truthPIX = &(*read_handle.back());
179  }
180 
181  if(m_useSCT) {
182  read_handle.emplace_back(m_truth_locationSCT,ctx);
183  if (not read_handle.back().isValid()) {
184  ATH_MSG_FATAL( "Could not find TruthSCT" );
185  return StatusCode::FAILURE;
186  }
187  event_data.m_truthSCT = &(*read_handle.back());
188  }
189 
190  if(m_clcutTRT > 0) {
191  read_handle.emplace_back(m_truth_locationTRT,ctx);
192  if (not read_handle.back().isValid()) {
193  ATH_MSG_FATAL( "Could not find TruthTRT" );
194  return StatusCode::FAILURE;
195  }
196  event_data.m_truthTRT = &(*read_handle.back());
197  }
198 
199  newClustersEvent (ctx,event_data);
200  newSpacePointsEvent (ctx,event_data);
201  event_data.m_nqtracks = qualityTracksSelection(event_data);
202  tracksComparison (ctx,event_data);
203  if(!event_data.m_particles[0].empty()) {
204 
205  efficiencyReconstruction(event_data);
206  if(msgLvl(MSG::DEBUG)) noReconstructedParticles(event_data);
207 
208  }
209 
210  {
211  std::lock_guard<std::mutex> lock(m_statMutex);
212  assert( event_data.m_trackCollectionStat.size() == m_trackCollectionStat.size());
213  for (unsigned int i=0; i< m_trackCollectionStat.size(); ++i ) {
214  m_trackCollectionStat[i] += event_data.m_trackCollectionStat[i];
215  }
216  m_eventStat += event_data.m_eventStat;
217  }
218 
219  if (msgLvl(MSG::DEBUG)) {
220  dumpevent(msg(),event_data);
221  msg() << endmsg;
222  }
223  return StatusCode::SUCCESS;
224 }
225 
227 // Finalize
229 
231  if(m_eventStat.m_events<=0) return StatusCode::SUCCESS;
232  const auto & w13 = std::setw(13);
233  const auto & p5 = std::setprecision(5);
234  const auto topNtail=[](const std::string & str){return "|" + str + "|";};
235  const std::string lineSeparator(83,'-');
236  const std::string spaceSeparator(83,' ');
237  std::stringstream out;
238  out<<std::fixed;
239  out<<topNtail(lineSeparator)<<"\n";
240  out<<"| TrackClusterAssValidation statistic for charge truth particles with |\n";
241  out<<topNtail(spaceSeparator)<<"\n";
242  out<<"| pT >="<<w13<<p5<<m_ptcut<<" MeV"<<" |\n";
243  if(m_ptcutmax < 1000000.) {
244  out<<"| pT <="<<w13<<p5<<m_ptcutmax<<" MeV"<<" |\n";
245  }
246  out<<"| |rapidity| <="<<w13<<p5<<m_rapcut<<" |\n";
247  out<<"| max vertex radius <="<<w13<<p5<<m_rmax<<" mm"<<" |\n";
248  out<<"| min vertex radius >="<<w13<<p5<<m_rmin<<" mm"<<" |\n";
249  out<<"| particles pdg ="<<std::setw(8)<<m_pdg<<" |\n";
250 
251  auto yesNo=[](const bool predicate){return predicate ? "yes" : "no ";};
252  out<<"| use Pixels information "<<yesNo(m_usePIX)<<" |\n";
253  out<<"| use SCT information "<<yesNo(m_useSCT)<<" |\n";
254  out<<"| use TRT information "<<yesNo(m_useTRT)<<" |\n";
255  out<<"| take into account outliers "<<yesNo(m_useOutliers)<<" |\n";
256  out<<topNtail(spaceSeparator)<<"\n";
257  if(!m_usePIX && !m_useSCT) return StatusCode::SUCCESS;
258  enum Regions{Barrel, Transition, Endcap, Forward, NRegions};
259  auto incrementArray=[](auto & array, const int idx){for (int j{};j!=NRegions;++j) array[idx][j] += array[idx+1][j];};
260  for(int i=48; i>=0; --i) {
261  m_eventStat.m_particleClusters [i] +=m_eventStat.m_particleClusters [i+1];
262  incrementArray(m_eventStat.m_particleClustersBTE, i);
263  m_eventStat.m_particleSpacePoints [i] +=m_eventStat.m_particleSpacePoints [i+1];
264  incrementArray(m_eventStat.m_particleSpacePointsBTE, i);
265  }
266  auto coerceToOne=[](const double & initialVal)->double {return (initialVal<1.) ? 1. : initialVal; };
267  /* Note that in all cases below, the dimension of the target (i.e. result) array is 10
268  * whereas the dimension of the source, or input, array is 50. The first index of the source
269  * is used for totals, to act as a denominator(coerced to 1 if necessary). Bin 49 is used for overflows.
270  * Stats for single clusters (index 1) appear to be unused in printout.
271  */
272  //all
273  double pa = coerceToOne(m_eventStat.m_particleClusters[0]);
274  std::array<double, 10> pc2ff{};
275  size_t clusterIdx = 2;
276  for (auto & thisCluster: pc2ff){
277  thisCluster = double(m_eventStat.m_particleClusters[ clusterIdx++ ])/ pa;
278  }
279  //barrel
280  pa = coerceToOne(m_eventStat.m_particleClustersBTE[0][Barrel]);
281  std::array<double, 10> pcBarrel2ff{};
282  size_t clusterBarrelIdx = 2;
283  for (auto & thisClusterB: pcBarrel2ff){
284  thisClusterB = double(m_eventStat.m_particleClustersBTE[ clusterBarrelIdx++ ][Barrel])/ pa;
285  }
286  //transition
287  pa = coerceToOne(m_eventStat.m_particleClustersBTE[0][Transition]);
288  std::array<double, 10> pcTransition2ff{};
289  size_t clusterTransitionIdx = 2;
290  for (auto & thisClusterT: pcTransition2ff){
291  thisClusterT = double(m_eventStat.m_particleClustersBTE[ clusterTransitionIdx++ ][Transition])/ pa;
292  }
293  //endcap
294  pa = coerceToOne(m_eventStat.m_particleClustersBTE[0][Endcap]);
295  std::array<double, 10> pcEndcap2ff{};
296  size_t clusterEndcapIdx = 2;
297  for (auto & thisClusterE: pcEndcap2ff){
298  thisClusterE = double(m_eventStat.m_particleClustersBTE[ clusterEndcapIdx++ ][Endcap])/ pa;
299  }
300  //fwd
301  pa = coerceToOne(m_eventStat.m_particleClustersBTE[0][3]);
302  std::array<double, 10> pcFwd2ff{};
303  size_t clusterFwdIdx = 2;
304  for (auto & thisClusterD: pcFwd2ff){
305  thisClusterD = double(m_eventStat.m_particleClustersBTE[ clusterFwdIdx++ ][Forward])/ pa;
306  }
307  //
308  //*** SPACE POINTS ***
309  //
310  //all
311  pa = coerceToOne(m_eventStat.m_particleSpacePoints[0]);
312  std::array<double, 10> sp2ff{};
313  size_t spacepointIdx = 2;
314  for (auto & thisSpacepoint: sp2ff){
315  thisSpacepoint = double(m_eventStat.m_particleSpacePoints[ spacepointIdx++ ])/ pa;
316  }
317  //barrel
318  pa = coerceToOne(m_eventStat.m_particleSpacePointsBTE[0][Barrel]);
319  std::array<double, 10> spBarrel2ff{};
320  size_t spacepointBarrelIdx = 2;
321  for (auto & thisSpacepoint: spBarrel2ff){
322  thisSpacepoint = double(m_eventStat.m_particleSpacePointsBTE[ spacepointBarrelIdx++ ][Barrel])/ pa;
323  }
324  //transition
325  pa = coerceToOne(m_eventStat.m_particleSpacePointsBTE[0][Transition]);
326  std::array<double, 10> spTransition2ff{};
327  size_t spacepointTransitionIdx = 2;
328  for (auto & thisSpacepoint: spTransition2ff){
329  thisSpacepoint = double(m_eventStat.m_particleSpacePointsBTE[ spacepointTransitionIdx++ ][Transition])/ pa;
330  }
331  //endcap
332  pa = coerceToOne(m_eventStat.m_particleSpacePointsBTE[0][Endcap]);
333  std::array<double, 10> spEndcap2ff{};
334  size_t spacepointEndcapIdx = 2;
335  for (auto & thisSpacepoint: spEndcap2ff){
336  thisSpacepoint = double(m_eventStat.m_particleSpacePointsBTE[ spacepointEndcapIdx++ ][Endcap])/ pa;
337  }
338  //Fwd
339  pa = coerceToOne(m_eventStat.m_particleSpacePointsBTE[0][Forward]);
340  std::array<double, 10> spFwd2ff{};
341  size_t spacepointFwdIdx = 2;
342  for (auto & thisSpacepoint: spFwd2ff){
343  thisSpacepoint = double(m_eventStat.m_particleSpacePointsBTE[ spacepointFwdIdx++ ][Forward])/ pa;
344  }
345  auto w8=std::setw(8);
346  out<<"| Probability for such charge particles to have some number silicon |\n";
347  out<<"| clusters | space points |\n";
348  out<<"| Total Barrel Transi Endcap Forward | Total Barrel Transi Endcap Forward |\n";
349 
350  for (size_t idx{0};idx != 10;++idx){
351  out<<"| >= "<<idx+2<< std::string((idx<8)?" ":" ")
352  <<w8<<p5<<pc2ff[idx]
353  <<w8<<p5<<pcBarrel2ff[idx]
354  <<w8<<p5<<pcTransition2ff[idx]
355  <<w8<<p5<<pcEndcap2ff[idx]
356  <<w8<<p5<<pcFwd2ff[idx]<<" | "
357 
358  <<w8<<p5<<sp2ff[idx]
359  <<w8<<p5<<spBarrel2ff[idx]
360  <<w8<<p5<<spTransition2ff[idx]
361  <<w8<<p5<<spEndcap2ff[idx]
362  <<w8<<p5<<spFwd2ff[idx]
363  <<" |\n";
364  }
365 
366  out<<topNtail(spaceSeparator)<<"\n";
367  out<<"| Additional cuts for truth particles are |\n";
368  out<<"| number silicon clusters >="<<w13<<m_clcut<<" |\n";
369  out<<"| number trt clusters >="<<w13<<m_clcutTRT<<" |\n";
370  out<<"| number space points >="<<w13<<m_spcut<<" |\n";
371 
372  pa = coerceToOne(m_eventStat.m_particleClusters[0]);
373  out<<"| Probability find truth particles with this cuts is "<<w8<<p5<<double(m_eventStat.m_events)/pa<<" |\n";
374  pa = coerceToOne(m_eventStat.m_particleClustersBTE[0][Barrel]);
375  out<<"| For barrel region "<<w8<<p5<<double(m_eventStat.m_eventsBTE[Barrel])/pa<<" |\n";
376  pa = coerceToOne(m_eventStat.m_particleClustersBTE[0][Transition]);
377  out<<"| For transition region "<<w8<<p5<<double(m_eventStat.m_eventsBTE[Transition])/pa<<" |\n";
378  pa = coerceToOne(m_eventStat.m_particleClustersBTE[0][Endcap]);
379  out<<"| For endcap region "<<w8<<p5<<double(m_eventStat.m_eventsBTE[Endcap])/pa<<" |\n";
380  pa = coerceToOne(m_eventStat.m_particleClustersBTE[0][Forward]);
381  out<<"| For forward region "<<w8<<p5<<double(m_eventStat.m_eventsBTE[Forward])/pa<<" |\n";
382  out<<"| |\n";
383  pa = coerceToOne(m_eventStat.m_nclustersNegBP);
384  double ratio = double(m_eventStat.m_nclustersPosBP)/pa;
385  double eratio = std::sqrt(ratio*(1.+ratio)/pa);
386  out<<"| Ratio barrel pixels clusters for +/- particles ="<<w8<<p5<<ratio<<" +-"<<w8<<p5<<eratio<<" |\n";
387  pa = coerceToOne(m_eventStat.m_nclustersNegEP);
388  ratio = double(m_eventStat.m_nclustersPosEP)/pa;
389  eratio = std::sqrt(ratio*(1.+ratio)/pa);
390  out<<"| Ratio endcap pixels clusters for +/- particles ="<<w8<<p5<<ratio<<" +-"<<w8<<p5<<eratio<<" |\n";
391  pa = coerceToOne(m_eventStat.m_nclustersNegBS);
392  ratio = double(m_eventStat.m_nclustersPosBS)/pa;
393  eratio = std::sqrt(ratio*(1.+ratio)/pa);
394  out<<"| Ratio barrel SCT clusters for +/- particles ="<<w8<<p5<<ratio<<" +-"<<w8<<p5<<eratio<<" |\n";
395  pa = coerceToOne(m_eventStat.m_nclustersNegES);
396  ratio = double(m_eventStat.m_nclustersPosES)/pa;
397  eratio = std::sqrt(ratio*(1.+ratio)/pa);
398  out<<"| Ratio endcap SCT clusters for +/- particles ="<<w8<<p5<<ratio<<" +-"<<w8<<p5<<eratio<<" |\n";
399  pa = double(m_eventStat.m_eventsNEG); if(pa < 1.) pa = 1.;
400  ratio = double(m_eventStat.m_eventsPOS)/pa;
401  eratio = std::sqrt(ratio*(1.+ratio)/pa);
402  out<<"| Number truth particles and +/- ratio ="<<std::setw(10)<<m_eventStat.m_events<<w8<<p5<<ratio<<" +-"<<w8<<p5<<eratio<<" |\n";
403  ratio = 0.;
404  if(m_eventStat.m_nclustersPTOT!=0) ratio = double(m_eventStat.m_nclustersPTOTt)/double(m_eventStat.m_nclustersPTOT);
405 
406  out<<"| Number pix clusters, truth clusters and ratio = "
407  <<std::setw(10)<<m_eventStat.m_nclustersPTOT
408  <<std::setw(10)<<m_eventStat.m_nclustersPTOTt
409  <<std::setw(12)<<std::setprecision(5)<<ratio<<" |\n";
410  ratio = 0.;
411  if(m_eventStat.m_nclustersSTOT!=0) ratio = double(m_eventStat.m_nclustersSTOTt)/double(m_eventStat.m_nclustersSTOT);
412  out<<"| Number sct clusters, truth clusters and ratio = "
413  <<std::setw(10)<<m_eventStat.m_nclustersSTOT
414  <<std::setw(10)<<m_eventStat.m_nclustersSTOTt
415  <<std::setw(12)<<std::setprecision(5)<<ratio<<" |\n";
416 
417  out<<"|-----------------------------------------------------------------------------------|\n\n";
418 
420  int nc = 0;
421  for(; t!=te; ++t) {
422  int n = 47-(t->key().size());
423  std::string s1; for(int i=0; i<n; ++i) s1.append(" "); s1.append("|");
424 
425 
426  out<<"|-----------------------------------------------------------------------------------|\n";
427  out<<"| Statistic for "<<(t->key())<<s1<<"\n";
428 
429  double ne = double(m_eventStat.m_events); if(ne < 1.) ne = 1.;
430  double ef [6]; for(int i=0; i!=6; ++i) ef [i] = double(m_trackCollectionStat[nc].m_efficiency [i]) /ne;
431  double ef0[6]; for(int i=0; i!=6; ++i) ef0[i] = double(m_trackCollectionStat[nc].m_efficiencyN [i][0])/ne;
432  double ef1[6]; for(int i=0; i!=6; ++i) ef1[i] = double(m_trackCollectionStat[nc].m_efficiencyN [i][1])/ne;
433  double ef2[6]; for(int i=0; i!=6; ++i) ef2[i] = double(m_trackCollectionStat[nc].m_efficiencyN [i][2])/ne;
434 
435 
436  using EffArray_t = std::array<double, 6>;
437  //
438  auto makeEffArray = [](const auto & threeDimArray, const size_t secondIdx, const size_t thirdIdx, const double denom){
439  EffArray_t result{};
440  size_t idx{0};
441  auto invDenom = 1./denom;
442  for (auto & entry: result){
443  entry = threeDimArray[idx++][secondIdx][thirdIdx]*invDenom;
444  }
445  return result;
446  };
447  //
448  const auto & efficiencyArrayInput = m_trackCollectionStat[nc].m_efficiencyBTE;
449  //
450  double neBTE = coerceToOne(m_eventStat.m_eventsBTE[Barrel]);
451  const EffArray_t efB0 = makeEffArray(efficiencyArrayInput,0,Barrel,neBTE);
452  const EffArray_t efB1 = makeEffArray(efficiencyArrayInput,1,Barrel,neBTE);
453  const EffArray_t efB2 = makeEffArray(efficiencyArrayInput,2,Barrel,neBTE);
454  const EffArray_t efB3 = makeEffArray(efficiencyArrayInput,3,Barrel,neBTE);
455  const EffArray_t efB4 = makeEffArray(efficiencyArrayInput,4,Barrel,neBTE);
456  //
457  neBTE = coerceToOne(m_eventStat.m_eventsBTE[Transition]);
458  const EffArray_t efT0 = makeEffArray(efficiencyArrayInput,0,Transition,neBTE);
459  const EffArray_t efT1 = makeEffArray(efficiencyArrayInput,1,Transition,neBTE);
460  const EffArray_t efT2 = makeEffArray(efficiencyArrayInput,2,Transition,neBTE);
461  const EffArray_t efT3 = makeEffArray(efficiencyArrayInput,3,Transition,neBTE);
462  const EffArray_t efT4 = makeEffArray(efficiencyArrayInput,4,Transition,neBTE);
463  //
464  neBTE = coerceToOne(m_eventStat.m_eventsBTE[Endcap]);
465  const EffArray_t efE0 = makeEffArray(efficiencyArrayInput,0,Endcap,neBTE);
466  const EffArray_t efE1 = makeEffArray(efficiencyArrayInput,1,Endcap,neBTE);
467  const EffArray_t efE2 = makeEffArray(efficiencyArrayInput,2,Endcap,neBTE);
468  const EffArray_t efE3 = makeEffArray(efficiencyArrayInput,3,Endcap,neBTE);
469  const EffArray_t efE4 = makeEffArray(efficiencyArrayInput,4,Endcap,neBTE);
470  //
471  neBTE = coerceToOne(m_eventStat.m_eventsBTE[Forward]);
472  const EffArray_t efD0 = makeEffArray(efficiencyArrayInput,0,Forward,neBTE);
473  const EffArray_t efD1 = makeEffArray(efficiencyArrayInput,1,Forward,neBTE);
474  const EffArray_t efD2 = makeEffArray(efficiencyArrayInput,2,Forward,neBTE);
475  const EffArray_t efD3 = makeEffArray(efficiencyArrayInput,3,Forward,neBTE);
476  const EffArray_t efD4 = makeEffArray(efficiencyArrayInput,4,Forward,neBTE);
477 
478 
479  double efrec = ef0[0]+ef0[1]+ef0[2]+ef1[0]+ef1[1]+ef2[0];
480  double efrecB = efB0[0]+efB0[1]+efB0[2]+efB1[0]+efB1[1]+efB2[0];
481  double efrecT = efT0[0]+efT0[1]+efT0[2]+efT1[0]+efT1[1]+efT2[0];
482  double efrecE = efE0[0]+efE0[1]+efE0[2]+efE1[0]+efE1[1]+efE2[0];
483  double efrecD = efD0[0]+efD0[1]+efD0[2]+efD1[0]+efD1[1]+efD2[0];
484 
485  ne = coerceToOne(m_eventStat.m_eventsPOS);
486  double efP[6]; for(int i=0; i!=6; ++i) efP[i] = double(m_trackCollectionStat[nc].m_efficiencyPOS[i])/ne;
487  ne = coerceToOne(m_eventStat.m_eventsNEG);
488  double efN[6]; for(int i=0; i!=6; ++i) efN[i] = double(m_trackCollectionStat[nc].m_efficiencyNEG[i])/ne;
489 
490  out<<"|-----------------------------------------------------------------------------------|\n";
491  out<<"| Probability to lose 0 1 2 3 4 >=5 clusters |\n";
492  out<<"|-----------------------------------------------------------------------------------|\n";
493 
494  auto formattedOutput=[&out](auto & effArray){
495  for (size_t i{};i!=6;++i){
496  out<<std::setw(9)<<std::setprecision(4)<<effArray[i];
497  }
498  out<<" |\n";
499  };
500 
501  out<<"| For all particles ";
502  formattedOutput(ef);
503  out<<"| For + particles ";
504  formattedOutput(efP);
505  out<<"| For - particles ";
506  formattedOutput(efN);
507  out<<"|-----------------------------------------------------------------------------------|\n";
508  out<<"| Barrel region |\n";
509  out<<"| 0 wrong clusters ";
510  formattedOutput(efB0);
511  out<<"| 1 wrong clusters ";
512  formattedOutput(efB1);
513  out<<"| 2 wrong clusters ";
514  formattedOutput(efB2);
515  out<<"| 3 wrong clusters ";
516  formattedOutput(efB3);
517  out<<"| >=4 wrong clusters ";
518  formattedOutput(efB4);
519  out<<"|-----------------------------------------------------------------------------------|\n";
520  out<<"| Transition region |\n";
521  out<<"| 0 wrong clusters ";
522  formattedOutput(efT0);
523  out<<"| 1 wrong clusters ";
524  formattedOutput(efT1);
525  out<<"| 2 wrong clusters ";
526  formattedOutput(efT2);
527  out<<"| 3 wrong clusters ";
528  formattedOutput(efT3);
529  out<<"| >=4 wrong clusters ";
530  formattedOutput(efT4);
531  out<<"|-----------------------------------------------------------------------------------|\n";
532  out<<"| Endcap region |\n";
533  out<<"| 0 wrong clusters ";
534  formattedOutput(efE0);
535  out<<"| 1 wrong clusters ";
536  formattedOutput(efE1);
537  out<<"| 2 wrong clusters ";
538  formattedOutput(efE2);
539  out<<"| 3 wrong clusters ";
540  formattedOutput(efE3);
541  out<<"| >=4 wrong clusters ";
542  formattedOutput(efE4);
543  out<<"|-----------------------------------------------------------------------------------|\n";
544  out<<"| Forward region |\n";
545  out<<"| 0 wrong clusters ";
546  formattedOutput(efD0);
547  out<<"| 1 wrong clusters ";
548  formattedOutput(efD1);
549  out<<"| 2 wrong clusters ";
550  formattedOutput(efD2);
551  out<<"| 3 wrong clusters ";
552  formattedOutput(efD3);
553  out<<"| >=4 wrong clusters ";
554  formattedOutput(efD4);
555 
556  out<<"|-----------------------------------------------------------------------------------|\n";
557  pa = coerceToOne(m_eventStat.m_particleClusters[0]);
558  out<<"| Efficiency reconstruction (number lose+wrong < 3) = "
559  <<std::setw(9)<<std::setprecision(5)<<efrec
560  <<" ("
561  <<std::setw(9)<<std::setprecision(5)<<efrec*double(m_eventStat.m_events)/pa
562  <<" ) "
563  <<" |\n";
564  pa = coerceToOne(m_eventStat.m_particleClustersBTE[0][Barrel]);
565  out<<"| For barrel region = "
566  <<std::setw(9)<<std::setprecision(5)<<efrecB
567  <<" ("
568  <<std::setw(9)<<std::setprecision(5)<<efrecB*double(m_eventStat.m_eventsBTE[Barrel])/pa
569  <<" ) "
570  <<" |\n";
571  pa = coerceToOne(m_eventStat.m_particleClustersBTE[0][Transition]);
572  out<<"| For transition region = "
573  <<std::setw(9)<<std::setprecision(5)<<efrecT
574  <<" ("
575  <<std::setw(9)<<std::setprecision(5)<<efrecT*double(m_eventStat.m_eventsBTE[Transition])/pa
576  <<" ) "
577  <<" |\n";
578  pa = coerceToOne(m_eventStat.m_particleClustersBTE[0][Endcap]);
579  out<<"| For endcap region = "
580  <<std::setw(9)<<std::setprecision(5)<<efrecE
581  <<" ("
582  <<std::setw(9)<<std::setprecision(5)<<efrecE*double(m_eventStat.m_eventsBTE[Endcap])/pa
583  <<" ) "
584  <<" |\n";
585  pa = coerceToOne(m_eventStat.m_particleClustersBTE[0][Forward]);
586  out<<"| For forward region = "
587  <<std::setw(9)<<std::setprecision(5)<<efrecD
588  <<" ("
589  <<std::setw(9)<<std::setprecision(5)<<efrecD*double(m_eventStat.m_eventsBTE[Forward])/pa
590  <<" ) "
591  <<" |\n";
592 
593  out<<"|-----------------------------------------------------------------------------------|\n";
594  out<<"| Reconstructed tracks + - +/-ratio error |\n";
595  out<<"|-----------------------------------------------------------------------------------|\n";
596 
597  pa = coerceToOne(m_trackCollectionStat[nc].m_ntracksNEGB);
598  ratio = double(m_trackCollectionStat[nc].m_ntracksPOSB)/pa;
599  eratio = std::sqrt(ratio*(1.+ratio)/pa);
600 
601  out<<"| Barrel "
602  <<std::setw(10)<<m_trackCollectionStat[nc].m_ntracksPOSB
603  <<std::setw(11)<<m_trackCollectionStat[nc].m_ntracksNEGB
604  <<std::setw(11)<<std::setprecision(5)<<ratio
605  <<std::setw(11)<<std::setprecision(5)<<eratio<<" |\n";
606  pa = coerceToOne(m_trackCollectionStat[nc].m_ntracksNEGE);
607  ratio = double(m_trackCollectionStat[nc].m_ntracksPOSE)/pa;
608  eratio = std::sqrt(ratio*(1.+ratio)/pa);
609 
610  out<<"| Endcap "
611  <<std::setw(10)<<m_trackCollectionStat[nc].m_ntracksPOSE
612  <<std::setw(11)<<m_trackCollectionStat[nc].m_ntracksNEGE
613  <<std::setw(11)<<std::setprecision(5)<<ratio
614  <<std::setw(11)<<std::setprecision(5)<<eratio<<" |\n";
615  pa = coerceToOne(m_trackCollectionStat[nc].m_ntracksNEGFWD);
616  ratio = double(m_trackCollectionStat[nc].m_ntracksPOSFWD)/pa;
617  eratio = std::sqrt(ratio*(1.+ratio)/pa);
618 
619  out<<"| Forward "
620  <<std::setw(10)<<m_trackCollectionStat[nc].m_ntracksPOSFWD
621  <<std::setw(11)<<m_trackCollectionStat[nc].m_ntracksNEGFWD
622  <<std::setw(11)<<std::setprecision(5)<<ratio
623  <<std::setw(11)<<std::setprecision(5)<<eratio<<" |\n";
624 
625 
626 
627  int nt=0;
628  int ft=0;
629  int kf=0;
630  for(int k = 0; k!=50; ++k) {
631  nt+=m_trackCollectionStat[nc].m_total[k];
632  ft+=m_trackCollectionStat[nc].m_fake [k];
633  if(!kf && nt) kf = k;
634  }
635 
636  if(kf) {
637 
638  out<<"|-----------------------------------------------------------------------------------|\n";
639  out<<"| Fake tracks rate for different number of clusters on track |\n";
640  out<<"|-----------------------------------------------------------------------------------|\n";
641 
642  for(int k = kf; k!=kf+6; ++k) {
643  out<<"| >= "<<std::setw(2)<<k<<" ";
644  }
645  out<<"|\n";
646 
647  for(int k = kf; k!=kf+6; ++k) {
648  double eff = 0.; if(nt>0) eff = double(ft)/double(nt);
649  out<<"|"<<std::setw(12)<<std::setprecision(5)<<eff<<" ";
650  nt-=m_trackCollectionStat[nc].m_total[k];
651  ft-=m_trackCollectionStat[nc].m_fake [k];
652  }
653  out<<"|\n";
654  out<<"|-----------------------------------------------------------------------------------|\n";
655  }
656  ++nc;
657  }
658  ATH_MSG_INFO("\n"<<out.str());
659  return StatusCode::SUCCESS;
660 }
661 
663 // Dumps conditions information into the MsgStream
665 
667 {
669 
670  int n;
671  out << level <<"\n";
672  out<<"|--------------------------------------------------------------------------------------------------------------------|\n";
673  for(; t!=te; ++t) {
674  n = 65-t->key().size();
675  std::string s1; for(int i=0; i<n; ++i) s1.append(" "); s1.append("|");
676  out<<"| Location of input tracks | "<<t->key()<<s1<<"\n";
677  }
678  auto padString = [](const std::string & s){
679  const int n = 65 - s.size();
680  return s + std::string(n, ' ');
681  };
682  std::string s2 = padString(m_spacepointsPixelname.key());
683  std::string s3 = padString(m_spacepointsSCTname.key());
684  std::string s4 = padString(m_spacepointsOverlapname.key());
685  std::string s5 = padString(m_clustersPixelname.key());
686  std::string s6 = padString(m_clustersSCTname.key());
687  //
688  std::string s9 = padString(m_clustersTRTname.key());
689  //
690  std::string s7 = padString(m_truth_locationPixel.key());
691  std::string s8 = padString(m_truth_locationSCT.key());
692  //
693  std::string s10 = padString(m_truth_locationTRT.key());
694 
695  out<<"| Pixel space points | "<<s2<<"|\n";
696  out<<"| SCT space points | "<<s3<<"|\n";
697  out<<"| Overlap space points | "<<s4<<"|\n";
698  out<<"| Pixel clusters | "<<s5<<"|\n";
699  out<<"| SCT clusters | "<<s6<<"|\n";
700  out<<"| TRT clusters | "<<s9<<"|\n";
701  out<<"| Truth location for pixels | "<<s7<<"|\n";
702  out<<"| Truth location for sct | "<<s8<<"|\n";
703  out<<"| Truth location for trt | "<<s10<<"|\n";
704  out<<"| pT cut | "
705  <<std::setw(14)<<std::setprecision(5)<<m_ptcut
706  <<" |\n";
707  out<<"| max pT cut | "
708  <<std::setw(14)<<std::setprecision(5)<<m_ptcutmax
709  <<" |\n";
710  out<<"| rapidity cut | "
711  <<std::setw(14)<<std::setprecision(5)<<m_rapcut
712  <<" |\n";
713  out<<"| min Radius | "
714  <<std::setw(14)<<std::setprecision(5)<<m_rmin
715  <<" |\n";
716  out<<"| max Radius | "
717  <<std::setw(14)<<std::setprecision(5)<<m_rmax
718  <<" |\n";
719  out<<"| Min. number clusters for generated track | "
720  <<std::setw(14)<<std::setprecision(5)<<m_clcut
721  <<" |\n";
722  out<<"| Min. number sp.points for generated track | "
723  <<std::setw(14)<<std::setprecision(5)<<m_spcut
724  <<" |\n";
725  out<<"| Min. number TRT clusters for generated track | "
726  <<std::setw(14)<<std::setprecision(5)<<m_clcutTRT
727  <<" |\n";
728  out<<"|--------------------------------------------------------------------------------------------------------------------|\n";
729  return out;
730 }
731 
733 // Dumps event information into the ostream
735 
737 {
738  out << MSG::DEBUG << "\n";
739  auto formatOutput = [&out](const auto val){
740  out<<std::setw(12)<<val
741  <<" |\n";
742  };
743  out<<"|---------------------------------------------------------------------|\n";
744  out<<"| m_nspacepoints | ";
745  formatOutput(event_data.m_nspacepoints);
746  out<<"| m_nclusters | ";
747  formatOutput(event_data.m_nclusters);
748  out<<"| Kine-Clusters size | ";
749  formatOutput(event_data.m_kinecluster.size());
750  out<<"| Kine-TRTClusters size | ";
751  formatOutput(event_data.m_kineclusterTRT.size());
752  out<<"| Kine-SpacePoints size | ";
753  formatOutput(event_data.m_kinespacepoint.size());
754  out<<"| Number good kine tracks | ";
755  formatOutput(event_data.m_nqtracks);
756  out<<"|---------------------------------------------------------------------|\n";
757  return out;
758 }
759 
760 
762 // New event for clusters information
764 
766 {
767  std::lock_guard<std::mutex> lock(m_statMutex);
768 
769  // Get pixel clusters container
770  //
771  std::unique_ptr<SG::ReadHandle<SiClusterContainer> > pixelcontainer;
772  std::unique_ptr<SG::ReadHandle<SiClusterContainer> > sctcontainer;
773  std::unique_ptr<SG::ReadHandle<TRT_DriftCircleContainer> > trtcontainer;
774 
775  if(m_usePIX) {
776  pixelcontainer = std::make_unique<SG::ReadHandle<SiClusterContainer> >(m_clustersPixelname,ctx);
777  if (!pixelcontainer->isValid()) ATH_MSG_DEBUG("Failed to create Pixel clusters container read handle with key " << m_clustersPixelname.key());
778  }
779 
780  // Get sct clusters container
781  //
782  if(m_useSCT) {
783  sctcontainer = std::make_unique<SG::ReadHandle<SiClusterContainer> >(m_clustersSCTname,ctx);
784  if (!sctcontainer->isValid()) ATH_MSG_DEBUG("Failed to create SCT clusters container read handle with key " << m_clustersSCTname.key());
785  }
786 
787  // Get trt cluster container
788  //
789  if(m_clcutTRT > 0) {
790  trtcontainer = std::make_unique<SG::ReadHandle<TRT_DriftCircleContainer> >(m_clustersTRTname,ctx);
791  if (!trtcontainer->isValid()) ATH_MSG_DEBUG("Failed to create TRT drift circle container read handle with key " << m_clustersTRTname.key());
792  }
793 
794  int Kine[1000];
795 
796  event_data.m_clusterHandles.reserve(3);
797  // Loop through all pixel clusters
798  //
799  if(pixelcontainer && pixelcontainer->isValid()) {
800  InDet::SiClusterContainer::const_iterator w = (*pixelcontainer)->begin();
801  InDet::SiClusterContainer::const_iterator we = (*pixelcontainer)->end ();
802 
803  for(; w!=we; ++w) {
804 
805  InDet::SiClusterCollection::const_iterator c = (*w)->begin();
806  InDet::SiClusterCollection::const_iterator ce = (*w)->end ();
807 
808  for(; c!=ce; ++c) {
809 
810  ++event_data.m_nclusters;
811  ++m_eventStat.m_nclustersPTOT;
812  if(isTruth(event_data,(*c))) ++m_eventStat.m_nclustersPTOTt;
813 
814 
815  int nk = kine(event_data,(*c),Kine,999);
816  for(int i=0; i!=nk; ++i) {
817  if(!isTheSameDetElement(event_data,Kine[i],(*c))) {
818  event_data.m_kinecluster.insert(std::make_pair(Kine[i],(*c)));
819  }
820  }
821  }
822  }
823  event_data.m_clusterHandles.push_back(std::move(pixelcontainer));
824 
825  }
826 
827  // Loop through all sct clusters
828  //
829  if(sctcontainer && sctcontainer->isValid()) {
830  InDet::SiClusterContainer::const_iterator w = (*sctcontainer)->begin();
831  InDet::SiClusterContainer::const_iterator we = (*sctcontainer)->end ();
832 
833  for(; w!=we; ++w) {
834 
835  InDet::SiClusterCollection::const_iterator c = (*w)->begin();
836  InDet::SiClusterCollection::const_iterator ce = (*w)->end ();
837 
838  for(; c!=ce; ++c) {
839 
840  ++event_data.m_nclusters;
841  ++m_eventStat.m_nclustersSTOT;
842  if(isTruth(event_data,(*c))) ++m_eventStat.m_nclustersSTOTt;
843 
844  int nk = kine(event_data,(*c),Kine,999);
845  for(int i=0; i!=nk; ++i) {
846  if(!isTheSameDetElement(event_data,Kine[i],(*c))) event_data.m_kinecluster.insert(std::make_pair(Kine[i],(*c)));
847  }
848  }
849  }
850  event_data.m_clusterHandles.push_back(std::move(sctcontainer));
851  }
852 
853  if(trtcontainer && trtcontainer->isValid()) {
854  // Loop through all trt clusters
855  //
856  InDet::TRT_DriftCircleContainer::const_iterator w = (*trtcontainer)->begin();
857  InDet::TRT_DriftCircleContainer::const_iterator we = (*trtcontainer)->end ();
858 
859  for(; w!=we; ++w) {
860 
861  InDet::TRT_DriftCircleCollection::const_iterator c = (*w)->begin();
862  InDet::TRT_DriftCircleCollection::const_iterator ce = (*w)->end ();
863 
864  for(; c!=ce; ++c) {
865 
866  ++event_data.m_nclustersTRT;
867  int nk = kine(event_data,(*c),Kine,999);
868  for(int i=0; i!=nk; ++i) event_data.m_kineclusterTRT.insert(std::make_pair(Kine[i],(*c)));
869  }
870  }
871  event_data.m_clusterHandles.push_back(std::move(trtcontainer));
872  }
873 }
874 
876 // New event for space points information
878 
880 {
881 
882  int Kine[1000];
883 
884  if(m_usePIX && !m_spacepointsPixelname.key().empty()) {
885  event_data.m_spacePointContainer.emplace_back(m_spacepointsPixelname,ctx);
886  if (!event_data.m_spacePointContainer.back().isValid()) {
887  ATH_MSG_DEBUG( "Invalid Pixels space points container read handle for key " << m_spacepointsPixelname.key() );
888  }
889  else {
890  SpacePointContainer::const_iterator spc = event_data.m_spacePointContainer.back()->begin();
891  SpacePointContainer::const_iterator spce = event_data.m_spacePointContainer.back()->end ();
892  for(; spc != spce; ++spc) {
893  SpacePointCollection::const_iterator sp = (*spc)->begin();
894  SpacePointCollection::const_iterator spe = (*spc)->end ();
895 
896  for(; sp != spe; ++sp) {
897 
898  ++event_data.m_nspacepoints;
899  int nk = kine(event_data,(*sp)->clusterList().first,Kine,999);
900  for(int i=0; i!=nk; ++i) {
901 
902  if(!isTheSameDetElement(event_data,Kine[i],(*sp))) {
903  event_data.m_kinespacepoint.insert(std::make_pair(Kine[i],(*sp)));
904  }
905  }
906  }
907  }
908  }
909  }
910  // Get sct space points containers from store gate
911  //
912  if(m_useSCT && !m_spacepointsSCTname.key().empty()) {
913  event_data.m_spacePointContainer.emplace_back(m_spacepointsSCTname,ctx);
914  if (!event_data.m_spacePointContainer.back().isValid()) {
915  ATH_MSG_DEBUG( "Invalid SCT space points container read handle for key " << m_spacepointsSCTname.key() );
916  }
917  else {
918  SpacePointContainer::const_iterator spc = event_data.m_spacePointContainer.back()->begin();
919  SpacePointContainer::const_iterator spce = event_data.m_spacePointContainer.back()->end ();
920 
921  for(; spc != spce; ++spc) {
922 
923  SpacePointCollection::const_iterator sp = (*spc)->begin();
924  SpacePointCollection::const_iterator spe = (*spc)->end ();
925 
926  for(; sp != spe; ++sp) {
927 
928 
929  ++event_data.m_nspacepoints;
930  int nk = kine(event_data,(*sp)->clusterList().first,(*sp)->clusterList().second,Kine,999);
931  for(int i=0; i!=nk; ++i) {
932  if(!isTheSameDetElement(event_data,Kine[i],(*sp))) {
933  event_data.m_kinespacepoint.insert(std::make_pair(Kine[i],(*sp)));
934  }
935  }
936  }
937  }
938  }
939  }
940  // Get sct overlap space points containers from store gate
941  //
942  if(m_useSCT && !m_spacepointsOverlapname.key().empty()) {
943  event_data.m_spacepointsOverlap=std::make_unique< SG::ReadHandle<SpacePointOverlapCollection> >(m_spacepointsOverlapname,ctx);
944  if (!event_data.m_spacepointsOverlap->isValid()) {
945  ATH_MSG_DEBUG( "Invalid overlap space points container read handle for key " << m_spacepointsOverlapname.key() );
946  }
947  else {
950 
951  for (; sp!=spe; ++sp) {
952 
953  ++event_data.m_nspacepoints;
954  int nk = kine(event_data,(*sp)->clusterList().first,(*sp)->clusterList().second,Kine,999);
955  for(int i=0; i!=nk; ++i) {
956  if(!isTheSameDetElement(event_data,Kine[i],(*sp))) {
957  event_data.m_kinespacepoint.insert(std::make_pair(Kine[i],(*sp)));
958  }
959  }
960  }
961  }
962  }
963 }
965 // Good kine tracks selection
967 
969 {
973 
974  if( c == event_data.m_kinecluster.end()) {
975  return 0;
976  }
977 
978  if( s == event_data.m_kinespacepoint.end()) {
979  return 0;
980  }
981 
982  if( m_clcutTRT > 0 && u == event_data.m_kineclusterTRT.end()) {
983  return 0;
984  }
985 
986  std::list<int> worskine;
987 
988  int rp = 0;
989  int t = 0;
990  int k0 = (*c).first;
991  int q0 = k0*charge(event_data,(*c),rp);
992  unsigned int nc = 1 ;
993 
994  auto coerceTo49 = [] (const size_t idx){
995  return (idx<50) ? idx : 49;
996  };
997 
998  for(++c; c!=event_data.m_kinecluster.end(); ++c) {
999 
1000  if((*c).first==k0) {++nc; continue;}
1001  q0 = charge(event_data,(*c),rp)*k0;
1002  //
1003  const size_t clusterIdx =coerceTo49(nc);
1004  ++event_data.m_eventStat.m_particleClusters [clusterIdx];
1005  ++event_data.m_eventStat.m_particleClustersBTE[clusterIdx][rp];
1006  //
1007  int ns = event_data.m_kinespacepoint.count(k0);
1008  const size_t spacepointIdx =coerceTo49(ns);
1009  ++event_data.m_eventStat.m_particleSpacePoints [spacepointIdx];
1010  ++event_data.m_eventStat.m_particleSpacePointsBTE[spacepointIdx][rp];
1011 
1012  if (nc < m_clcut ) worskine.push_back(k0);
1013  else if(event_data.m_kinespacepoint.count(k0)< m_spcut ) worskine.push_back(k0);
1014  else if(event_data.m_kineclusterTRT.count(k0)< m_clcutTRT) worskine.push_back(k0);
1015  else {
1016  event_data.m_particles[0].push_back(InDet::PartPropCache(q0,rp)); ++t;
1017  }
1018 
1019  k0 = (*c).first;
1020  q0 =charge(event_data,(*c),rp)*k0;
1021  nc = 1 ;
1022  }
1023 
1024  ++event_data.m_eventStat.m_particleClusters [coerceTo49(nc)];
1025  ++event_data.m_eventStat.m_particleClustersBTE[coerceTo49(nc)][rp];
1026  int ns = event_data.m_kinespacepoint.count(k0);
1027  ++event_data.m_eventStat.m_particleSpacePoints [coerceTo49(ns)];
1028  ++event_data.m_eventStat.m_particleSpacePointsBTE[coerceTo49(ns)][rp];
1029 
1030  if (nc < m_clcut ) worskine.push_back(k0);
1031  else if(event_data.m_kinespacepoint.count(k0)< m_spcut ) worskine.push_back(k0);
1032  else if(event_data.m_kineclusterTRT.count(k0)< m_clcutTRT) worskine.push_back(k0);
1033  else {
1034  event_data.m_particles[0].push_back(InDet::PartPropCache(q0,rp)); ++t;
1035  }
1036  for(auto & pThisCluster: worskine) {
1037  event_data.m_kinecluster .erase(pThisCluster);
1038  event_data.m_kineclusterTRT.erase(pThisCluster);
1039  event_data.m_kinespacepoint.erase(pThisCluster);
1040  }
1041 
1042  for(c = event_data.m_kinecluster.begin(); c!= event_data.m_kinecluster.end(); ++c) {
1043  const Trk::PrepRawData*
1044  d = (*c).second;
1046  de= dynamic_cast<const InDetDD::SiDetectorElement*>(d->detectorElement());
1047  if (not de) continue;
1048  int q = charge(event_data,*c,rp);
1049 
1050  if (q<0) {
1051  if(de->isBarrel()) {
1052  de->isPixel() ? ++event_data.m_eventStat.m_nclustersNegBP : ++event_data.m_eventStat.m_nclustersNegBS;
1053  }
1054  else {
1055  de->isPixel() ? ++event_data.m_eventStat.m_nclustersNegEP : ++event_data.m_eventStat.m_nclustersNegES;
1056  }
1057 
1058  }
1059  else if(q>0) {
1060  if(de->isBarrel()) {
1061  de->isPixel() ? ++event_data.m_eventStat.m_nclustersPosBP : ++event_data.m_eventStat.m_nclustersPosBS;
1062  }
1063  else {
1064  de->isPixel() ? ++event_data.m_eventStat.m_nclustersPosEP : ++event_data.m_eventStat.m_nclustersPosES;
1065  }
1066  }
1067  }
1068 
1069 
1070 
1071  for(const auto& p: event_data.m_particles[0]) {
1072  for(SG::ReadHandleKeyArray<TrackCollection>::size_type nc=1; nc<m_tracklocation.size(); ++nc) event_data.m_particles[nc].push_back(p);
1073  }
1074  return t;
1075 }
1076 
1078 // Recontructed track comparison with kine information
1080 
1082 {
1083  if(!event_data.m_nqtracks) return;
1084 
1085 
1086  int nc = -1;
1087  event_data.m_trackcontainer.reserve(m_tracklocation.size());
1088  for(const SG::ReadHandleKey<TrackCollection> &track_key : m_tracklocation ) {
1089  if(++nc >= 100) return;
1090  event_data.m_tracks[nc].clear();
1091 
1092  event_data.m_trackcontainer.emplace_back(track_key,ctx );
1093  if (!event_data.m_trackcontainer.back().isValid()) {
1094  continue;
1095  }
1096 
1097  // Loop through all found tracks
1098  //
1099  TrackCollection::const_iterator t,te = event_data.m_trackcontainer.back()->end();
1100 
1101  int KINE[200],NKINE[200];
1102 
1103  for (t=event_data.m_trackcontainer.back()->begin(); t!=te; ++t) {
1104 
1105  Trk::TrackStates::const_iterator s = (*t)->trackStateOnSurfaces()->begin(),
1106  se = (*t)->trackStateOnSurfaces()->end();
1107 
1108  int NK = 0;
1109  int NC = 0;
1110  int N0 = 0;
1111  int nkm = 0;
1112  bool qp = false;
1113 
1114  const Trk::TrackParameters* tpf = (*s)->trackParameters(); if(!tpf) continue;
1115  const AmgVector(5)& Vpf = tpf ->parameters ();
1116  double pTf = std::abs(std::sin(Vpf[3])/Vpf[4]);
1117  bool qTf = pTf > m_ptcut;
1118  for(; s!=se; ++s) {
1119 
1120  if(!qp) {
1121 
1122  const Trk::TrackParameters* tp = (*s)->trackParameters();
1123 
1124  if(tp) {
1125  qp = true;
1126  const AmgVector(5)& Vp = tp->parameters();
1127  double pT = std::sin(Vp[3])/Vp[4] ;
1128  double rap = std::abs(std::log(std::tan(.5*Vp[3])));
1129  if (pT > m_ptcut && pT < m_ptcutmax) {
1130  if (rap < 1. ) ++event_data.m_trackCollectionStat[nc].m_ntracksPOSB;
1131  else if(rap < 3.0) ++event_data.m_trackCollectionStat[nc].m_ntracksPOSE;
1132  else if(rap < m_rapcut) ++event_data.m_trackCollectionStat[nc].m_ntracksPOSFWD;
1133  }
1134  else if(pT < -m_ptcut && pT > -m_ptcutmax) {
1135  if (rap < 1. ) ++event_data.m_trackCollectionStat[nc].m_ntracksNEGB;
1136  else if(rap < 3.0) ++event_data.m_trackCollectionStat[nc].m_ntracksNEGE;
1137  else if(rap < m_rapcut) ++event_data.m_trackCollectionStat[nc].m_ntracksNEGFWD;
1138  }
1139  }
1140  }
1141 
1142  if(!m_useOutliers && !(*s)->type(Trk::TrackStateOnSurface::Measurement)) continue;
1143 
1144  const Trk::MeasurementBase* mb = (*s)->measurementOnTrack();
1145  if(!mb) continue;
1146 
1147  const Trk::RIO_OnTrack* ri = dynamic_cast<const Trk::RIO_OnTrack*>(mb);
1148  if(!ri) continue;
1149 
1150  const Trk::PrepRawData* rd = ri->prepRawData();
1151  if(!rd) continue;
1152 
1153  const InDet::SiCluster* si = dynamic_cast<const InDet::SiCluster*>(rd);
1154  if(!si) continue;
1155 
1156  if(!m_usePIX && dynamic_cast<const InDet::PixelCluster*>(si)) continue;
1157  if(!m_useSCT && dynamic_cast<const InDet::SCT_Cluster*> (si)) continue;
1158 
1159 
1160  int Kine[1000], nk=kine0(event_data,rd,Kine,999); ++NC; if(!nk) ++N0;
1161 
1162  for(int k = 0; k!=nk; ++k) {
1163 
1164  int n = 0;
1165  for(; n!=NK; ++n) {if(Kine[k]==KINE[n]) {++NKINE[n]; break;}}
1166  if(n==NK) {KINE[NK] = Kine[k]; NKINE[NK] = 1; if (NK < 200) ++NK;}
1167  }
1168  for(int n=0; n!=NK; ++n) {if(NKINE[n]>nkm) nkm = NKINE[n];}
1169  }
1170 
1171  for(int n=0; n!=NK; ++n) {
1172  if(NKINE[n]==nkm) {
1173  int NQ = 1000*NKINE[n]+(NC-NKINE[n]);
1174 
1175  event_data.m_tracks[nc].insert(std::make_pair(KINE[n],NQ));
1176  if(qTf) {
1177  if(NC-N0 > 2) {
1178  ++event_data.m_trackCollectionStat[nc].m_total[NC]; if(NC-NKINE[n] > 2) {++event_data.m_trackCollectionStat[nc].m_fake[NC];}
1179  }
1180  }
1181  }
1182  }
1183  }
1184 
1185  }
1186 }
1187 
1189 // Particles and reconstructed tracks comparision
1191 
1193 {
1195 
1196  event_data.m_difference[nc].clear();
1197  auto p = event_data.m_particles[nc].begin(), pe =event_data.m_particles[nc].end();
1198  if(p==pe) return;
1199  std::multimap<int,int>::const_iterator t, te = event_data.m_tracks[nc].end();
1200 
1201  while (p!=pe) {
1202 
1203  int k = (*p).barcode();
1204  int n = event_data.m_kinecluster.count(k);
1205  int m = 0;
1206  int w = 0;
1207  t = event_data.m_tracks[nc].find(k);
1208  for(; t!=te; ++t) {
1209  if((*t).first!=k) break;
1210  int ts = (*t).second/1000;
1211  int ws = (*t).second%1000;
1212  if (ts > m ) {m = ts; w = ws;}
1213  else if(ts==m && w > ws) { w = ws;}
1214  }
1215  int d = n-m; if(d<0) d = 0; else if(d > 5) d=5; if(w>4) w = 4;
1216  if(m) {
1217  ++event_data.m_trackCollectionStat[nc].m_efficiency [d];
1218  ++event_data.m_trackCollectionStat[nc].m_efficiencyN[d][w];
1219  }
1220  int ch = (*p).charge();
1221  if(m) {
1222  ++event_data.m_trackCollectionStat[nc].m_efficiencyBTE[d][w][(*p).rapidity()];
1223  ch > 0 ? ++event_data.m_trackCollectionStat[nc].m_efficiencyPOS[d] : ++event_data.m_trackCollectionStat[nc].m_efficiencyNEG[d];
1224  }
1225  if(nc==0) {
1226  ++event_data.m_eventStat.m_events; ch > 0 ? ++event_data.m_eventStat.m_eventsPOS : ++event_data.m_eventStat.m_eventsNEG;
1227  ++event_data.m_eventStat.m_eventsBTE[(*p).rapidity()];
1228  }
1229  if(d==0) event_data.m_particles[nc].erase(p++);
1230  else {event_data.m_difference[nc].push_back(n-m); ++p;}
1231  }
1232  }
1233 }
1234 
1236 // Pointer to particle production for space point
1238 
1240 (const InDet::TrackClusterAssValidation::EventData_t &event_data,const Trk::PrepRawData* d1,const Trk::PrepRawData* d2,int* Kine,int nmax) const
1241 {
1242  int nkine = 0;
1243  int Kine1[1000],Kine2[1000];
1244  int n1 = kine(event_data,d1,Kine1,nmax); if(!n1) return nkine;
1245  int n2 = kine(event_data,d2,Kine2,nmax); if(!n2) return nkine;
1246 
1247  for(int i = 0; i!=n1; ++i) {
1248  for(int j = 0; j!=n2; ++j) {
1249  if(Kine1[i]==Kine2[j]) {Kine[nkine++] = Kine1[i]; break;}
1250  }
1251  }
1252  return nkine;
1253 }
1254 
1256 // Pointer to particle production for cluster
1258 
1260 (const InDet::TrackClusterAssValidation::EventData_t &event_data,const Trk::PrepRawData* d,int* Kine,int nmax) const
1261 {
1262 
1263  PRD_MultiTruthCollection::const_iterator mce;
1264  PRD_MultiTruthCollection::const_iterator mc = findTruth(event_data,d,mce);
1265 
1266  Identifier ID = d->identify();
1267  int nkine = 0;
1268 
1269  for(; mc!=mce; ++mc) {
1270 
1271  if( (*mc).first != ID ) return nkine;
1272 
1273  int k = (*mc).second.barcode(); if(k<=0) continue;
1274 
1275  const HepMC::ConstGenParticlePtr pa = (*mc).second.cptr();
1276  if(!pa or !pa->production_vertex()) continue;
1277 
1278  int pdg = std::abs(pa->pdg_id()); if(m_pdg && m_pdg != pdg ) continue;
1279 
1280  const HepPDT::ParticleData* pd = m_particleDataTable->particle(pdg);
1281  if(!pd or std::abs(pd->charge()) < .5) continue;
1282 
1283  // pT cut
1284  //
1285  double px = pa->momentum().px();
1286  double py = pa->momentum().py();
1287  double pz = pa->momentum().pz();
1288  double pt = std::sqrt(px*px+py*py);
1289  if( pt < m_ptcut || pt > m_ptcutmax) continue;
1290 
1291  // Rapidity cut
1292  //
1293  double t = std::abs(pz)/pt;
1294  if( t > m_tcut ) continue;
1295 
1296  // Radius cut
1297  //
1298  double vx = pa->production_vertex()->position().x();
1299  double vy = pa->production_vertex()->position().y();
1300  double r = std::sqrt(vx*vx+vy*vy);
1301  if( r < m_rmin || r > m_rmax) continue;
1302 
1303  Kine[nkine] = k; if(++nkine >= nmax) break;
1304  }
1305  return nkine;
1306 }
1307 
1309 // Pointer to particle production for cluster
1311 
1313 (const InDet::TrackClusterAssValidation::EventData_t &event_data,const Trk::PrepRawData* d,int* Kine,int nmax)
1314 {
1315 
1316  PRD_MultiTruthCollection::const_iterator mce;
1317  PRD_MultiTruthCollection::const_iterator mc = findTruth(event_data, d,mce);
1318 
1319  Identifier ID = d->identify();
1320  int nkine = 0;
1321 
1322  for(; mc!=mce; ++mc) {
1323 
1324  if( (*mc).first != ID ) return nkine;
1325 
1326  int k = (*mc).second.barcode(); if(k<=0) continue;
1327  Kine[nkine] = k; if(++nkine >= nmax) break;
1328  }
1329  return nkine;
1330 }
1331 
1333 // Test for cluster association with truth particles
1335 
1338 {
1339  PRD_MultiTruthCollection::const_iterator mce;
1340  PRD_MultiTruthCollection::const_iterator mc = findTruth(event_data,d,mce);
1341  return mc!=mce;
1342 }
1343 
1345 // Test detector element
1347 
1349 (const InDet::TrackClusterAssValidation::EventData_t &event_data, int K,const Trk::PrepRawData* d)
1350 {
1351  std::multimap<int,const Trk::PrepRawData*>::const_iterator k = event_data.m_kinecluster.find(K);
1352  for(; k!=event_data.m_kinecluster.end(); ++k) {
1353 
1354  if((*k).first!= K) return false;
1355  if(d->detectorElement()==(*k).second->detectorElement()) return true;
1356  }
1357  return false;
1358 }
1359 
1361 // Test detector element
1363 
1365 (const InDet::TrackClusterAssValidation::EventData_t &event_data, int K,const Trk::SpacePoint* sp)
1366 {
1367  const Trk::PrepRawData* p1 = sp->clusterList().first;
1368  const Trk::PrepRawData* p2 = sp->clusterList().second;
1369 
1370  std::multimap<int,const Trk::SpacePoint*>::const_iterator k = event_data.m_kinespacepoint.find(K);
1371 
1372  if(!p2) {
1373 
1374  for(; k!=event_data.m_kinespacepoint.end(); ++k) {
1375  if((*k).first!= K) return false;
1376 
1377  const Trk::PrepRawData* n1 = (*k).second->clusterList().first ;
1378  const Trk::PrepRawData* n2 = (*k).second->clusterList().second;
1379 
1380  if(p1->detectorElement() == n1->detectorElement()) return true;
1381  if(!n2) continue;
1382  if(p1->detectorElement() == n2->detectorElement()) return true;
1383  }
1384  return false;
1385  }
1386 
1387  for(; k!=event_data.m_kinespacepoint.end(); ++k) {
1388  if((*k).first!= K) return false;
1389 
1390  const Trk::PrepRawData* n1 = (*k).second->clusterList().first ;
1391  const Trk::PrepRawData* n2 = (*k).second->clusterList().second;
1392 
1393  if(p1->detectorElement() == n1->detectorElement()) return true;
1394  if(p2->detectorElement() == n1->detectorElement()) return true;
1395  if(!n2) continue;
1396  if(p1->detectorElement() == n2->detectorElement()) return true;
1397  if(p2->detectorElement() == n2->detectorElement()) return true;
1398  }
1399  return false;
1400 }
1401 
1403 // Dump information about no recontructed particles
1405 
1407 {
1408 
1410 
1411  auto p = event_data.m_particles[nc].begin(), pe =event_data.m_particles[nc].end();
1412  if(p==pe) continue;
1413 
1414  std::list<int>::const_iterator dif = event_data.m_difference[nc].begin();
1415 
1416  std::multimap<int,const Trk::PrepRawData*>::const_iterator c,ce = event_data.m_kinecluster.end();
1417 
1418  int n = 69-m_tracklocation[nc].key().size();
1419  std::string s1; for(int i=0; i<n; ++i) s1.append(" "); s1.append("|");
1420  std::stringstream out;
1421  out<<"|----------------------------------------------------------------------------------------|\n";
1422  out<<"| "<<m_tracklocation[nc]<<s1<<"\n";
1423  out<<"|----------------------------------------------------------------------------------------|\n";
1424  out<<"| # pdg kine Ncl Ntr Nsp Lose pT(MeV) rapidity radius z |\n";
1425  out<<"|----------------------------------------------------------------------------------------|\n";
1426 
1427  n = 0;
1428  for(; p!=pe; ++p) {
1429 
1430  int k = (*p).barcode();
1431 
1432  c = event_data.m_kinecluster.find(k); if(c==ce) continue;
1433  const Trk::PrepRawData* d = (*c).second;
1434 
1435  PRD_MultiTruthCollection::const_iterator mce;
1436  PRD_MultiTruthCollection::const_iterator mc = findTruth(event_data,d,mce);
1437 
1438  Identifier ID = d->identify();
1439  bool Q = false;
1440  for(; mc!=mce; ++mc) {
1441  if((*mc).first != ID) break;
1442  if((*mc).second.barcode()==k) {Q=true; break;}
1443  }
1444 
1445  if(!Q) continue;
1446 
1447  const HepMC::ConstGenParticlePtr pa = (*mc).second.cptr();
1448 
1449  double px = pa->momentum().px();
1450  double py = pa->momentum().py();
1451  double pz = pa->momentum().pz();
1452  double vx = pa->production_vertex()->position().x();
1453  double vy = pa->production_vertex()->position().y();
1454  double vz = pa->production_vertex()->position().z();
1455  double pt = std::sqrt(px*px+py*py);
1456  double t = std::atan2(pt,pz);
1457  double ra =-std::log(std::tan(.5*t));
1458  double r = std::sqrt(vx*vx+vy*vy);
1459  ++n;
1460  out<<"| "
1461  <<std::setw(4)<<n
1462  <<std::setw(6)<<pa->pdg_id()
1463  <<std::setw(10)<<HepMC::barcode(pa)
1464  <<std::setw(4)<<event_data.m_kinecluster .count(k)
1465  <<std::setw(4)<<event_data.m_kineclusterTRT.count(k)
1466  <<std::setw(4)<<event_data.m_kinespacepoint.count(k)
1467  <<std::setw(4)<<(*dif)
1468  <<std::setw(12)<<std::setprecision(5)<<pt
1469  <<std::setw(12)<<std::setprecision(5)<<ra
1470  <<std::setw(12)<<std::setprecision(5)<<r
1471  <<std::setw(12)<<std::setprecision(5)<<vz
1472  <<" |\n";
1473  ++dif;
1474 
1475  }
1476  out<<"|----------------------------------------------------------------------------------------|\n";
1477  ATH_MSG_INFO("\n"<<out.str());
1478  }
1479  return true;
1480 }
1481 
1483 // Cluster truth information
1485 
1486 PRD_MultiTruthCollection::const_iterator
1489  const Trk::PrepRawData* d,
1490  PRD_MultiTruthCollection::const_iterator& mce)
1491 {
1492  const InDet::SCT_Cluster * si = dynamic_cast<const InDet::SCT_Cluster*> (d);
1493  const InDet::PixelCluster * px = dynamic_cast<const InDet::PixelCluster*> (d);
1494  const InDet::TRT_DriftCircle* tr = dynamic_cast<const InDet::TRT_DriftCircle*>(d);
1495 
1496  PRD_MultiTruthCollection::const_iterator mc;
1497 
1498  if (px && event_data.m_truthPIX) {mc=event_data.m_truthPIX->find(d->identify()); mce=event_data.m_truthPIX->end();}
1499  else if(si && event_data.m_truthSCT) {mc=event_data.m_truthSCT->find(d->identify()); mce=event_data.m_truthSCT->end();}
1500  else if(tr && event_data.m_truthTRT) {mc=event_data.m_truthTRT->find(d->identify()); mce=event_data.m_truthTRT->end();}
1501  else {
1502  const PRD_MultiTruthCollection *truth[] {event_data. m_truthPIX,event_data.m_truthSCT, event_data.m_truthTRT};
1503  for (int i=0; i<3; i++) {
1504  if (truth[i]) {
1505  mce=truth[i]->end();
1506  return truth[i]->end();
1507  }
1508  }
1509  throw std::runtime_error("Neither Pixel, SCT nor TRT truth.");
1510  }
1511  return mc;
1512 }
1513 
1515 // Cluster truth information
1517 
1518 int InDet::TrackClusterAssValidation::charge(const InDet::TrackClusterAssValidation::EventData_t &event_data,std::pair<int,const Trk::PrepRawData*> pa,int& rap) const
1519 {
1520  int k = pa.first;
1521  const Trk::PrepRawData* d = pa.second;
1522  PRD_MultiTruthCollection::const_iterator mce;
1523  PRD_MultiTruthCollection::const_iterator mc = findTruth(event_data,d,mce);
1524 
1525  for(; mc!=mce; ++mc) {
1526  if((*mc).second.barcode()==k) {
1527 
1528  const HepMC::ConstGenParticlePtr pat = (*mc).second.cptr();
1529 
1530  rap = 0;
1531  double px = pat->momentum().px();
1532  double py = pat->momentum().py();
1533  double pz = pat->momentum().pz();
1534  double pt = std::sqrt(px*px+py*py) ;
1535  double t = std::atan2(pt,pz) ;
1536  double ra = std::abs(std::log(std::tan(.5*t)));
1537  // Forward
1538  if (ra > 3.0)
1539  rap = 3;
1540  else
1541  // other regions
1542  ra > 1.6 ? rap = 2 : ra > .8 ? rap = 1 : rap = 0;
1543 
1544  int pdg = pat->pdg_id();
1545  const HepPDT::ParticleData* pd = m_particleDataTable->particle(std::abs(pdg));
1546  if(!pd) return 0;
1547  double ch = pd->charge(); if(pdg < 0) ch = -ch;
1548  if(ch > .5) return 1;
1549  if(ch < -.5) return -1;
1550  return 0;
1551  }
1552  }
1553  return 0;
1554 }
Trk::SpacePoint::clusterList
const std::pair< const PrepRawData *, const PrepRawData * > & clusterList() const
return the pair of cluster pointers by reference
Definition: Tracking/TrkEvent/TrkSpacePoint/TrkSpacePoint/SpacePoint.h:127
xAOD::iterator
JetConstituentVector::iterator iterator
Definition: JetConstituentVector.cxx:68
Trk::SpacePoint
Definition: Tracking/TrkEvent/TrkSpacePoint/TrkSpacePoint/SpacePoint.h:35
beamspotman.r
def r
Definition: beamspotman.py:676
CalculateHighPtTerm.pT
pT
Definition: ICHEP2016/CalculateHighPtTerm.py:57
xAOD::eratio
setCharge setNTRTHiThresholdHits eratio
Definition: TrigElectron_v1.cxx:96
InDet::TrackClusterAssValidation::m_clcutTRT
unsigned int m_clcutTRT
Definition: TrackClusterAssValidation.h:68
ATH_MSG_FATAL
#define ATH_MSG_FATAL(x)
Definition: AthMsgStreamMacros.h:34
InDet::EventStat_t::m_nclustersNegEP
int m_nclustersNegEP
Definition: TrackClusterAssValidationUtils.h:111
InDet::TrackClusterAssValidation::m_useTRT
bool m_useTRT
Definition: TrackClusterAssValidation.h:59
InDet::TrackClusterAssValidation::initialize
StatusCode initialize()
Definition: TrackClusterAssValidation.cxx:87
DataModel_detail::const_iterator
Const iterator class for DataVector/DataList.
Definition: DVLIterator.h:82
sendEI_SPB.ch
ch
Definition: sendEI_SPB.py:35
PlotCalibFromCool.ft
ft
Definition: PlotCalibFromCool.py:329
python.SystemOfUnits.s
int s
Definition: SystemOfUnits.py:131
InDet::TrackClusterAssValidation::m_spacepointsOverlapname
SG::ReadHandleKey< SpacePointOverlapCollection > m_spacepointsOverlapname
Definition: TrackClusterAssValidation.h:79
InDet::TrackClusterAssValidation::EventData_t::m_nqtracks
int m_nqtracks
Definition: TrackClusterAssValidation.h:111
get_generator_info.result
result
Definition: get_generator_info.py:21
test_pyathena.px
px
Definition: test_pyathena.py:18
python.SystemOfUnits.m
int m
Definition: SystemOfUnits.py:91
python.PerfMonSerializer.p
def p
Definition: PerfMonSerializer.py:743
max
#define max(a, b)
Definition: cfImp.cxx:41
python.Constants.FATAL
int FATAL
Definition: Control/AthenaCommon/python/Constants.py:19
ATH_MSG_INFO
#define ATH_MSG_INFO(x)
Definition: AthMsgStreamMacros.h:31
ID
std::vector< Identifier > ID
Definition: CalibHitIDCheck.h:24
InDet::EventStat_t::m_nclustersNegBP
int m_nclustersNegBP
Definition: TrackClusterAssValidationUtils.h:109
AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::declareProperty
Gaudi::Details::PropertyBase & declareProperty(Gaudi::Property< T > &t)
Definition: AthCommonDataStore.h:145
InDet::TrackClusterAssValidation::m_useOutliers
bool m_useOutliers
Definition: TrackClusterAssValidation.h:60
hist_file_dump.d
d
Definition: hist_file_dump.py:137
InDet::TrackClusterAssValidation::EventData_t::m_nclusters
int m_nclusters
Definition: TrackClusterAssValidation.h:110
PRD_MultiTruthCollection
A PRD is mapped onto all contributing particles.
Definition: PRD_MultiTruthCollection.h:24
GenVertex.h
InDet::TrackClusterAssValidation::EventData_t::m_particles
std::vector< std::list< PartPropCache > > m_particles
Definition: TrackClusterAssValidation.h:124
test_pyathena.pt
pt
Definition: test_pyathena.py:11
python.AthDsoLogger.out
out
Definition: AthDsoLogger.py:71
ParticleTest.tp
tp
Definition: ParticleTest.py:25
TRT_PAI_gasdata::NC
const int NC
Number of levels for Carbon.
Definition: TRT_PAI_gasdata.h:237
InDet::TrackClusterAssValidation::charge
int charge(const InDet::TrackClusterAssValidation::EventData_t &event_data, std::pair< int, const Trk::PrepRawData * >, int &) const
Definition: TrackClusterAssValidation.cxx:1518
InDet::EventStat_t::m_nclustersNegBS
int m_nclustersNegBS
Definition: TrackClusterAssValidationUtils.h:110
AthCommonMsg< Gaudi::Algorithm >::msgLvl
bool msgLvl(const MSG::Level lvl) const
Definition: AthCommonMsg.h:30
dq_defect_virtual_defect_validation.d1
d1
Definition: dq_defect_virtual_defect_validation.py:79
SG::HandleKeyArray
Definition: StoreGate/StoreGate/HandleKeyArray.h:38
Trk::RIO_OnTrack
Definition: RIO_OnTrack.h:70
InDet::TrackClusterAssValidation::m_truth_locationSCT
SG::ReadHandleKey< PRD_MultiTruthCollection > m_truth_locationSCT
Definition: TrackClusterAssValidation.h:84
InDet::TrackClusterAssValidation::newSpacePointsEvent
void newSpacePointsEvent(const EventContext &ctx, InDet::TrackClusterAssValidation::EventData_t &event_data) const
Definition: TrackClusterAssValidation.cxx:879
ReadCellNoiseFromCoolCompare.s4
s4
Definition: ReadCellNoiseFromCoolCompare.py:381
InDet::TrackClusterAssValidation::m_clustersPixelname
SG::ReadHandleKey< SiClusterContainer > m_clustersPixelname
Definition: TrackClusterAssValidation.h:81
read_hist_ntuple.t
t
Definition: read_hist_ntuple.py:5
TrackClusterAssValidation.h
SG::VarHandleKey::key
const std::string & key() const
Return the StoreGate ID for the referenced object.
Definition: AthToolSupport/AsgDataHandles/Root/VarHandleKey.cxx:141
SG::ReadHandleKey
Property holding a SG store/key/clid from which a ReadHandle is made.
Definition: StoreGate/StoreGate/ReadHandleKey.h:39
InDet::TrackClusterAssValidation::m_statMutex
std::mutex m_statMutex
Definition: TrackClusterAssValidation.h:63
drawFromPickle.exp
exp
Definition: drawFromPickle.py:36
keylayer_zslicemap.se
se
Definition: keylayer_zslicemap.py:194
InDet::TrackClusterAssValidation::m_spacepointsSCTname
SG::ReadHandleKey< SpacePointContainer > m_spacepointsSCTname
Definition: TrackClusterAssValidation.h:77
LArG4AODNtuplePlotter.pe
pe
Definition: LArG4AODNtuplePlotter.py:116
GenParticle.h
InDet::TrackClusterAssValidation::EventData_t::m_kineclusterTRT
std::multimap< int, const Trk::PrepRawData * > m_kineclusterTRT
Definition: TrackClusterAssValidation.h:122
Trk::u
@ u
Enums for curvilinear frames.
Definition: ParamDefs.h:83
InDet::TRT_DriftCircle
Definition: TRT_DriftCircle.h:32
InDet::TrackClusterAssValidation::EventData_t::m_tracks
std::vector< std::multimap< int, int > > m_tracks
Definition: TrackClusterAssValidation.h:126
InDet::EventStat_t::m_nclustersNegES
int m_nclustersNegES
Definition: TrackClusterAssValidationUtils.h:112
drawFromPickle.atan
atan
Definition: drawFromPickle.py:36
InDet::EventStat_t::m_particleClusters
int m_particleClusters[50]
Definition: TrackClusterAssValidationUtils.h:100
AthenaPoolTestRead.sc
sc
Definition: AthenaPoolTestRead.py:27
python.iconfTool.models.loaders.level
level
Definition: loaders.py:20
AthReentrantAlgorithm
An algorithm that can be simultaneously executed in multiple threads.
Definition: AthReentrantAlgorithm.h:83
InDet::TrackClusterAssValidation::kine
int kine(const InDet::TrackClusterAssValidation::EventData_t &event_data, const Trk::PrepRawData *, const Trk::PrepRawData *, int *, int) const
Definition: TrackClusterAssValidation.cxx:1240
InDet::TrackClusterAssValidation::EventData_t::m_truthPIX
const PRD_MultiTruthCollection * m_truthPIX
Definition: TrackClusterAssValidation.h:118
InDet::EventStat_t::m_particleSpacePointsBTE
int m_particleSpacePointsBTE[50][4]
Definition: TrackClusterAssValidationUtils.h:103
InDet::TrackClusterAssValidation::EventData_t::m_nspacepoints
int m_nspacepoints
Definition: TrackClusterAssValidation.h:109
InDet::TrackClusterAssValidation::m_clcut
unsigned int m_clcut
Definition: TrackClusterAssValidation.h:67
InDet::TrackClusterAssValidation::m_spcut
unsigned int m_spcut
Definition: TrackClusterAssValidation.h:69
mc
Definition: mc.PG_single_nu_valid.py:1
TrigConf::MSGTC::Level
Level
Definition: Trigger/TrigConfiguration/TrigConfBase/TrigConfBase/MsgStream.h:21
InDet::TrackClusterAssValidation::m_SCTDetEleCollKey
SG::ReadCondHandleKey< InDetDD::SiDetectorElementCollection > m_SCTDetEleCollKey
Definition: TrackClusterAssValidation.h:89
InDet::TrackClusterAssValidation::kine0
static int kine0(const InDet::TrackClusterAssValidation::EventData_t &event_data, const Trk::PrepRawData *, int *, int)
Definition: TrackClusterAssValidation.cxx:1313
InDet::TrackClusterAssValidation::m_truth_locationTRT
SG::ReadHandleKey< PRD_MultiTruthCollection > m_truth_locationTRT
Definition: TrackClusterAssValidation.h:85
lumiFormat.i
int i
Definition: lumiFormat.py:92
InDet::TrackClusterAssValidation::m_truth_locationPixel
SG::ReadHandleKey< PRD_MultiTruthCollection > m_truth_locationPixel
Definition: TrackClusterAssValidation.h:83
InDet::TrackClusterAssValidation::noReconstructedParticles
bool noReconstructedParticles(const InDet::TrackClusterAssValidation::EventData_t &event_data) const
Definition: TrackClusterAssValidation.cxx:1406
Identifier
Definition: DetectorDescription/Identifier/Identifier/Identifier.h:32
beamspotman.n
n
Definition: beamspotman.py:731
InDet::TrackClusterAssValidation::EventData_t::m_trackCollectionStat
std::vector< TrackCollectionStat_t > m_trackCollectionStat
Definition: TrackClusterAssValidation.h:127
InDet::PartPropCache
Definition: TrackClusterAssValidationUtils.h:15
InDet::TrackClusterAssValidation::dumptools
MsgStream & dumptools(MsgStream &out, MSG::Level level) const
Definition: TrackClusterAssValidation.cxx:666
endmsg
#define endmsg
Definition: AnalysisConfig_Ntuple.cxx:63
EL::StatusCode
::StatusCode StatusCode
StatusCode definition for legacy code.
Definition: PhysicsAnalysis/D3PDTools/EventLoop/EventLoop/StatusCode.h:22
ATH_MSG_DEBUG
#define ATH_MSG_DEBUG(x)
Definition: AthMsgStreamMacros.h:29
InDet::TrackClusterAssValidation::~TrackClusterAssValidation
virtual ~TrackClusterAssValidation()
HepMC::barcode
int barcode(const T *p)
Definition: Barcode.h:16
AmgVector
AmgVector(4) T2BSTrackFilterTool
Definition: T2BSTrackFilterTool.cxx:114
InDet::TrackClusterAssValidation::m_useSCT
bool m_useSCT
Definition: TrackClusterAssValidation.h:58
TrackCollection.h
InDet::TrackClusterAssValidation::m_spacepointsPixelname
SG::ReadHandleKey< SpacePointContainer > m_spacepointsPixelname
Definition: TrackClusterAssValidation.h:78
Amg::pz
@ pz
Definition: GeoPrimitives.h:40
InDet::TrackClusterAssValidation::EventData_t::m_eventStat
EventStat_t m_eventStat
Definition: TrackClusterAssValidation.h:128
ATH_CHECK
#define ATH_CHECK
Definition: AthCheckMacros.h:40
InDet::TrackClusterAssValidation::EventData_t::m_clusterHandles
std::vector< std::unique_ptr< SG::VarHandleBase > > m_clusterHandles
Definition: TrackClusterAssValidation.h:114
InDet::TrackClusterAssValidation::m_pdg
int m_pdg
Definition: TrackClusterAssValidation.h:61
InDet::TrackClusterAssValidation::isTruth
static bool isTruth(const InDet::TrackClusterAssValidation::EventData_t &, const Trk::PrepRawData *)
Definition: TrackClusterAssValidation.cxx:1337
IdentifiableContainerMT::const_iterator
Definition: IdentifiableContainerMT.h:82
drawFromPickle.tan
tan
Definition: drawFromPickle.py:36
Trk::ParametersBase
Definition: ParametersBase.h:55
SG::VarHandleKey::initialize
StatusCode initialize(bool used=true)
If this object is used as a property, then this should be called during the initialize phase.
Definition: AthToolSupport/AsgDataHandles/Root/VarHandleKey.cxx:103
xAOD::double
double
Definition: CompositeParticle_v1.cxx:159
InDet::TrackClusterAssValidation::EventData_t::m_truthSCT
const PRD_MultiTruthCollection * m_truthSCT
Definition: TrackClusterAssValidation.h:119
InDet::TrackClusterAssValidation::m_ptcutmax
double m_ptcutmax
Definition: TrackClusterAssValidation.h:71
InDet::SCT_Cluster
Definition: InnerDetector/InDetRecEvent/InDetPrepRawData/InDetPrepRawData/SCT_Cluster.h:34
compute_lumi.denom
denom
Definition: compute_lumi.py:76
InDet::TrackClusterAssValidation::EventData_t::m_difference
std::vector< std::list< int > > m_difference
Definition: TrackClusterAssValidation.h:125
InDet::TrackClusterAssValidation::findTruth
static PRD_MultiTruthCollection::const_iterator findTruth(const InDet::TrackClusterAssValidation::EventData_t &event_data, const Trk::PrepRawData *, PRD_MultiTruthCollection::const_iterator &)
Definition: TrackClusterAssValidation.cxx:1488
InDet::TrackClusterAssValidation::tracksComparison
void tracksComparison(const EventContext &ctx, InDet::TrackClusterAssValidation::EventData_t &event_data) const
Definition: TrackClusterAssValidation.cxx:1081
lumiFormat.array
array
Definition: lumiFormat.py:98
Trk::PrepRawData
Definition: PrepRawData.h:62
InDet::TrackClusterAssValidation::m_tracklocation
SG::ReadHandleKeyArray< TrackCollection > m_tracklocation
Definition: TrackClusterAssValidation.h:76
Trk::MeasurementBase
Definition: MeasurementBase.h:58
dso-stats.pat
pat
Definition: dso-stats.py:39
InDet::EventStat_t::m_nclustersPosES
int m_nclustersPosES
Definition: TrackClusterAssValidationUtils.h:108
GetAllXsec.entry
list entry
Definition: GetAllXsec.py:132
InDet::TrackClusterAssValidation::EventData_t::m_trackcontainer
std::vector< SG::ReadHandle< TrackCollection > > m_trackcontainer
Definition: TrackClusterAssValidation.h:115
InDet::TrackClusterAssValidation::EventData_t::m_truthTRT
const PRD_MultiTruthCollection * m_truthTRT
Definition: TrackClusterAssValidation.h:120
HepMC::ConstGenParticlePtr
const GenParticle * ConstGenParticlePtr
Definition: GenParticle.h:38
Amg::py
@ py
Definition: GeoPrimitives.h:39
ReadCellNoiseFromCoolCompare.s3
s3
Definition: ReadCellNoiseFromCoolCompare.py:380
name
std::string name
Definition: Control/AthContainers/Root/debug.cxx:195
RIO_OnTrack.h
PixelClusterContainer.h
InDet::TrackClusterAssValidation::m_tcut
double m_tcut
Definition: TrackClusterAssValidation.h:73
InDet::TrackClusterAssValidation::EventData_t
Definition: TrackClusterAssValidation.h:91
InDetDD::SiDetectorElement
Definition: SiDetectorElement.h:109
SG::CondHandleKey::initialize
StatusCode initialize(bool used=true)
InDet::TrackClusterAssValidation::m_clustersTRTname
SG::ReadHandleKey< TRT_DriftCircleContainer > m_clustersTRTname
Definition: TrackClusterAssValidation.h:82
InDet::TrackClusterAssValidation::TrackClusterAssValidation
TrackClusterAssValidation(const std::string &name, ISvcLocator *pSvcLocator)
Definition: TrackClusterAssValidation.cxx:24
InDet::TrackClusterAssValidation::m_pixelDetEleCollKey
SG::ReadCondHandleKey< InDetDD::SiDetectorElementCollection > m_pixelDetEleCollKey
Definition: TrackClusterAssValidation.h:87
Trk::RIO_OnTrack::prepRawData
virtual const Trk::PrepRawData * prepRawData() const =0
returns the PrepRawData (also known as RIO) object to which this RIO_OnTrack is associated.
InDet::TrackClusterAssValidation::m_particleDataTable
const HepPDT::ParticleDataTable * m_particleDataTable
Definition: TrackClusterAssValidation.h:131
xAOD::k0
@ k0
for Fatras usage
Definition: TrackingPrimitives.h:202
Rtt_histogram.n1
n1
Definition: Rtt_histogram.py:21
InDet::PixelCluster
Definition: InnerDetector/InDetRecEvent/InDetPrepRawData/InDetPrepRawData/PixelCluster.h:49
InDet::EventStat_t::m_eventsBTE
int m_eventsBTE[4]
Definition: TrackClusterAssValidationUtils.h:98
InDet::TrackClusterAssValidation::efficiencyReconstruction
void efficiencyReconstruction(InDet::TrackClusterAssValidation::EventData_t &event_data) const
Definition: TrackClusterAssValidation.cxx:1192
python.compareTCTs.ratio
ratio
Definition: compareTCTs.py:295
python.CaloScaleNoiseConfig.str
str
Definition: CaloScaleNoiseConfig.py:78
InDet::EventStat_t::m_nclustersPosBP
int m_nclustersPosBP
Definition: TrackClusterAssValidationUtils.h:105
InDet::TrackClusterAssValidation::isTheSameDetElement
static bool isTheSameDetElement(const InDet::TrackClusterAssValidation::EventData_t &event_data, int, const Trk::PrepRawData *)
Definition: TrackClusterAssValidation.cxx:1349
TRT_PAI_physicsConstants::mb
const double mb
1mb to cm2
Definition: TRT_PAI_physicsConstants.h:15
InDet::EventStat_t::m_events
int m_events
Definition: TrackClusterAssValidationUtils.h:95
Pythia8_RapidityOrderMPI.val
val
Definition: Pythia8_RapidityOrderMPI.py:14
InDet::TrackClusterAssValidation::m_rapcut
double m_rapcut
Definition: TrackClusterAssValidation.h:72
InDet::EventStat_t
Definition: TrackClusterAssValidationUtils.h:93
InDet::TrackClusterAssValidation::finalize
StatusCode finalize()
Definition: TrackClusterAssValidation.cxx:230
dq_defect_virtual_defect_validation.d2
d2
Definition: dq_defect_virtual_defect_validation.py:81
DEBUG
#define DEBUG
Definition: page_access.h:11
ReadCellNoiseFromCoolCompare.s2
s2
Definition: ReadCellNoiseFromCoolCompare.py:379
InDet::TrackClusterAssValidation::EventData_t::m_nclustersTRT
int m_nclustersTRT
Definition: TrackClusterAssValidation.h:112
InDet::EventStat_t::m_eventsPOS
int m_eventsPOS
Definition: TrackClusterAssValidationUtils.h:96
AthCommonMsg< Gaudi::Algorithm >::msg
MsgStream & msg() const
Definition: AthCommonMsg.h:24
InDet::TrackClusterAssValidation::m_usePIX
bool m_usePIX
Definition: TrackClusterAssValidation.h:57
python.CaloCondTools.log
log
Definition: CaloCondTools.py:20
InDet::EventStat_t::m_eventsNEG
int m_eventsNEG
Definition: TrackClusterAssValidationUtils.h:97
SCT_ClusterContainer.h
InDet::TrackClusterAssValidation::EventData_t::m_spacepointsOverlap
std::unique_ptr< SG::ReadHandle< SpacePointOverlapCollection > > m_spacepointsOverlap
Definition: TrackClusterAssValidation.h:117
LArNewCalib_DelayDump_OFC_Cali.idx
idx
Definition: LArNewCalib_DelayDump_OFC_Cali.py:69
extractSporadic.q
list q
Definition: extractSporadic.py:98
python.SystemOfUnits.ns
int ns
Definition: SystemOfUnits.py:130
InDet::TrackClusterAssValidation::qualityTracksSelection
int qualityTracksSelection(InDet::TrackClusterAssValidation::EventData_t &event_data) const
Definition: TrackClusterAssValidation.cxx:968
dqt_zlumi_alleff_HIST.eff
int eff
Definition: dqt_zlumi_alleff_HIST.py:113
str
Definition: BTagTrackIpAccessor.cxx:11
InDet::TrackClusterAssValidation::execute
StatusCode execute(const EventContext &ctx) const
Definition: TrackClusterAssValidation.cxx:164
InDet::TrackClusterAssValidation::dumpevent
static MsgStream & dumpevent(MsgStream &out, const InDet::TrackClusterAssValidation::EventData_t &event_data)
Definition: TrackClusterAssValidation.cxx:736
InDet::EventStat_t::m_nclustersPosEP
int m_nclustersPosEP
Definition: TrackClusterAssValidationUtils.h:107
drawFromPickle.sin
sin
Definition: drawFromPickle.py:36
beamspotnt.nt
def nt
Definition: bin/beamspotnt.py:1063
ReadHandle.h
Handle class for reading from StoreGate.
rp
ReadCards * rp
Definition: IReadCards.cxx:26
python.IoTestsLib.w
def w
Definition: IoTestsLib.py:200
covarianceTool.mc
mc
Definition: covarianceTool.py:554
InDet::TrackClusterAssValidation::m_rmax
double m_rmax
Definition: TrackClusterAssValidation.h:75
python.CaloScaleNoiseConfig.ts
ts
Definition: CaloScaleNoiseConfig.py:86
InDet::TrackClusterAssValidation::EventData_t::m_kinespacepoint
std::multimap< int, const Trk::SpacePoint * > m_kinespacepoint
Definition: TrackClusterAssValidation.h:123
InDet::EventStat_t::m_nclustersPosBS
int m_nclustersPosBS
Definition: TrackClusterAssValidationUtils.h:106
python.compressB64.c
def c
Definition: compressB64.py:93
InDet::TrackClusterAssValidation::EventData_t::m_kinecluster
std::multimap< int, const Trk::PrepRawData * > m_kinecluster
Definition: TrackClusterAssValidation.h:121
nmax
const int nmax(200)
InDet::TrackClusterAssValidation::m_rmin
double m_rmin
Definition: TrackClusterAssValidation.h:74
InDet::SiCluster
Definition: InnerDetector/InDetRecEvent/InDetPrepRawData/InDetPrepRawData/SiCluster.h:40
LheEventFiller_Common.ef
ef
Definition: SFGen_i/share/common/LheEventFiller_Common.py:7
plotBeamSpotMon.nc
int nc
Definition: plotBeamSpotMon.py:83
Trk::TrackStateOnSurface::Measurement
@ Measurement
This is a measurement, and will at least contain a Trk::MeasurementBase.
Definition: TrackStateOnSurface.h:101
InDet::TrackClusterAssValidation::m_ptcut
double m_ptcut
Definition: TrackClusterAssValidation.h:70
InDet::EventStat_t::m_particleClustersBTE
int m_particleClustersBTE[50][4]
Definition: TrackClusterAssValidationUtils.h:102
InDet::TrackClusterAssValidation::newClustersEvent
void newClustersEvent(const EventContext &ctx, InDet::TrackClusterAssValidation::EventData_t &event_data) const
Definition: TrackClusterAssValidation.cxx:765
InDet::TrackClusterAssValidation::m_clustersSCTname
SG::ReadHandleKey< SiClusterContainer > m_clustersSCTname
Definition: TrackClusterAssValidation.h:80
fitman.k
k
Definition: fitman.py:528
Trk::PrepRawData::detectorElement
virtual const TrkDetElementBase * detectorElement() const =0
return the detector element corresponding to this PRD The pointer will be zero if the det el is not d...
InDet::TrackClusterAssValidation::EventData_t::m_spacePointContainer
std::vector< SG::ReadHandle< SpacePointContainer > > m_spacePointContainer
Definition: TrackClusterAssValidation.h:116
InDet::EventStat_t::m_particleSpacePoints
int m_particleSpacePoints[50]
Definition: TrackClusterAssValidationUtils.h:101