ATLAS Offline Software
PurityAnalysis.cxx
Go to the documentation of this file.
1 
11 // #include "TrigInDetAnalysisExample/PurityAnalysis.h"
12 #include "PurityAnalysis.h"
13 
14 #include "TrigInDetAnalysis/TrackFilter.h"
15 
16 
17 #include "TF1.h"
18 
19 extern int r;
20 extern int ev;
21 extern int lb;
22 extern int ts;
23 
24 extern double a0;
25 
26 extern bool hipt;
27 
28 
29 
30 
31 
33 
34 void PurityAnalysis::initialise() {
35 
36  std::cout << "PurityAnalysis::initialise() " << name() << std::endl;
37 
38  //+++ pT ranges
39  // double tmp_maxPt = 50000.;
40  // double tmp_absResPt = 0.0005;
41  // double tmp_maxPt = 50.;
42  // double tmp_absResPt = 0.5;
43 
44  const int pTResBins = 100;
45 
46  //+++ Eta ranges
47  double tmp_maxEta = 3.;
48  double tmp_absResEta = 0.04; // 0.0005;
49 
50  //+++ Phi ranges
51  double tmp_maxPhi = 3.142;
52  double tmp_absResPhi = 0.02; // 0.0001;
53 
54 
55  const int etaBins = 60;
56  const int etaResBins = 300;
57 
58  const int phiBins = 36;
59  const int phiResBins = 100;
60 
61  const int zBins = 100;
62  const double zMax = 400;
63 
64  const int zresBins = 100;
65  const double zresMax = 10;
66 
67  const int d0Bins = 100;
68  const double d0Max = 10;
69 
70  const int d0resBins = 100;
71  const double d0resMax = 10;
72 
73  // beamspot corrected position
74 
75  const int a0Bins = 100;
76  const double a0Max = 10;
77 
78  const int a0resBins = 100;
79  const double a0resMax = 5;
80 
81  //+++ Book histograms
82 
83  // calculate a logarithmic binning in pt
84 
85  const int ptnbins = 20;
86  double ptbinlims[ptnbins+1];
87  for ( int i=0 ; i<=ptnbins ; i++ ) { ptbinlims[i] = std::pow(10, 3.0*i/ptnbins+2)/1000; }
88 
89 
90  TDirectory* dir = gDirectory;
91 
92  std::cout << "PurityAnalysis::initialize() Directory " << gDirectory->GetName() << " " << name() << std::endl;
93 
94  m_dir = new TIDDirectory(name());
95  m_dir->push();
96 
97  // TIDDirectory d("histos");
98  // d.push();
99 
100  std::cout << "PurityAnalysis::initialize() Directory " << gDirectory->GetName() << " package directory, " << name() << std::endl;
101 
102  // int Nptbins = 6;
103  // double _ptlims[7] = { 0, 500, 1000, 1500, 2000, 5000, 10000 };
104 
105  // TH2F* effpt2d = new TH2F("pteta2d", "pteta", Nptbins, _ptlims, 40, -tmp_maxEta, tmp_maxEta );
106  // TH2F* effeta2d = new TH2F("etapt2d", "pteta", ptnbins, ptbinlims, 6, -tmp_maxEta, tmp_maxEta );
107 
108  // eff_pteta = new Efficiency2D( effpt2d, "pteta" );
109  // eff_etapt = new Efficiency2D( effeta2d, "etapt" );
110 
111  // effpt2d->SetDirectory(0);
112  // effeta2d->SetDirectory(0);
113 
114  // delete effpt2d;
115  // delete effeta2d;
116 
117  Efficiency1D* heff[8];
118  Efficiency1D* hpurity[6];
119 
120  addHistogram( m_hchi2=new TH1F("chi2", "chi2", 100, 0, 20) );
121 
122  // "reference" quantities
123  addHistogram( new TH1F( "pT", "pT", ptnbins, ptbinlims ) );
124  addHistogram( new TH1F( "eta", "eta", etaBins, -tmp_maxEta, tmp_maxEta ) );
125  addHistogram( new TH1F( "phi", "phi", phiBins, -tmp_maxPhi, tmp_maxPhi ) );
126  addHistogram( new TH1F( "z0", "z0", zBins, -zMax, zMax ) );
127  addHistogram( new TH1F( "d0", "d0", d0Bins, -d0Max, d0Max ) );
128  addHistogram( new TH1F( "a0", "a0", a0Bins, -a0Max, a0Max ) );
129 
130  // efficienies and purities
131  heff[0] = new Efficiency1D( find("pT"), "pT_eff" );
132  heff[1] = new Efficiency1D( find("eta"), "eta_eff" );
133  heff[2] = new Efficiency1D( find("phi"), "phi_eff" );
134  heff[3] = new Efficiency1D( find("z0"), "z0_eff" );
135  heff[4] = new Efficiency1D( find("d0"), "d0_eff" );
136  heff[5] = new Efficiency1D( find("a0"), "a0_eff" );
137 
138  heff[6] = new Efficiency1D( find("pT"), "pTm_eff" );
139  heff[7] = new Efficiency1D( find("pT"), "pTp_eff" );
140 
141  m_eff_pt = heff[0];
142  m_eff_eta = heff[1];
143  m_eff_phi = heff[2];
144  m_eff_z0 = heff[3];
145  m_eff_d0 = heff[4];
146  m_eff_a0 = heff[5];
147 
148  m_eff_ptm = heff[6];
149  m_eff_ptp = heff[7];
150 
151  // addHistogram ( hDeltaR = new TH1F("DeltaR", "DeltaR", 100, 0, 0.1 ) );
152  addHistogram ( m_hDeltaR = new TH1F("DeltaR", "DeltaR", 100, 0, 0.2 ) );
153 
154  hpurity[0] = new Efficiency1D( find("pT"), "pT_pur" );
155  hpurity[1] = new Efficiency1D( find("eta"), "eta_pur" );
156  hpurity[2] = new Efficiency1D( find("phi"), "phi_pur" );
157  hpurity[3] = new Efficiency1D( find("z0"), "z0_pur" );
158  hpurity[4] = new Efficiency1D( find("d0"), "d0_pur" );
159  hpurity[5] = new Efficiency1D( find("a0"), "a0_pur" );
160 
161  m_purity_pt = hpurity[0];
162  m_purity_eta = hpurity[1];
163  m_purity_phi = hpurity[2];
164  m_purity_z0 = hpurity[3];
165  m_purity_d0 = hpurity[4];
166  m_purity_a0 = hpurity[5];
167 
168  // "test" quantities
169  addHistogram( new TH1F( "pT_rec", "pT_rec", ptnbins, ptbinlims ) );
170  addHistogram( new TH1F( "eta_rec", "eta_rec", etaBins, -tmp_maxEta, tmp_maxEta ) );
171  addHistogram( new TH1F( "phi_rec", "phi_rec", phiBins, -tmp_maxPhi, tmp_maxPhi ) );
172  addHistogram( new TH1F( "z0_rec", "z0_rec", zBins, -zMax, zMax ) );
173  addHistogram( new TH1F( "d0_rec", "d0_rec", d0Bins, -d0Max, d0Max ) );
174  addHistogram( new TH1F( "a0_rec", "a0_rec", a0Bins, -a0Max, a0Max ) );
175 
176  // resolutions
177  // addHistogram( new TH1F( "pT_res", "pT_res", 2*pTResBins, -2*tmp_absResPt, 2*tmp_absResPt ) );
178  addHistogram( new TH1F( "pT_res", "pT_res", 4*pTResBins, -0.1, 0.1 ) );
179  addHistogram( new TH1F( "spT_res", "spT_res", 4*pTResBins, -0.1, 0.1 ) );
180  addHistogram( new TH1F( "ipT_res", "pT_res", 4*pTResBins, -0.4, 0.4 ) );
181  addHistogram( new TH1F( "eta_res", "eta_res", etaResBins, -2*tmp_absResEta, 2*tmp_absResEta ) );
182  addHistogram( new TH1F("etai_res", "etai_res", 1000, -0.04, 0.04 ) );
183  addHistogram( new TH1F( "phi_res", "phi_res", 2*phiResBins, -2*tmp_absResPhi, 2*tmp_absResPhi ) );
184  addHistogram( new TH1F( "z0_res", "z0_res", zresBins, -zresMax, zresMax ) );
185  addHistogram( new TH1F( "d0_res", "d0_res", 4*d0resBins, -0.5*d0resMax, 0.5*d0resMax ) );
186  addHistogram( new TH1F( "a0_res", "a0_res", 4*a0resBins, -a0resMax, a0resMax ) );
187 
188  std::cout << "booked" << std::endl;
189 
190 
191  // hit occupancies
192 
193  int NHits = 40;
194  int Ntracks = 50;
195 
196  addHistogram( new TH1F( "nsct", "nsct", NHits, -0.5, float(NHits-0.5) ) );
197  addHistogram( new TH1F( "nsct_rec", "nsct_rec", NHits, -0.5, float(NHits-0.5) ) );
198 
199  addHistogram( new TH1F( "npix", "npix", NHits, -0.5, float(NHits-0.5) ) );
200  addHistogram( new TH1F( "npix_rec", "npix_rec", NHits, -0.5, float(NHits-0.5) ) );
201 
202  addHistogram( new TH1F( "ntrt", "ntrt", NHits, -0.5, float(NHits-0.5) ) );
203  addHistogram( new TH1F( "ntrt_rec", "ntrt_rec", NHits, -0.5, float(NHits-0.5) ) );
204 
205  addHistogram( new TH1F( "nstraw", "nstraw", NHits*4, -0.5, float(4*NHits-0.5) ) );
206  addHistogram( new TH1F( "nstraw_rec", "nstraw_rec", NHits*4, -0.5, float(4*NHits-0.5) ) );
207 
208  addHistogram( new TH1F( "ntracks", "ntracks", Ntracks, -0.5, float(Ntracks+0.5) ) );
209  addHistogram( new TH1F( "ntracks_rec", "ntracks_rec", Ntracks, -0.5, float(Ntracks+0.5) ) );
210 
211 
212  m_dir->pop();
213 
214  dir->cd();
215 
216  // std::cout << "initialize() Directory " << gDirectory->GetName() << " on leaving" << std::endl;
217 
218 }
219 
220 
221 
222 
223 
224 
226 
227 void PurityAnalysis::finalise() {
228 
229  std::cout << "PurityAnalysis::finalise() " << name()
230  << "\tNreco " << m_Nreco
231  << "\tNref " << m_Nref
232  << "\tNmatched " << m_Nmatched << " tracks"
233  << std::endl;
234 
235  // if ( m_Nreco==0 ) return;
236 
237  // TIDDirectory d( name() );
238  // d.push();
239 
240  m_dir->push();
241 
242  //std::map<std::string, TH1F*>::iterator hitr=m_histos.begin();
243  //std::map<std::string, TH1F*>::iterator hend=m_histos.end();
244  // for ( ; hitr!=hend ; hitr++ ) hitr->second->Write();
245 
246  // std::cout << "DBG >" << m_eff_pt->Hist()->GetName() << "< DBG" << std::endl;
247 
248  Efficiency1D* heff[8] = { m_eff_pt, m_eff_eta, m_eff_phi, m_eff_z0, m_eff_d0, m_eff_a0, m_eff_ptm, m_eff_ptp };
249  for ( int i=8 ; i-- ; ) { heff[i]->finalise(); } // heff[i]->Hist()->Write(); }
250 
251  // std::cout << "DBG >" << m_purity_pt->Hist()->GetName() << "< DBG" << std::endl;
252 
253  // m_eff_pteta->finalise(); m_eff_pteta->Write("eta_efficiency_binned_pt", "x");
254  // for ( int i=1 ; i<=m_eff_pteta->GetNbinsX() ; i++ ) {
255  // TH1F* h = m_eff_pteta->SliceX(i);
256  // }
257 
258  // m_eff_etapt->finalise(); m_eff_etapt->Write("pt_efficieny_binned_eta", "y");
259  // for ( int i=1 ; i<=m_eff_etapt->GetNbinsY() ; i++ ) {
260  // TH1F* h = m_eff_etapt->SliceY(i);
261  // }
262 
263  Efficiency1D* hpurity[6] = { m_purity_pt, m_purity_eta, m_purity_phi, m_purity_z0, m_purity_d0, m_purity_a0 };
264  for ( int i=6 ; i-- ; ) { hpurity[i]->finalise(); } // hpurity[i]->Hist()->Write(); }
265 
266  // d.pop();
267  m_dir->pop();
268 
269 }
270 
271 #if 0
272 int duff_counter = 0;
273 
274 void duff(const std::string& label="") {
275  std::cout << "counter " << duff_counter << label << std::endl;
276  duff_counter++;
277 }
278 #endif
279 
281 
282 void PurityAnalysis::execute(const std::vector<TIDA::Track*>& reftracks,
283  const std::vector<TIDA::Track*>& testtracks,
284  TrackAssociator* matcher )
285 {
286  if ( m_print ) std::cout << "PurityAnalysis::execute() \t " << name()
287  << "\tref " << reftracks.size()
288  << "\ttest " << testtracks.size() << std::endl;
289 
290  // std::cout << "\tx " << m_xBeamReference << "\ty " << m_yBeamReference
291  // << "\tx " << m_xBeamTest << "\ty " << m_yBeamTest << std::endl;
292 
293 
294  // std::cout << "PurityAnalysis (resolutions really) filling " << std::endl;
295 
296  // should have these as a class variable
297  static std::string varName[10] = { "pT", "eta", "phi", "z0", "d0", "a0", "nsct", "npix", "ntrt", "nstraw" };
298 
299  // std::cout << "PurityAnalysis ref size " << reftracks.size() << "\ttest size " << testtracks.size() << std::endl;
300 
301  std::map<std::string, TH1F*>::iterator hmitr = m_histos.find("ntracks");
302  if ( hmitr!=m_histos.end() ) hmitr->second->Fill( reftracks.size() );
303 
304  hmitr = m_histos.find("ntracks_rec");
305  if ( hmitr!=m_histos.end() ) hmitr->second->Fill( testtracks.size() );
306 
307  bool dump = false;
308 
309  m_Nreco += testtracks.size();
310  m_Nref += reftracks.size();
311 
312  // std::cout << "PurityAnalysis ref tracks " << std::endl;
313 
314  m_Nref = 0;
315  for ( int i=reftracks.size() ; i-- ; ) {
316  double phit = reftracks[i]->phi();
317  double a0t = reftracks[i]->a0() + sin(phit)*m_xBeamReference - cos(phit)*m_yBeamReference;
318  if ( std::fabs(a0t)<a0 ) m_Nref++;
319  }
320 
321  // if ( testtracks.size() ) std::cout << "NTRACKS " << testtracks.size() << std::endl;
322 
323  for ( int i=reftracks.size() ; i-- ; ) {
324 
325  double pTt = reftracks[i]->pT()/1000;
326  double etat = reftracks[i]->eta();
327  double phit = reftracks[i]->phi();
328  double z0t = reftracks[i]->z0() + std::cos(phit)*m_xBeamReference + std::sin(phit)*m_yBeamReference;
329  double d0t = reftracks[i]->a0();
330  // this will be changed when we know the beam spot position
331  // double a0t = reftracks[i]->a0() + sin(phit)*m_xBeam - cos(phit)*m_yBeam;
332  double a0t = reftracks[i]->a0() + std::sin(phit)*m_xBeamReference - std::cos(phit)*m_yBeamReference;
333 
334  if ( std::fabs(a0t)>a0 ) continue;
335 
336  double chi2t = reftracks[i]->chi2();
337 
338  m_hchi2->Fill( chi2t );
339 
340  double nsctt = reftracks[i]->sctHits();
341  double npixt = reftracks[i]->pixelHits();
342 
343  double ntrtt = reftracks[i]->trHits();
344  double nstrawt = reftracks[i]->strawHits();
345 
346  // double ts_scale = (ts-1260400000)*3000.0/(1260700000-1260400000);
347 
348  // std::cout << "Fill h2 " << " " << h2m << " " << *reftracks[i] << std::endl;
349 
350  const TIDA::Track* matchedreco = matcher->matched(reftracks[i]);
351 
352  // std::cout << "\t\tPurityAnalysis " << name() << "\t" << i << " " << *reftracks[i] << " -> ";
353 
354  // raw reference track distributions
355  double vpart[10] = { std::fabs(pTt), etat, phit, z0t, d0t, a0t, nsctt, npixt, ntrtt, nstrawt };
356  for ( int it=0 ; it<10 ; it++ ) {
357  // std::string hname = varName[it];
358  // std::map<std::string, TH1F*>::iterator hmitr = m_histos.find(hname);
359  // if ( hmitr!=m_histos.end() ) hmitr->second->Fill( vpart[it] );
360 
361  if ( TH1F* hptr = find( varName[it] ) ) hptr->Fill( vpart[it] );
362  else std::cerr << "hmmm histo " << varName[it] << " not found" << std::endl;
363 
364  }
365 
366 
367  if ( matchedreco ) {
368 
369  // efficiency histos
370  m_eff_pt->Fill(std::fabs(pTt));
371  m_eff_z0->Fill(z0t);
372  m_eff_eta->Fill(etat);
373  m_eff_phi->Fill(phit);
374  m_eff_d0->Fill(d0t);
375  m_eff_a0->Fill(a0t);
376 
377  // signed pT
378  if ( pTt<0 ) m_eff_ptm->Fill(std::fabs(pTt));
379  else m_eff_ptp->Fill(std::fabs(pTt));
380 
381  m_Nmatched++;
382 
383  // residual histos
384  double pTr = matchedreco->pT()/1000;
385  double etar = matchedreco->eta();
386  double phir = matchedreco->phi();
387  double z0r = matchedreco->z0() + std::cos(phir)*m_xBeamTest + std::sin(phir)*m_yBeamTest; ;
388  double d0r = matchedreco->a0();
389  double a0r = matchedreco->a0() + sin(phir)*m_xBeamTest - cos(phir)*m_yBeamTest; // this will be changed when we know the beam spot position
390 
391  // double nsctr = matchedreco->sctHits();
392  // double npixr = matchedreco->pixelHits();
393 
394  // double ntrtr = matchedreco->trHits();
395  // double nstrawr = matchedreco->strawHits();
396 
397 #if 0
398  // if ( m_print ) std::cout << "PurityAnalysis::execute() \t " << name() << "\t" << i << " "
399  if ( m_print ) std::cout << "PurityAnalysis::execute() \t\t" << i << " "
400  << *reftracks[i] << " -> " << *matchedreco << "\t"
401  << pTr << " " << pTt << " " << d0r << " " << d0t << std::endl;
402 #endif
403 
404  double vres[6] = { 1.0/std::fabs(pTr)-1.0/std::fabs(pTt), etar-etat, phir-phit, z0r-z0t, d0r-d0t, a0r-a0t };
405  for ( int it=0 ; it<6 ; it++ ) {
406  if ( it==0 ) {
407  find("ipT_res")->Fill( vres[0] );
408  find("spT_res")->Fill( 1.0/pTr-1.0/pTt );
409  }
410  if ( TH1F* hptr = find(varName[it]+"_res") ) hptr->Fill( vres[it] );
411  else std::cerr << "hmmm histo " << varName[it]+"_res" << " not found" << std::endl;
412  }
413 
414  if ( TH1F* hptr = find("etai_res") ) hptr->Fill( etat-etar );
415 
416 
417  double Delphi = phit-phir;
418  double Deleta = etat-etar;
419 
420  if ( Delphi<-M_PI ) Delphi+=2*M_PI;
421  if ( Delphi>M_PI ) Delphi -=2*M_PI;
422 
423  double DeltaR = std::sqrt(Delphi*Delphi+Deleta*Deleta);
424 
425  m_hDeltaR->Fill(DeltaR);
426 
427 
428  // in this loop over the reference tracks, could fill efficiency
429  // histograms
430 
431  // m_eff_pteta->Fill( pTt, etat );
432  // m_eff_etapt->Fill( pTt, etat );
433 
435 
436 #if 0
437  // raw test track distributions
438  double vpart[10] = { std::fabs(pTr), etar, phir, z0r, d0r, a0r, nsctr, npixr, ntrtr, nstrawr };
439  for ( int it=0 ; it<10 ; it++ ) {
440  // std::string hname = name()+"_"+varName[it]+"_rec";
441  // std::string hname = varName[it]+"_rec";
442  // std::map<std::string, TH1F*>::iterator hmitr = m_histos.find(hname);
443  // if ( hmitr!=m_histos.end() ) hmitr->second->Fill( vpar[it] );
444  // else std::cerr << "hmmm histo " << hname << " not found" << std::endl;
445  if ( TH1F* hptr = find(varName[it]+"_rec") ) hptr->Fill( vpart[it] );
446  else std::cerr << "hmmm histo " << varName[it]+"_rec" << " not found" << std::endl;
447  }
448 #endif
449 
450  }
451  else {
452  // fill efficiencies with unmatched histos
453  // std::cout << "NULL" << std::endl;
454  m_eff_pt->FillDenom(std::fabs(pTt));
455  m_eff_z0->FillDenom(z0t);
456  m_eff_eta->FillDenom(etat);
457  m_eff_phi->FillDenom(phit);
458  m_eff_d0->FillDenom(d0t);
459  m_eff_a0->FillDenom(a0t);
460 
461  // signed pT
462  if ( pTt<0 ) m_eff_ptm->FillDenom(std::fabs(pTt));
463  else m_eff_ptp->FillDenom(std::fabs(pTt));
464 
465  dump = false;
466 
467 #if 0
468  if ( std::fabs(pTt)>4 ) {
469  dump = true;
470 
471  hipt = true;
472  std::cout << mname << "\tMISSING TRACK run " << r << "\tevent " << ev
473  << "\tlb " << lb << "\t" << *reftracks[i];
474 
475  if ( std::fabs(pTt)>=30 ) std::cout << "\tvery high pt";
476  if ( std::fabs(pTt)>4 &&
477  std::fabs(pTt)<30 ) std::cout << "\t high pt";
478 
479  std::cout << std::endl;
480 
481  }
482 #endif
483 
484 
485  // m_eff_pteta->FillDenom( pTt, etat );
486  // m_eff_etapt->FillDenom( pTt, etat );
487 
488  }
489 
490  }
491 
492  // for fake/purity histograms, loop over the test tracks
493  // and get the corresponding matched reference tracks from the
494  // reverse map in the TrackAscociator class - revmatched()
495 
496  static int icount = 0;
497 
498  // if ( icount%1000 ) std::cout << "chain " << name() << "\t " << m_Nreco << " tracks" << std::endl;
499  // if ( icount%1000 )
500  if ( m_print ) std::cout << "PurityAnalysis::execute() \t " << name() << "\t " << icount << " events\t " << testtracks.size() << " tracks (" << m_Nreco << ")" << "\n---------------" << std::endl;
501 
502  icount++;
503 
504  for ( int i=testtracks.size() ; i-- ; ) {
505 
506  // std::cout << "\t\tPurityAnalysis purity " << name() << "\t" << i << " " << *testtracks[i] << " -> ";
507 
508  // double pTr = std::fabs(testtracks[i]->pT());
509  double pTr = testtracks[i]->pT()/1000;
510  double etar = testtracks[i]->eta();
511  double phir = testtracks[i]->phi();
512  double z0r = testtracks[i]->z0();
513  double d0r = testtracks[i]->a0();
514  double a0r = testtracks[i]->a0() + sin(phir)*m_xBeamTest - cos(phir)*m_yBeamTest; // this will be changed when we know the beam spot position
515  // double a0rp = testtracks[i]->a0() - sin(phir)*m_xBeam - cos(phir)*m_yBeam; // this will be changed when we know the beam spot position
516 
517  // std::cout << "d0 " << d0r << "\tphi " << phir << "\tx " << m_xBeamTest << "\ty " << m_yBeamTest << std::endl;
518 
519  double nsctr = testtracks[i]->sctHits();
520  double npixr = testtracks[i]->pixelHits();
521 
522  double ntrtr = testtracks[i]->trHits();
523  double nstrawr = testtracks[i]->strawHits();
524 
525 
526  const TIDA::Track* matchedref = matcher->revmatched(testtracks[i]);
527 
528  // if ( matchedref ) std::cout << *matchedref << std::endl;
529  // else std::cout << "NULL" << std::endl;
530 
531 #if 1
532  // raw test track distributions
533  double vpart[10] = { std::fabs(pTr), etar, phir, z0r, d0r, a0r, nsctr, npixr, ntrtr, nstrawr };
534  for ( int it=0 ; it<10 ; it++ ) {
535  // std::string hname = name()+"_"+varName[it]+"_rec";
536  // std::string hname = varName[it]+"_rec";
537  // std::map<std::string, TH1F*>::iterator hmitr = m_histos.find(hname);
538  // if ( hmitr!=m_histos.end() ) hmitr->second->Fill( vpar[it] );
539  // else std::cerr << "hmmm histo " << hname << " not found" << std::endl;
540  if ( TH1F* hptr = find(varName[it]+"_rec") ) hptr->Fill( vpart[it] );
541  else std::cerr << "hmmm histo " << varName[it]+"_rec" << " not found" << std::endl;
542  }
543 #endif
544 
545 
546  // purities
547  if ( matchedref ) {
548 
549  // std::cout << *matchedref << std::endl;
550 
551  m_purity_pt->Fill(std::fabs(pTr));
552  m_purity_z0->Fill(z0r);
553  m_purity_eta->Fill(etar);
554  m_purity_phi->Fill(phir);
555  m_purity_d0->Fill(d0r);
556  m_purity_a0->Fill(a0r);
557 
558  // hnpix_v_sct_match->Fill( nsctr*0.5, npixr*0.5 );
559 
560  }
561  else {
562  // std::cout << "NULL" << std::endl;
563  m_purity_pt->FillDenom(std::fabs(pTr));
564  m_purity_z0->FillDenom(z0r);
565  m_purity_eta->FillDenom(etar);
566  m_purity_phi->FillDenom(phir);
567  m_purity_d0->FillDenom(d0r);
568  m_purity_a0->FillDenom(a0r);
569  }
570 
571  }
572 
573  if ( dump && m_print ) {
574 
575  std::cout << "PurityAnalysis::execute() missed a high pT track - dumping tracks" << std::endl;
576 
577  for ( int i=reftracks.size() ; i-- ; ) {
578 
579  if ( std::fabs( reftracks[i]->pT() ) > 1000 ) {
580  std::cout << "\t dump " << *reftracks[i];
581  const TIDA::Track* matchedreco = matcher->matched(reftracks[i]);
582  if ( matchedreco ) std::cout << " <--> " << *matchedreco << std::endl;
583  else std::cout << std::endl;
584  }
585 
586  }
587 
588  for ( int i=testtracks.size() ; i-- ; ) {
589  const TIDA::Track* matchedref = matcher->revmatched(testtracks[i]);
590  if ( matchedref==0 ) std::cout << "\t\t\t\t\t " << *testtracks[i] << std::endl;
591  }
592 
593  }
594 
595  if ( m_print ) std::cout << "PurityAnalysis::execute() exiting" << std::endl;
596 
597 }
598 
599 
600 
xAOD::iterator
JetConstituentVector::iterator iterator
Definition: JetConstituentVector.cxx:68
TIDA::Associator
Definition: TIDAAssociator.h:24
CalculateHighPtTerm.pT
pT
Definition: ICHEP2016/CalculateHighPtTerm.py:57
PurityAnalysis::m_purity_pt
Efficiency1D * m_purity_pt
Definition: PurityAnalysis.h:103
ev
int ev
Definition: globals.cxx:25
PurityAnalysis::m_eff_phi
Efficiency1D * m_eff_phi
Definition: PurityAnalysis.h:98
TIDDirectory::pop
void pop()
Definition: TIDDirectory.h:79
r
int r
Definition: globals.cxx:22
TIDA::Track::a0
double a0() const
Definition: Trigger/TrigAnalysis/TrigInDetAnalysis/TrigInDetAnalysis/Track.h:49
met::DeltaR
@ DeltaR
Definition: METRecoCommon.h:11
TrackFilter.h
base class for a single track selection filter allowing parameter setting for complex track selection
PurityAnalysis::execute
virtual void execute(const std::vector< TIDA::Track * > &reftracks, const std::vector< TIDA::Track * > &testtracks, TrackAssociator *matcher)
fill all the histograms - matched histograms, efficiencies etc
Definition: PurityAnalysis.cxx:282
ConvertOldUJHistosToNewHistos.etaBins
list etaBins
Definition: ConvertOldUJHistosToNewHistos.py:145
skel.it
it
Definition: skel.GENtoEVGEN.py:407
M_PI
#define M_PI
Definition: ActiveFraction.h:11
T_Efficiency::finalise
void finalise(double scale=100)
actually calculate the efficiencies
Definition: T_Efficiency.h:71
TIDA::Track::pT
double pT() const
Definition: Trigger/TrigAnalysis/TrigInDetAnalysis/TrigInDetAnalysis/Track.h:50
PurityAnalysis::addHistogram
void addHistogram(TH1F *h)
Definition: PurityAnalysis.h:77
PurityAnalysis::m_eff_pt
Efficiency1D * m_eff_pt
Definition: PurityAnalysis.h:93
drawFromPickle.cos
cos
Definition: drawFromPickle.py:36
TIDA::Track::z0
double z0() const
Definition: Trigger/TrigAnalysis/TrigInDetAnalysis/TrigInDetAnalysis/Track.h:48
TIDA::Associator::matched
virtual const S * matched(T *t)
Definition: TIDAAssociator.h:45
PurityAnalysis::m_eff_d0
Efficiency1D * m_eff_d0
Definition: PurityAnalysis.h:100
run_Egamma1_LArStrip_Fex.dump
dump
Definition: run_Egamma1_LArStrip_Fex.py:87
PurityAnalysis::m_purity_phi
Efficiency1D * m_purity_phi
Definition: PurityAnalysis.h:105
TIDDirectory
Definition: TIDDirectory.h:25
TrackAnalysis::name
const std::string & name() const
return identifier
Definition: TrackAnalysis.h:52
PurityAnalysis::finalise
virtual void finalise()
calculate the efficiencies and write them out with all the histograms
Definition: PurityAnalysis.cxx:227
Efficiency1D
Definition: Efficiency1D.h:19
ParseInputs.gDirectory
gDirectory
Definition: Final2012/ParseInputs.py:133
lumiFormat.i
int i
Definition: lumiFormat.py:85
ts
int ts
Definition: globals.cxx:24
PurityAnalysis::m_purity_eta
Efficiency1D * m_purity_eta
Definition: PurityAnalysis.h:104
PixelAthClusterMonAlgCfg.varName
string varName
end cluster ToT and charge
Definition: PixelAthClusterMonAlgCfg.py:125
Efficiency1D::FillDenom
virtual void FillDenom(double x, float w=1)
Definition: Efficiency1D.h:42
lb
int lb
Definition: globals.cxx:23
PurityAnalysis::m_Nref
int m_Nref
Definition: PurityAnalysis.h:117
add-xsec-uncert-quadrature-N.label
label
Definition: add-xsec-uncert-quadrature-N.py:104
TIDA::Associator::revmatched
virtual const T * revmatched(S *t)
Definition: TIDAAssociator.h:52
a0
double a0
Definition: globals.cxx:27
PurityAnalysis::initialise
virtual void initialise()
book all the histograms
Definition: PurityAnalysis.cxx:34
TIDA::Track::phi
double phi() const
Definition: Trigger/TrigAnalysis/TrigInDetAnalysis/TrigInDetAnalysis/Track.h:47
Efficiency1D::Fill
virtual void Fill(double x, double w=1)
Definition: Efficiency1D.h:37
PurityAnalysis::find
TH1F * find(const std::string &n)
Definition: PurityAnalysis.h:82
MakeTH3DFromTH2Ds.zBins
list zBins
Definition: MakeTH3DFromTH2Ds.py:86
TrackAnalysis::m_yBeamReference
double m_yBeamReference
Definition: TrackAnalysis.h:154
PurityAnalysis::m_hchi2
TH1F * m_hchi2
Definition: PurityAnalysis.h:120
beamspotman.dir
string dir
Definition: beamspotman.py:619
PurityAnalysis::m_print
bool m_print
flag to print out the matched tracks etc
Definition: PurityAnalysis.h:123
PurityAnalysis.h
PurityAnalysis::m_dir
TIDDirectory * m_dir
Definition: PurityAnalysis.h:89
TIDA::Track::eta
double eta() const
Definition: Trigger/TrigAnalysis/TrigInDetAnalysis/TrigInDetAnalysis/Track.h:46
TIDDirectory::push
void push()
Definition: TIDDirectory.h:78
PurityAnalysis::m_eff_ptm
Efficiency1D * m_eff_ptm
Definition: PurityAnalysis.h:95
PurityAnalysis::m_histos
std::map< std::string, TH1F * > m_histos
Definition: PurityAnalysis.h:91
PurityAnalysis::m_eff_eta
Efficiency1D * m_eff_eta
Definition: PurityAnalysis.h:97
TrackAnalysis::m_xBeamTest
double m_xBeamTest
test sample
Definition: TrackAnalysis.h:158
PurityAnalysis::m_hDeltaR
TH1F * m_hDeltaR
Definition: PurityAnalysis.h:113
PurityAnalysis::m_eff_z0
Efficiency1D * m_eff_z0
Definition: PurityAnalysis.h:99
TIDA::Track
Definition: Trigger/TrigAnalysis/TrigInDetAnalysis/TrigInDetAnalysis/Track.h:26
PurityAnalysis::m_purity_a0
Efficiency1D * m_purity_a0
Definition: PurityAnalysis.h:108
TrackAnalysis::m_yBeamTest
double m_yBeamTest
Definition: TrackAnalysis.h:159
python.TrigEgammaMonitorHelper.TH1F
def TH1F(name, title, nxbins, bins_par2, bins_par3=None, path='', **kwargs)
Definition: TrigEgammaMonitorHelper.py:24
PurityAnalysis::m_Nmatched
int m_Nmatched
Definition: PurityAnalysis.h:118
PurityAnalysis::m_purity_z0
Efficiency1D * m_purity_z0
Definition: PurityAnalysis.h:106
drawFromPickle.sin
sin
Definition: drawFromPickle.py:36
PurityAnalysis::m_eff_ptp
Efficiency1D * m_eff_ptp
Definition: PurityAnalysis.h:94
PurityAnalysis::m_purity_d0
Efficiency1D * m_purity_d0
Definition: PurityAnalysis.h:107
pow
constexpr int pow(int base, int exp) noexcept
Definition: ap_fixedTest.cxx:15
PurityAnalysis::m_Nreco
int m_Nreco
number of reconstructed tracks
Definition: PurityAnalysis.h:116
hipt
bool hipt
Definition: globals.cxx:29
PurityAnalysis::m_eff_a0
Efficiency1D * m_eff_a0
Definition: PurityAnalysis.h:101
TrackAnalysis::m_xBeamReference
double m_xBeamReference
beamline positions reference sample
Definition: TrackAnalysis.h:153
Generated on Thu Oct 30 2025 21:16:37 for ATLAS Offline Software by   doxygen 1.8.18