Belle II Software  release-05-02-19
PerformanceEvaluationBaseClass.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2010 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Giulia Casarosa *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include <tracking/modules/trackingPerformanceEvaluation/PerformanceEvaluationBaseClass.h>
12 
13 #include <framework/logging/Logger.h>
14 
15 #include <root/TAxis.h>
16 
17 using namespace Belle2;
18 
19 PerformanceEvaluationBaseClass:: PerformanceEvaluationBaseClass()
20 {
21 }
22 
23 PerformanceEvaluationBaseClass::~ PerformanceEvaluationBaseClass()
24 {
25 
26 }
27 
28 TH1F* PerformanceEvaluationBaseClass::createHistogram1D(const char* name, const char* title,
29  Int_t nbins, Double_t min, Double_t max,
30  const char* xtitle, TList* histoList)
31 {
32 
33  TH1F* h = new TH1F(name, title, nbins, min, max);
34 
35  h->GetXaxis()->SetTitle(xtitle);
36 
37  if (histoList)
38  histoList->Add(h);
39 
40  return h;
41 }
42 
43 TH1F* PerformanceEvaluationBaseClass::createHistogram1D(const char* name, const char* title,
44  Int_t nbins, Double_t* bins,
45  const char* xtitle, TList* histoList)
46 {
47 
48  TH1F* h = new TH1F(name, title, nbins, bins);
49 
50  h->GetXaxis()->SetTitle(xtitle);
51 
52  if (histoList)
53  histoList->Add(h);
54 
55  return h;
56 }
57 
58 TH2F* PerformanceEvaluationBaseClass::createHistogram2D(const char* name, const char* title,
59  Int_t nbinsX, Double_t minX, Double_t maxX,
60  const char* titleX,
61  Int_t nbinsY, Double_t minY, Double_t maxY,
62  const char* titleY, TList* histoList)
63 {
64 
65  TH2F* h = new TH2F(name, title, nbinsX, minX, maxX, nbinsY, minY, maxY);
66 
67  h->GetXaxis()->SetTitle(titleX);
68  h->GetYaxis()->SetTitle(titleY);
69 
70  if (histoList)
71  histoList->Add(h);
72 
73  return h;
74 }
75 
76 TH2F* PerformanceEvaluationBaseClass::createHistogram2D(const char* name, const char* title,
77  Int_t nbinsX, Double_t* binsX,
78  const char* titleX,
79  Int_t nbinsY, Double_t* binsY,
80  const char* titleY,
81  TList* histoList)
82 {
83 
84  TH2F* h = new TH2F(name, title, nbinsX, binsX, nbinsY, binsY);
85 
86  h->GetXaxis()->SetTitle(titleX);
87  h->GetYaxis()->SetTitle(titleY);
88 
89  if (histoList)
90  histoList->Add(h);
91 
92  return h;
93 }
94 
95 TH3F* PerformanceEvaluationBaseClass::createHistogram3D(const char* name, const char* title,
96  Int_t nbinsX, Double_t minX, Double_t maxX,
97  const char* titleX,
98  Int_t nbinsY, Double_t minY, Double_t maxY,
99  const char* titleY,
100  Int_t nbinsZ, Double_t minZ, Double_t maxZ,
101  const char* titleZ,
102  TList* histoList)
103 {
104 
105  TH3F* h = new TH3F(name, title, nbinsX, minX, maxX, nbinsY, minY, maxY, nbinsZ, minZ, maxZ);
106 
107  h->GetXaxis()->SetTitle(titleX);
108  h->GetYaxis()->SetTitle(titleY);
109  h->GetZaxis()->SetTitle(titleZ);
110 
111  if (histoList)
112  histoList->Add(h);
113 
114  return h;
115 }
116 
117 TH3F* PerformanceEvaluationBaseClass::createHistogram3D(const char* name, const char* title,
118  Int_t nbinsX, Double_t* binsX,
119  const char* titleX,
120  Int_t nbinsY, Double_t* binsY,
121  const char* titleY,
122  Int_t nbinsZ, Double_t* binsZ,
123  const char* titleZ,
124  TList* histoList)
125 {
126 
127  TH3F* h = new TH3F(name, title, nbinsX, binsX, nbinsY, binsY, nbinsZ, binsZ);
128 
129  h->GetXaxis()->SetTitle(titleX);
130  h->GetYaxis()->SetTitle(titleY);
131  h->GetZaxis()->SetTitle(titleZ);
132 
133  if (histoList)
134  histoList->Add(h);
135 
136  return h;
137 }
138 
139 TH1* PerformanceEvaluationBaseClass::duplicateHistogram(const char* newname, const char* newtitle,
140  TH1* h, TList* histoList)
141 {
142 
143  TH1F* h1 = dynamic_cast<TH1F*>(h);
144  TH2F* h2 = dynamic_cast<TH2F*>(h);
145  TH3F* h3 = dynamic_cast<TH3F*>(h);
146 
147  TH1* newh = 0;
148 
149  if (h1)
150  newh = new TH1F(*h1);
151  if (h2)
152  newh = new TH2F(*h2);
153  if (h3)
154  newh = new TH3F(*h3);
155 
156  newh->SetName(newname);
157  newh->SetTitle(newtitle);
158 
159  if (histoList)
160  histoList->Add(newh);
161 
162 
163  return newh;
164 }
165 
166 TH1F* PerformanceEvaluationBaseClass::createHistogramsRatio(const char* name, const char* title,
167  TH1* hNum, TH1* hDen, bool isEffPlot,
168  int axisRef)
169 {
170 
171  TH1F* h1den = dynamic_cast<TH1F*>(hDen);
172  TH1F* h1num = dynamic_cast<TH1F*>(hNum);
173  TH2F* h2den = dynamic_cast<TH2F*>(hDen);
174  TH2F* h2num = dynamic_cast<TH2F*>(hNum);
175  TH3F* h3den = dynamic_cast<TH3F*>(hDen);
176  TH3F* h3num = dynamic_cast<TH3F*>(hNum);
177 
178  TH1* hden = 0;
179  TH1* hnum = 0;
180 
181  if (h1den) {
182  hden = new TH1F(*h1den);
183  hnum = new TH1F(*h1num);
184  }
185  if (h2den) {
186  hden = new TH2F(*h2den);
187  hnum = new TH2F(*h2num);
188  }
189  if (h3den) {
190  hden = new TH3F(*h3den);
191  hnum = new TH3F(*h3num);
192  }
193 
194  TAxis* the_axis;
195  TAxis* the_other1;
196  TAxis* the_other2;
197 
198  if (axisRef == 0) {
199  the_axis = hden->GetXaxis();
200  the_other1 = hden->GetYaxis();
201  the_other2 = hden->GetZaxis();
202  } else if (axisRef == 1) {
203  the_axis = hden->GetYaxis();
204  the_other1 = hden->GetXaxis();
205  the_other2 = hden->GetZaxis();
206  } else if (axisRef == 2) {
207  the_axis = hden->GetZaxis();
208  the_other1 = hden->GetXaxis();
209  the_other2 = hden->GetYaxis();
210  } else
211  return NULL;
212 
213 
214  TH1F* h;
215  if (the_axis->GetXbins()->GetSize())
216  h = new TH1F(name, title, the_axis->GetNbins(), (the_axis->GetXbins())->GetArray());
217  else
218  h = new TH1F(name, title, the_axis->GetNbins(), the_axis->GetXmin(), the_axis->GetXmax());
219  h->GetXaxis()->SetTitle(the_axis->GetTitle());
220 
221  h->GetYaxis()->SetRangeUser(0.00001, 1);
222 
223  Int_t bin = 0;
224  Int_t nBins = 0;
225 
226  for (int the_bin = 1; the_bin < the_axis->GetNbins() + 1; the_bin++) {
227 
228  double num = 0;
229  double den = 0;
230 
231  for (int other1_bin = 1; other1_bin < the_other1->GetNbins() + 1; other1_bin++)
232  for (int other2_bin = 1; other2_bin < the_other2->GetNbins() + 1; other2_bin++) {
233 
234  if (axisRef == 0) bin = hden->GetBin(the_bin, other1_bin, other2_bin);
235  else if (axisRef == 1) bin = hden->GetBin(other1_bin, the_bin, other2_bin);
236  else if (axisRef == 2) bin = hden->GetBin(other1_bin, other2_bin, the_bin);
237 
238  if (hden->IsBinUnderflow(bin))
239  B2INFO(" bin = " << bin << "(" << the_bin << "," << other1_bin << "," << other2_bin << "), UNDERFLOW");
240  if (hden->IsBinOverflow(bin))
241  B2INFO(" bin = " << bin << "(" << the_bin << "," << other1_bin << "," << other2_bin << "), OVERFLOW");
242 
243  num += hnum->GetBinContent(bin);
244  den += hden->GetBinContent(bin);
245 
246  nBins++;
247  }
248 
249  double eff = 0;
250  double err = 0;
251 
252  if (den > 0) {
253  eff = (double)num / den;
254  err = sqrt(eff * (1 - eff)) / sqrt(den);
255  }
256 
257  if (isEffPlot) {
258  h->SetBinContent(the_bin, eff);
259  h->SetBinError(the_bin, err);
260  } else {
261  h->SetBinContent(the_bin, 1 - eff);
262  h->SetBinError(the_bin, err);
263  }
264  }
265 
266  return h;
267 
268 }
269 
270 
271 void PerformanceEvaluationBaseClass::addInefficiencyPlots(TList* histoList, TH3F* h3_xPerMCParticle, TH3F* h3_MCParticle)
272 {
273 
274  if ((h3_xPerMCParticle == NULL) || (h3_MCParticle == NULL))
275  return;
276 
277  //normalized to MCParticles
278  TH1F* h_ineff_pt = createHistogramsRatio("hineffpt", "inefficiency VS pt, normalized to MCParticles", h3_xPerMCParticle,
279  h3_MCParticle, false, 0);
280  histoList->Add(h_ineff_pt);
281 
282  TH1F* h_ineff_theta = createHistogramsRatio("hinefftheta", "inefficiency VS #theta, normalized to MCParticles",
283  h3_xPerMCParticle, h3_MCParticle, false, 1);
284  histoList->Add(h_ineff_theta);
285 
286  TH1F* h_ineff_phi = createHistogramsRatio("hineffphi", "inefficiency VS #phi, normalized to MCParticles", h3_xPerMCParticle,
287  h3_MCParticle, false, 2);
288  histoList->Add(h_ineff_phi);
289 
290 }
291 
292 void PerformanceEvaluationBaseClass::addEfficiencyPlots(TList* histoList, TH3F* h3_xPerMCParticle, TH3F* h3_MCParticle)
293 {
294  if ((h3_xPerMCParticle == NULL) || (h3_MCParticle == NULL))
295  return;
296 
297  //normalized to MCParticles
298  TH1F* h_eff_pt = createHistogramsRatio("heffpt", "efficiency VS pt, normalized to MCParticles", h3_xPerMCParticle,
299  h3_MCParticle, true, 0);
300  histoList->Add(h_eff_pt);
301 
302  TH1F* h_eff_theta = createHistogramsRatio("hefftheta", "efficiency VS #theta, normalized to MCParticles", h3_xPerMCParticle,
303  h3_MCParticle, true, 1);
304  histoList->Add(h_eff_theta);
305 
306  TH1F* h_eff_phi = createHistogramsRatio("heffphi", "efficiency VS #phi, normalized to MCParticles", h3_xPerMCParticle,
307  h3_MCParticle, true, 2);
308  histoList->Add(h_eff_phi);
309 
310 }
311 
312 
313 
314 void PerformanceEvaluationBaseClass::addPurityPlots(TList* histoList, TH3F* h3_MCParticlesPerX, TH3F* h3_X)
315 {
316  if ((h3_X == NULL) || (h3_MCParticlesPerX == NULL))
317  return;
318 
319  //purity histograms
320  TH1F* h_pur_pt = createHistogramsRatio("hpurpt", "purity VS pt", h3_MCParticlesPerX, h3_X, true, 0);
321  histoList->Add(h_pur_pt);
322 
323  TH1F* h_pur_theta = createHistogramsRatio("hpurtheta", "purity VS #theta", h3_MCParticlesPerX, h3_X, true, 1);
324  histoList->Add(h_pur_theta);
325 
326  TH1F* h_pur_phi = createHistogramsRatio("hpurphi", "purity VS #phi", h3_MCParticlesPerX, h3_X, true, 2);
327  histoList->Add(h_pur_phi);
328 
329 }
330 
331 TH1F* PerformanceEvaluationBaseClass::effPlot1D(TH1F* h1_den, TH1F* h1_num, const char* name, const char* title,
332  bool geo_accettance, TList* histoList)
333 {
334 
335  std::string total;
336  std::string trueTitle;
337 
338  if (geo_accettance == false) {
339  std::string name1 = "_noGeoAcc";
340  total = std::string(name) + name1;
341  trueTitle = std::string(title) + name1;
342  } else {
343  std::string name2 = "_withGeoAcc";
344  total = std::string(name) + name2;
345  trueTitle = std::string(title) + name2;
346  }
347 
348  TH1F* h = (TH1F*)duplicateHistogram(total.c_str(), trueTitle.c_str(), h1_den, histoList);
349  h->GetYaxis()->SetRangeUser(0., 1);
350 
351  for (int bin = 0; bin < h->GetXaxis()->GetNbins(); bin++) {
352  float num = h1_num->GetBinContent(bin + 1);
353  float den = h1_den->GetBinContent(bin + 1);
354  double eff = 0.;
355  double err = 0.;
356 
357  if (den > 0) {
358  eff = (double)num / den;
359  err = sqrt(eff * (1 - eff)) / sqrt(den);
360  }
361  h->SetBinContent(bin + 1, eff);
362  h->SetBinError(bin + 1, err);
363  }
364 
365 
366  if (histoList)
367  histoList->Add(h);
368 
369  return h;
370 }
371 
372 TH1F* PerformanceEvaluationBaseClass::effPlot1D(TH1F* h1_MC, TH1F* h1_RecoTrack, TH1F* h1_Track, const char* name,
373  const char* title, TList* histoList)
374 {
375  if (h1_Track == NULL) B2INFO("h_Track missing");
376 
377  std::string name1 = "_noGeoAcc";
378  std::string name2 = "_withGeoAcc";
379 
380  std::string total1 = std::string(name) + name1;
381  std::string total2 = std::string(name) + name2;
382 
383  std::string title1 = std::string(title) + name1;
384  std::string title2 = std::string(title) + name2;
385 
386  TH1F* h[2];
387 
388  h[0] = (TH1F*)duplicateHistogram(total2.c_str(), title2.c_str(), h1_RecoTrack, histoList);
389  h[0]->GetYaxis()->SetRangeUser(0., 1);
390 
391  for (int bin = 0; bin < h[0]->GetXaxis()->GetNbins(); bin++) {
392  float num = h1_Track->GetBinContent(bin + 1);
393  float den = h1_RecoTrack->GetBinContent(bin + 1);
394  double eff = 0.;
395  double err = 0.;
396 
397  if (den > 0) {
398  eff = (double)num / den;
399  err = sqrt(eff * (1 - eff)) / sqrt(den);
400  }
401  h[0]->SetBinContent(bin + 1, eff);
402  h[0]->SetBinError(bin + 1, err);
403  }
404 
405  h[1] = (TH1F*)duplicateHistogram(total1.c_str(), title1.c_str(), h1_MC, histoList);
406  h[1]->GetYaxis()->SetRangeUser(0., 1);
407 
408  for (int bin = 0; bin < h[1]->GetXaxis()->GetNbins(); bin++) {
409  float num = h1_Track->GetBinContent(bin + 1);
410  float den = h1_MC->GetBinContent(bin + 1);
411  double eff = 0.;
412  double err = 0.;
413 
414  if (den > 0) {
415  eff = (double)num / den;
416  err = sqrt(eff * (1 - eff)) / sqrt(den);
417  }
418  h[1]->SetBinContent(bin + 1, eff);
419  h[1]->SetBinError(bin + 1, err);
420  }
421 
422  if (histoList) {
423  histoList->Add(h[0]);
424  histoList->Add(h[1]);
425  }
426 
427  return *h;
428 }
429 
430 TH2F* PerformanceEvaluationBaseClass::effPlot2D(TH2F* h2_den, TH2F* h2_num,
431  const char* name, const char* title, bool geo_accettance, TList* histoList)
432 {
433  std::string err_char = "_error_";
434  std::string addTitle = "Errors, ";
435 
436  std::string total;
437  std::string error;
438  std::string trueTitle;
439 
440  std::string titleErr = addTitle + std::string(title);
441 
442  if (geo_accettance == false) {
443  std::string name1 = "_noGeoAcc";
444  total = std::string(name) + name1;
445  trueTitle = std::string(title) + name1;
446  error = std::string(name) + err_char + std::string(name1);
447  } else {
448  std::string name2 = "_withGeoAcc";
449  total = std::string(name) + name2;
450  trueTitle = std::string(title) + name2;
451  error = std::string(name) + err_char + std::string(name2);
452  }
453 
454  TH2F* h2[2];
455  h2[0] = (TH2F*)duplicateHistogram(total.c_str(), trueTitle.c_str(), h2_den, histoList);
456  h2[1] = (TH2F*)duplicateHistogram(error.c_str(), titleErr.c_str(), h2_den, histoList);
457 
458  for (int binX = 0; binX < h2[0]->GetXaxis()->GetNbins(); binX++) {
459  for (int binY = 0; binY < h2[0]->GetYaxis()->GetNbins(); binY++) {
460  float num = h2_num->GetBinContent(binX + 1, binY + 1);
461  float den = h2_den->GetBinContent(binX + 1, binY + 1);
462  double eff = 0.;
463  double err = 0.;
464 
465  if (den > 0) {
466  eff = (double)num / den;
467  err = sqrt(eff * (1 - eff)) / sqrt(den);
468  }
469 
470  h2[0]->SetBinContent(binX + 1, binY + 1, eff);
471  h2[0]->SetBinError(binX + 1, binY + 1, err);
472  h2[1]->SetBinContent(binX + 1, binY + 1, err);
473  }
474  }
475 
476 
477  if (histoList) {
478  histoList->Add(h2[0]);
479  histoList->Add(h2[1]);
480  }
481 
482  return *h2;
483 
484 }
485 
486 TH2F* PerformanceEvaluationBaseClass::effPlot2D(TH2F* h2_MC, TH2F* h2_RecoTrack, TH2F* h2_Track, const char* name,
487  const char* title, TList* histoList)
488 {
489  if (h2_Track == NULL) B2INFO("h_Track missing");
490 
491  std::string name1 = "_noGeoAcc";
492  std::string name2 = "_withGeoAcc";
493  std::string err_char = "_error_";
494  std::string addTitle = "Errors, ";
495 
496  std::string total1 = std::string(name) + name1;
497  std::string total2 = std::string(name) + name2;
498 
499  std::string title1 = std::string(title) + name1;
500  std::string title2 = std::string(title) + name2;
501 
502  std::string error1 = std::string(name) + err_char + name1;
503  std::string error2 = std::string(name) + err_char + name2;
504  std::string titleErr = addTitle + std::string(title);
505 
506  TH2F* h2[4];
507 
508  h2[0] = (TH2F*)duplicateHistogram(total2.c_str(), title2.c_str(), h2_RecoTrack, histoList);
509  h2[1] = (TH2F*)duplicateHistogram(error2.c_str(), titleErr.c_str(), h2_RecoTrack, histoList);
510 
511  for (int binX = 0; binX < h2[0]->GetXaxis()->GetNbins(); binX++) {
512  for (int binY = 0; binY < h2[0]->GetYaxis()->GetNbins(); binY++) {
513  float num = h2_Track->GetBinContent(binX + 1, binY + 1);
514  float den = h2_RecoTrack->GetBinContent(binX + 1, binY + 1);
515  double eff = 0.;
516  double err = 0.;
517 
518  if (den > 0) {
519  eff = num / den;
520  err = sqrt(eff * (1 - eff)) / sqrt(den);
521  }
522 
523  h2[0]->SetBinContent(binX + 1, binY + 1, eff);
524  h2[0]->SetBinError(binX + 1, binY + 1, err);
525  h2[1]->SetBinContent(binX + 1, binY + 1, err);
526  }
527  }
528 
529  h2[2] = (TH2F*)duplicateHistogram(total1.c_str(), title1.c_str(), h2_MC, histoList);
530  h2[3] = (TH2F*)duplicateHistogram(error1.c_str(), titleErr.c_str(), h2_MC, histoList);
531 
532  for (int binX = 0; binX < h2[2]->GetXaxis()->GetNbins(); binX++) {
533  for (int binY = 0; binY < h2[2]->GetYaxis()->GetNbins(); binY++) {
534  float num = h2_Track->GetBinContent(binX + 1, binY + 1);
535  float den = h2_MC->GetBinContent(binX + 1, binY + 1);
536  double eff = 0.;
537  double err = 0.;
538 
539  if (den > 0) {
540  eff = num / den;
541  err = sqrt(eff * (1 - eff)) / sqrt(den);
542  }
543 
544  h2[2]->SetBinContent(binX + 1, binY + 1, eff);
545  h2[2]->SetBinError(binX + 1, binY + 1, err);
546  h2[3]->SetBinContent(binX + 1, binY + 1, err);
547  }
548  }
549 
550  if (histoList) {
551  histoList->Add(h2[0]);
552  histoList->Add(h2[1]);
553  histoList->Add(h2[2]);
554  histoList->Add(h2[3]);
555  }
556 
557  return *h2;
558 }
559 
560 TH1F* PerformanceEvaluationBaseClass::geoAcc1D(TH1F* h1_den, TH1F* h1_num, const char* name, const char* title, TList* histoList)
561 {
562  TH1F* h = (TH1F*)duplicateHistogram(name, title, h1_den, histoList);
563  h->GetYaxis()->SetRangeUser(0., 1);
564 
565  for (int bin = 0; bin < h->GetXaxis()->GetNbins(); bin++) {
566  float num = h1_num->GetBinContent(bin + 1);
567  float den = h1_den->GetBinContent(bin + 1);
568  double eff = 0.;
569  double err = 0.;
570 
571  if (den > 0) {
572  eff = (double)num / den;
573  err = sqrt(eff * (1 - eff)) / sqrt(den);
574  }
575  h->SetBinContent(bin + 1, eff);
576  h->SetBinError(bin + 1, err);
577  }
578 
579  if (histoList)
580  histoList->Add(h);
581 
582  return h;
583 }
584 
585 TH2F* PerformanceEvaluationBaseClass::geoAcc2D(TH2F* h2_den, TH2F* h2_num,
586  const char* name, const char* title, TList* histoList)
587 {
588  std::string err_char = "_err";
589  std::string addTitle = "Errors, ";
590 
591  std::string error = std::string(name) + err_char;
592  std::string titleErr = addTitle + std::string(title);
593 
594  TH2F* h2[2];
595  h2[0] = (TH2F*)duplicateHistogram(name, title, h2_den, histoList);
596  h2[1] = (TH2F*)duplicateHistogram(error.c_str(), titleErr.c_str(), h2_den, histoList);
597 
598  for (int binX = 0; binX < h2[0]->GetXaxis()->GetNbins(); binX++) {
599  for (int binY = 0; binY < h2[0]->GetYaxis()->GetNbins(); binY++) {
600  float num = h2_num->GetBinContent(binX + 1, binY + 1);
601  float den = h2_den->GetBinContent(binX + 1, binY + 1);
602  double eff = 0.;
603  double err = 0.;
604 
605  if (den > 0) {
606  eff = (double)num / den;
607  err = sqrt(eff * (1 - eff)) / sqrt(den);
608  }
609  h2[0]->SetBinContent(binX + 1, binY + 1, eff);
610  h2[0]->SetBinError(binX + 1, binY + 1, err);
611  h2[1]->SetBinContent(binX + 1, binY + 1, err);
612  }
613  }
614 
615  if (histoList) {
616  histoList->Add(h2[0]);
617  histoList->Add(h2[1]);
618  }
619 
620  return *h2;
621 
622 }
623 
624 TH1F* PerformanceEvaluationBaseClass::V0FinderEff(TH1F* h1_dau0, TH1F* h1_dau1, TH1F* h1_Mother,
625  const char* name, const char* title, TList* histoList)
626 {
627 
628  TH1F* h = (TH1F*)duplicateHistogram(name, title, h1_Mother, histoList);
629  h->GetYaxis()->SetRangeUser(0., 1);
630 
631  for (int bin = 0; bin < h->GetXaxis()->GetNbins(); bin++) {
632  double dau0 = h1_dau0->GetBinContent(bin);
633  double dau1 = h1_dau1->GetBinContent(bin);
634  double Mother = h1_Mother->GetBinContent(bin);
635  double dau0Err = h1_dau0->GetBinError(bin);
636  double dau1Err = h1_dau1->GetBinError(bin);
637  double MotherErr = h1_Mother->GetBinError(bin);
638 
639  double binCont = 1. * Mother / dau0 / dau1;
640  double binErr = binCont * sqrt((dau0Err / dau0) * (dau0Err / dau0) + (dau1Err / dau1) * (dau1Err / dau1) * (MotherErr / Mother) *
641  (MotherErr / Mother));
642 
643  h->SetBinContent(bin, binCont);
644  h->SetBinError(bin, binErr);
645  }
646 
647  if (histoList)
648  histoList->Add(h);
649 
650  return h;
651 }
652 
Belle2::PerformanceEvaluationBaseClass::effPlot1D
TH1F * effPlot1D(TH1F *h1_den, TH1F *h1_num, const char *name, const char *title, bool geo_accettance, TList *histoList=NULL)
Create a 1D efficiency histogram and add it to the TList of 1D-histograms.
Definition: PerformanceEvaluationBaseClass.cc:331
Belle2::PerformanceEvaluationBaseClass::geoAcc1D
TH1F * geoAcc1D(TH1F *h1_den, TH1F *h1_num, const char *name, const char *title, TList *histoList=NULL)
Create a 1D efficiency histogram for geometric acceptance and add it to the TList of 1D-histograms.
Definition: PerformanceEvaluationBaseClass.cc:560
Belle2::PerformanceEvaluationBaseClass::createHistogram2D
TH2F * createHistogram2D(const char *name, const char *title, Int_t nbinsX, Double_t minX, Double_t maxX, const char *titleX, Int_t nbinsY, Double_t minY, Double_t maxY, const char *titleY, TList *histoList=NULL)
Create a 2D histogram and add it to the TList of 2D-histograms.
Definition: PerformanceEvaluationBaseClass.cc:58
Belle2::PerformanceEvaluationBaseClass::addInefficiencyPlots
void addInefficiencyPlots(TList *graphList=NULL, TH3F *h3_xPerMCParticle=NULL, TH3F *h3_MCParticle=NULL)
Create pt-, theta- and phi-inefficiency 1D histograms and add them to the TList of 1D-histograms.
Definition: PerformanceEvaluationBaseClass.cc:271
Belle2::PerformanceEvaluationBaseClass::createHistogramsRatio
TH1F * createHistogramsRatio(const char *name, const char *title, TH1 *hNum, TH1 *hDen, bool isEffPlot, int axisRef)
Make a new 1D histogram from the ratio of two others and add it to the TList of 1D-histograms.
Definition: PerformanceEvaluationBaseClass.cc:166
Belle2::PerformanceEvaluationBaseClass::effPlot2D
TH2F * effPlot2D(TH2F *h2_den, TH2F *h2_num, const char *name, const char *title, bool geo_accettance, TList *histoList=NULL)
Create a 2D efficiency histogram and add it to the TList of 2D-histograms.
Definition: PerformanceEvaluationBaseClass.cc:430
Belle2::PerformanceEvaluationBaseClass::V0FinderEff
TH1F * V0FinderEff(TH1F *h1_dau0, TH1F *h1_dau1, TH1F *h1_Mother, const char *name, const char *title, TList *histoList=NULL)
Create a 1D efficiency histogram for V0 finding and add it to the TList of 1D-histograms.
Definition: PerformanceEvaluationBaseClass.cc:624
Belle2::PerformanceEvaluationBaseClass::addEfficiencyPlots
void addEfficiencyPlots(TList *graphList=NULL, TH3F *h3_xPerMCParticle=NULL, TH3F *h3_MCParticle=NULL)
Create pt-, theta- and phi-efficiency 1D histograms and add them to the TList of 1D-histograms.
Definition: PerformanceEvaluationBaseClass.cc:292
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::PerformanceEvaluationBaseClass::addPurityPlots
void addPurityPlots(TList *graphList=NULL, TH3F *h3_xPerMCParticle=NULL, TH3F *h3_MCParticle=NULL)
Create pt-, theta- and phi-purity 1D histograms and add them to the TList of 1D-histograms.
Definition: PerformanceEvaluationBaseClass.cc:314
Belle2::PerformanceEvaluationBaseClass::createHistogram3D
TH3F * createHistogram3D(const char *name, const char *title, Int_t nbinsX, Double_t minX, Double_t maxX, const char *titleX, Int_t nbinsY, Double_t minY, Double_t maxY, const char *titleY, Int_t nbinsZ, Double_t minZ, Double_t maxZ, const char *titleZ, TList *histoList=NULL)
Create a 3D histogram and add it to the TList of 3D-histograms.
Definition: PerformanceEvaluationBaseClass.cc:95
Belle2::PerformanceEvaluationBaseClass::createHistogram1D
TH1F * createHistogram1D(const char *name, const char *title, Int_t nbins, Double_t min, Double_t max, const char *xtitle, TList *histoList=NULL)
Create a 1D histogram and add it to the TList of 1D-histograms.
Definition: PerformanceEvaluationBaseClass.cc:28
Belle2::PerformanceEvaluationBaseClass::duplicateHistogram
TH1 * duplicateHistogram(const char *newname, const char *newtitle, TH1 *h, TList *histoList=NULL)
Make a copy of a 1D histogram and add it to the TList of 1D-histograms.
Definition: PerformanceEvaluationBaseClass.cc:139
Belle2::PerformanceEvaluationBaseClass::geoAcc2D
TH2F * geoAcc2D(TH2F *h2_den, TH2F *h2_num, const char *name, const char *title, TList *histoList=NULL)
Create a 2D efficiency histogram for geometric acceptance and add it to the TList of 2D-histograms.
Definition: PerformanceEvaluationBaseClass.cc:585