FairRoot/PandaRoot
Macros | Functions | Variables
cutfinder_fullev.C File Reference
#include <algorithm>
#include "TFile.h"
#include "TTree.h"
#include "TString.h"
#include "TH1F.h"
#include "TCanvas.h"
#include "TROOT.h"
#include "TEventList.h"
#include "TDirectory.h"
#include <iostream>
#include "TLine.h"
#include "TLatex.h"
#include "TStyle.h"
#include "TObjArray.h"
#include "TPRegexp.h"
#include <map>
#include <utility>

Go to the source code of this file.

Macros

#define MAX   1000
 
#define BINS   500
 
#define NCAN   3
 

Functions

bool mycompare (int i, int j)
 
double bestSignificance (TH1 *h1, TH1 *h2, double &bestcut, double fac2=1.)
 
int countEvents (TTree *t, TString ccut)
 
int countEvents (TTree *t, TEventList &el)
 
double diffQA (TH1F &ht1, TH1F &ht2)
 
double bestEff (TH1 *h1, TH1 *h2, double &bestcut, double supr2=0.8, double fac2=1.)
 
double bestSuppression (TH1 *h1, TH1 *h2, double &bestcut, double eff=0.9)
 
double bestSuppressionEvt (TTree *t, TString precut, TString varname, TString sigcut, TString bgcut, double &bestcut, double eff, int idx)
 
double bestSeparation (TTree *t, TString varname, TEventList &els, TEventList &elb, double &bestcut, double eff, int idx)
 
double bestSuppressionEvt (TTree *t, TString varname, TEventList &els, TEventList &elb, double &bestcut, double eff, int idx)
 
void findLimits (TTree *t, TString var, TString ccut, double &low, double &high, double frac=0.98)
 
int cutfinder_fullev (TString pre, TString ntp, TString sigcut, TString precut="", double n0s=50., double n0b=500., double supr=0.9)
 

Variables

std::map< int, int > evcnt
 
double signi [MAX]
 
int idx [MAX]
 
TString vars [MAX]
 
double cut [MAX]
 
double minh [MAX]
 
double maxh [MAX]
 
double N0_sig
 
double N0_bg
 
double Nsigev
 
double Nbgev
 

Macro Definition Documentation

#define BINS   500

Definition at line 21 of file cutfinder_fullev.C.

#define MAX   1000

Definition at line 20 of file cutfinder_fullev.C.

Referenced by cutfinder_fullev().

#define NCAN   3

Definition at line 22 of file cutfinder_fullev.C.

Function Documentation

double bestEff ( TH1 *  h1,
TH1 *  h2,
double &  bestcut,
double  supr2 = 0.8,
double  fac2 = 1. 
)

Definition at line 147 of file cutfinder_fullev.C.

References i, and n.

148 {
149  int n=h1->GetNbinsX();
150  double sum1=0, sum2=0;
151  double bestr=0, bestl=0;
152  double int1 = h1->Integral();
153  double int2 = h2->Integral();
154  double bestcutl=-999.;
155  double bestcutr=999.;
156 
157  bestcut=-999.;
158  int i=1;
159 
160  while (sum2<supr2)
161  {
162  sum1 += h1->GetBinContent(i)/int1;
163  sum2 += h2->GetBinContent(i++)/int2;
164  }
165  bestr=(1.0-sum1);
166  bestcutr = h1->GetBinCenter(i-1);
167 
168  sum1=0;
169  sum2=0;
170  i=n;
171 
172  while (sum2<supr2)
173  {
174  sum1 += h1->GetBinContent(i)/int1;
175  sum2 += h2->GetBinContent(i--)/int2;
176  }
177 
178  bestl=(1.0-sum1);
179  bestcutl = h1->GetBinCenter(i-1);
180 
181  if (bestr>bestl)
182  {
183  bestcut = bestcutr;
184  return bestr;
185  }
186 
187  bestcut = bestcutl;
188  return bestl;
189 }
i
Int_t i
Definition: run_full.C:25
n
int n
Definition: HypStatDecay.cxx:21
h1
TH1F * h1
Definition: anal_hit_digi_cluster_fwendcap.C:121
h2
TH1F * h2
Definition: anal_hit_digi_cluster_fwendcap.C:122
double bestSeparation ( TTree *  t,
TString  varname,
TEventList &  els,
TEventList &  elb,
double &  bestcut,
double  eff,
int  idx 
)

Definition at line 326 of file cutfinder_fullev.C.

References ev, i, idx, max(), maxh, min(), minh, Nbgev, Nsigev, and run.

327 {
328  int i;
329  t->SetBranchStatus("*",0);
330  t->SetBranchStatus("ev",1);
331  t->SetBranchStatus("run",1);
332  t->SetBranchStatus(varname,1);
333 
334  Float_t ev, var, run;
335  t->SetBranchAddress("ev",&ev);
336  t->SetBranchAddress("run",&run);
337  t->SetBranchAddress(varname,&var);
338 
339  std::map<int, int> sigcnt;
340  std::map<int, int> bgcnt;
341  std::map<int, int> bgcnt2;
342 
343  std::vector<pair<double, int> > sigvals;
344  std::vector<pair<double, int> > bgvals;
345 
346  int Nsigval = els.GetN(), Nbgval = elb.GetN();
347  double min = 1e9, max=-1e9;
348 
349  // prepare the pairs of variable value, eventnumber for signal and count signals
350  for (i=0;i<els.GetN();++i)
351  {
352  t->GetEntry(els.GetEntry(i));
353  sigvals.push_back(std::make_pair(var, ev+run*2000));
354  if (var>max) max = var;
355  if (var<min) min = var;
356  }
357 
358  // prepare the pairs of variable value, eventnumber for background and count background
359  for (i=0;i<elb.GetN();++i)
360  {
361  t->GetEntry(elb.GetEntry(i));
362  bgvals.push_back(std::make_pair(var, ev+run*2000));
363  }
364 
365  // sort signals by first element = variable value
366  sort(sigvals.begin(), sigvals.end());
367 
368  i=0;
369  sigcnt.clear();
370  // find cut for current efficiency requirement left cut
371  while ( (sigcnt.size()/(double)Nsigev)<eff && i<Nsigval ) sigcnt[sigvals[i++].second]+=1;
372  double leftcut = sigvals[i].first;
373 
374  sigcnt.clear();
375  i=sigvals.size()-1;
376  // find cut for current efficiency requirement right cut
377  while ( (sigcnt.size()/(double)Nsigev)<eff && i>=0) sigcnt[sigvals[i--].second]+=1;
378  double rightcut = sigvals[i].first;
379 
380  bgcnt.clear();
381  bgcnt2.clear();
382 
383  for (i=0;i<bgvals.size();++i)
384  {
385  if (bgvals[i].first<=leftcut) bgcnt[bgvals[i].second];
386  if (bgvals[i].first>=rightcut) bgcnt2[bgvals[i].second];
387  if (bgvals[i].first>max) max = bgvals[i].first;
388  if (bgvals[i].first<min) min = bgvals[i].first;
389  }
390 
391  minh[idx]=min; maxh[idx]=max;
392 
393  t->SetBranchStatus("*",1);
394 
395 
396  double leftsupr = (Nbgev-bgcnt.size())/Nbgev;
397  double rightsupr = (Nbgev-bgcnt2.size())/Nbgev;
398 
399  cout <<varname<<" l:"<<leftcut<<"("<<leftsupr<<") r:"<<rightcut<<"("<<rightsupr<<")"<<endl;
400 
401  if (leftcut!=leftcut || rightcut!=rightcut) return 0;
402 
403  bestcut = rightcut;
404  if (leftsupr>rightsupr)
405  {
406  bestcut = leftcut;
407  return leftsupr;
408  }
409  return rightsupr;
410 
411 
412 }
els
TEventList * els
Definition: findcuts.C:52
run
Int_t run
Definition: autocutx.C:47
i
Int_t i
Definition: run_full.C:25
minh
double minh[MAX]
Definition: cutfinder_fullev.C:31
ev
int ev
Definition: GammaSpectraAnalysis.C:44
max
friend F32vec4 max(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:26
idx
int idx[MAX]
Definition: cutfinder_fullev.C:27
maxh
double maxh[MAX]
Definition: cutfinder_fullev.C:31
min
friend F32vec4 min(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:25
Nbgev
double Nbgev
Definition: cutfinder_fullev.C:34
elb
TEventList * elb
Definition: findcuts.C:52
Nsigev
double Nsigev
Definition: cutfinder_fullev.C:34
t
TTree * t
Definition: bump_analys.C:13
double bestSignificance ( TH1 *  h1,
TH1 *  h2,
double &  bestcut,
double  fac2 = 1. 
)

Definition at line 41 of file cutfinder_fullev.C.

References i, n, and sqrt().

42 {
43  int n=h1->GetNbinsX();
44  double sum1=0, sum2=0;
45  double bestr=0, bestl=0;
46  double int1 = h1->Integral();
47  double int2 = h2->Integral();
48  double bestcutl=-999.;
49  double bestcutr=999.;
50 
51  bestcut=-999.;
52 
53  for (int i=1;i<n;++i)
54  {
55  double cur1 = h1->GetBinContent(i)/int1;
56  double cur2 = h2->GetBinContent(i)/int2;
57  sum1 += cur1;
58  sum2 += cur2;
59 
60  double vall = 0, valr=0;
61 
62  if ((sum1+fac2*sum2)>0) vall=sum1/sqrt(fac2*sum2+sum1);
63  if ((1.-sum1+fac2*(1.-sum2))>0) valr=(1.-sum1)/sqrt(1.-sum1+fac2*(1.-sum2));
64 
65  if (valr>bestr) { bestr=valr; bestcutr = h1->GetBinCenter(i);}
66  if (vall>bestl) { bestl=vall; bestcutl = h1->GetBinCenter(i);}
67  }
68 
69  if (bestr>bestl)
70  {
71  bestcut = bestcutr;
72  return bestr;
73  }
74 
75  bestcut = bestcutl;
76  return bestl;
77 }
i
Int_t i
Definition: run_full.C:25
sqrt
friend F32vec4 sqrt(const F32vec4 &a)
Definition: P4_F32vec4.h:29
n
int n
Definition: HypStatDecay.cxx:21
h1
TH1F * h1
Definition: anal_hit_digi_cluster_fwendcap.C:121
h2
TH1F * h2
Definition: anal_hit_digi_cluster_fwendcap.C:122
double bestSuppression ( TH1 *  h1,
TH1 *  h2,
double &  bestcut,
double  eff = 0.9 
)

Definition at line 191 of file cutfinder_fullev.C.

References i, and n.

192 {
193  int n=h1->GetNbinsX();
194  double sum1=0, sum2=0;
195  double bestr=0, bestl=0;
196  double int1 = h1->Integral();
197  double int2 = h2->Integral();
198  double bestcutl=-999.;
199  double bestcutr=999.;
200 
201  bestcut=-999.;
202  int i=1;
203 
204  while (sum1<eff)
205  {
206  sum1 += h1->GetBinContent(i)/int1;
207  sum2 += h2->GetBinContent(i++)/int2;
208  }
209  bestr=(1.0-sum2);
210  bestcutr = h1->GetBinCenter(i-1);
211 
212  sum1=0;
213  sum2=0;
214  i=n;
215 
216  while (sum1<eff)
217  {
218  sum1 += h1->GetBinContent(i)/int1;
219  sum2 += h2->GetBinContent(i--)/int2;
220  }
221 
222  bestl=(1.0-sum2);
223  bestcutl = h1->GetBinCenter(i-1);
224 
225  if (bestr>bestl)
226  {
227  bestcut = bestcutr;
228  return bestr;
229  }
230 
231  bestcut = bestcutl;
232  return bestl;
233 }
i
Int_t i
Definition: run_full.C:25
n
int n
Definition: HypStatDecay.cxx:21
h1
TH1F * h1
Definition: anal_hit_digi_cluster_fwendcap.C:121
h2
TH1F * h2
Definition: anal_hit_digi_cluster_fwendcap.C:122
double bestSuppressionEvt ( TTree *  t,
TString  precut,
TString  varname,
TString  sigcut,
TString  bgcut,
double &  bestcut,
double  eff,
int  idx 
)

Definition at line 235 of file cutfinder_fullev.C.

References elb, els, ev, i, idx, max(), maxh, min(), minh, Nbgev, Nsigev, and run.

Referenced by cutfinder_fullev().

236 {
237  int i;
238  t->SetEventList(0);
239  t->SetBranchStatus("*",1);
240 
241  TEventList els("els"), elb("elb");
242  t->Draw(">>els",precut+"&&"+sigcut);
243  t->Draw(">>elb",precut+"&&"+bgcut);
244 
245  t->SetBranchStatus("*",0);
246  t->SetBranchStatus("ev",1);
247  t->SetBranchStatus("run",1);
248  t->SetBranchStatus(varname,1);
249 
250  Float_t ev, var, run;
251  t->SetBranchAddress("ev",&ev);
252  t->SetBranchAddress("run",&run);
253  t->SetBranchAddress(varname,&var);
254 
255  std::map<int, int> sigcnt;
256  std::map<int, int> bgcnt;
257  std::map<int, int> bgcnt2;
258 
259  std::vector<pair<double, int> > sigvals;
260  std::vector<pair<double, int> > bgvals;
261 
262  int Nsigval = els.GetN(), Nbgval = elb.GetN();
263  double min = 1e9, max=-1e9;
264 
265  // prepare the pairs of variable value, eventnumber for signal and count signals
266  for (i=0;i<els.GetN();++i)
267  {
268  t->GetEntry(els.GetEntry(i));
269  sigvals.push_back(std::make_pair(var, ev+run*2000));
270  if (var>max) max = var;
271  if (var<min) min = var;
272  }
273 
274  // prepare the pairs of variable value, eventnumber for background and count background
275  for (i=0;i<elb.GetN();++i)
276  {
277  t->GetEntry(elb.GetEntry(i));
278  bgvals.push_back(std::make_pair(var, ev+run*2000));
279  }
280 
281  // sort signals by first element = variable value
282  sort(sigvals.begin(), sigvals.end());
283 
284  i=0;
285  sigcnt.clear();
286  // find cut for current efficiency requirement left cut
287  while ( (sigcnt.size()/(double)Nsigev)<eff && i<Nsigval ) sigcnt[sigvals[i++].second]+=1;
288  double leftcut = sigvals[i].first;
289 
290  sigcnt.clear();
291  i=sigvals.size()-1;
292  // find cut for current efficiency requirement right cut
293  while ( (sigcnt.size()/(double)Nsigev)<eff && i>=0) sigcnt[sigvals[i--].second]+=1;
294  double rightcut = sigvals[i].first;
295 
296  bgcnt.clear();
297  bgcnt2.clear();
298 
299  for (i=0;i<bgvals.size();++i)
300  {
301  if (bgvals[i].first<=leftcut) bgcnt[bgvals[i].second];
302  if (bgvals[i].first>=rightcut) bgcnt2[bgvals[i].second];
303  if (bgvals[i].first>max) max = bgvals[i].first;
304  if (bgvals[i].first<min) min = bgvals[i].first;
305  }
306 
307  minh[idx]=min; maxh[idx]=max;
308 
309  t->SetBranchStatus("*",1);
310 
311 
312  double leftsupr = (Nbgev-bgcnt.size())/Nbgev;
313  double rightsupr = (Nbgev-bgcnt2.size())/Nbgev;
314 
315  cout <<varname<<" l:"<<leftcut<<"("<<leftsupr<<") r:"<<rightcut<<"("<<rightsupr<<")"<<endl;
316 
317  bestcut = rightcut;
318  if (leftsupr>rightsupr)
319  {
320  bestcut = leftcut;
321  return leftsupr;
322  }
323  return rightsupr;
324 }
els
TEventList * els
Definition: findcuts.C:52
run
Int_t run
Definition: autocutx.C:47
i
Int_t i
Definition: run_full.C:25
minh
double minh[MAX]
Definition: cutfinder_fullev.C:31
ev
int ev
Definition: GammaSpectraAnalysis.C:44
max
friend F32vec4 max(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:26
idx
int idx[MAX]
Definition: cutfinder_fullev.C:27
maxh
double maxh[MAX]
Definition: cutfinder_fullev.C:31
min
friend F32vec4 min(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:25
Nbgev
double Nbgev
Definition: cutfinder_fullev.C:34
elb
TEventList * elb
Definition: findcuts.C:52
Nsigev
double Nsigev
Definition: cutfinder_fullev.C:34
t
TTree * t
Definition: bump_analys.C:13
double bestSuppressionEvt ( TTree *  t,
TString  varname,
TEventList &  els,
TEventList &  elb,
double &  bestcut,
double  eff,
int  idx 
)

Definition at line 414 of file cutfinder_fullev.C.

References ev, i, idx, max(), maxh, min(), minh, Nbgev, Nsigev, and run.

415 {
416  int i;
417  t->SetBranchStatus("*",0);
418  t->SetBranchStatus("ev",1);
419  t->SetBranchStatus("run",1);
420  t->SetBranchStatus(varname,1);
421 
422  Float_t ev, var, run;
423  t->SetBranchAddress("ev",&ev);
424  t->SetBranchAddress("run",&run);
425  t->SetBranchAddress(varname,&var);
426 
427  std::map<int, int> sigcnt;
428  std::map<int, int> bgcnt;
429  std::map<int, int> bgcnt2;
430 
431  std::vector<pair<double, int> > sigvals;
432  std::vector<pair<double, int> > bgvals;
433 
434  int Nsigval = els.GetN(), Nbgval = elb.GetN();
435  double min = 1e9, max=-1e9;
436 
437  // prepare the pairs of variable value, eventnumber for signal and count signals
438  for (i=0;i<els.GetN();++i)
439  {
440  t->GetEntry(els.GetEntry(i));
441  sigvals.push_back(std::make_pair(var, ev+run*2000));
442  if (var>max) max = var;
443  if (var<min) min = var;
444  }
445 
446  // prepare the pairs of variable value, eventnumber for background and count background
447  for (i=0;i<elb.GetN();++i)
448  {
449  t->GetEntry(elb.GetEntry(i));
450  bgvals.push_back(std::make_pair(var, ev+run*2000));
451  }
452 
453  // sort signals by first element = variable value
454  sort(sigvals.begin(), sigvals.end());
455 
456  i=0;
457  sigcnt.clear();
458  // find cut for current efficiency requirement left cut
459  while ( (sigcnt.size()/(double)Nsigev)<eff && i<Nsigval ) sigcnt[sigvals[i++].second]+=1;
460  double leftcut = sigvals[i].first;
461 
462  sigcnt.clear();
463  i=sigvals.size()-1;
464  // find cut for current efficiency requirement right cut
465  while ( (sigcnt.size()/(double)Nsigev)<eff && i>=0) sigcnt[sigvals[i--].second]+=1;
466  double rightcut = sigvals[i].first;
467 
468  bgcnt.clear();
469  bgcnt2.clear();
470 
471  for (i=0;i<bgvals.size();++i)
472  {
473  if (bgvals[i].first<=leftcut) bgcnt[bgvals[i].second];
474  if (bgvals[i].first>=rightcut) bgcnt2[bgvals[i].second];
475  if (bgvals[i].first>max) max = bgvals[i].first;
476  if (bgvals[i].first<min) min = bgvals[i].first;
477  }
478 
479  minh[idx]=min; maxh[idx]=max;
480 
481  t->SetBranchStatus("*",1);
482 
483 
484  double leftsupr = (Nbgev-bgcnt.size())/Nbgev;
485  double rightsupr = (Nbgev-bgcnt2.size())/Nbgev;
486 
487  cout <<varname<<" l:"<<leftcut<<"("<<leftsupr<<") r:"<<rightcut<<"("<<rightsupr<<")"<<endl;
488 
489  if (leftcut!=leftcut || rightcut!=rightcut) return 0;
490 
491  bestcut = rightcut;
492  if (leftsupr>rightsupr)
493  {
494  bestcut = leftcut;
495  return leftsupr;
496  }
497  return rightsupr;
498 }
els
TEventList * els
Definition: findcuts.C:52
run
Int_t run
Definition: autocutx.C:47
i
Int_t i
Definition: run_full.C:25
minh
double minh[MAX]
Definition: cutfinder_fullev.C:31
ev
int ev
Definition: GammaSpectraAnalysis.C:44
max
friend F32vec4 max(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:26
idx
int idx[MAX]
Definition: cutfinder_fullev.C:27
maxh
double maxh[MAX]
Definition: cutfinder_fullev.C:31
min
friend F32vec4 min(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:25
Nbgev
double Nbgev
Definition: cutfinder_fullev.C:34
elb
TEventList * elb
Definition: findcuts.C:52
Nsigev
double Nsigev
Definition: cutfinder_fullev.C:34
t
TTree * t
Definition: bump_analys.C:13
int countEvents ( TTree *  t,
TString  ccut 
)

Definition at line 79 of file cutfinder_fullev.C.

References ev, evcnt, i, and run.

Referenced by cutfinder_fullev().

80 {
81  t->SetEventList(0);
82  t->SetBranchStatus("*",1);
83  t->Draw(">>el",ccut);
84  t->SetBranchStatus("*",0);
85  t->SetBranchStatus("evt",1);
86  t->SetBranchStatus("run",1);
87 
88  TEventList *el=(TEventList*)gDirectory->Get("el");
89 
90  Float_t ev,run;
91  t->SetBranchAddress("ev",&ev);
92  t->SetBranchAddress("run",&run);
93 
94  evcnt.clear();
95 
96  for (int i=0;i<el->GetN();++i)
97  {
98  t->GetEntry(el->GetEntry(i));
99  evcnt[ev+run*2000]+=1;
100  }
101  t->SetBranchStatus("*",1);
102 
103  return evcnt.size();
104 }
run
Int_t run
Definition: autocutx.C:47
i
Int_t i
Definition: run_full.C:25
ev
int ev
Definition: GammaSpectraAnalysis.C:44
evcnt
std::map< int, int > evcnt
Definition: cutfinder_fullev.C:24
t
TTree * t
Definition: bump_analys.C:13
int countEvents ( TTree *  t,
TEventList &  el 
)

Definition at line 105 of file cutfinder_fullev.C.

References ev, evcnt, i, and run.

106 {
107  t->SetEventList(&el);
108  t->SetBranchStatus("*",0);
109  t->SetBranchStatus("ev",1);
110  t->SetBranchStatus("run",1);
111 
112  Float_t ev, run;
113  t->SetBranchAddress("ev",&ev);
114  t->SetBranchAddress("run",&run);
115 
116  evcnt.clear();
117 
118  for (int i=0;i<el.GetN();++i)
119  {
120  t->GetEntry(el.GetEntry(i));
121  evcnt[ev+2000*run]+=1;
122  }
123  t->SetBranchStatus("*",1);
124 
125  return evcnt.size();
126 }
run
Int_t run
Definition: autocutx.C:47
i
Int_t i
Definition: run_full.C:25
ev
int ev
Definition: GammaSpectraAnalysis.C:44
evcnt
std::map< int, int > evcnt
Definition: cutfinder_fullev.C:24
t
TTree * t
Definition: bump_analys.C:13
int cutfinder_fullev ( TString  pre,
TString  ntp,
TString  sigcut,
TString  precut = "",
double  n0s = 50.,
double  n0b = 500.,
double  supr = 0.9 
)

Definition at line 534 of file cutfinder_fullev.C.

References bestSuppressionEvt(), c1, cnt, countEvents(), cut, elb, els, f, h1, h2, i, idx, MAX, maxh, minh, mycompare(), N0_bg, N0_sig, Nbgev, nsig, Nsigev, signi, t, TString, and vars.

535 {
536  gStyle->SetTitleY(0.993);
537  gStyle->SetTitleH(0.07);
538  gStyle->SetPadTopMargin(0.08);
539  gStyle->SetTitleBorderSize(0);
540  gStyle->SetOptStat(0);
541 
542  N0_sig = n0s*1000.;
543  N0_bg = n0b*1000.;
544 
545 
546  int i,j;
547  double sigl, sigh, bgl, bgh;
548 
549  TFile *f=new TFile("datafull/M"+pre+"_"+ntp+".root","READ");
550  //TFile *f=new TFile(fname,"READ");
551  TTree *t=(TTree*)f->Get(ntp);
552 
553  TEventList els("els");
554  TEventList elb("elb");
555 
556  TCanvas *c1=new TCanvas("c1","c1",10,10,1800,600);
557  c1->Divide(6,2);
558 
559  TObjArray* branches = t->GetListOfBranches();
560 
561  TString bgcut = "!("+sigcut+")";
562  //sigcut+="&&mct";
563  TString bstring = precut+" "+sigcut;
564  if (precut!="")
565  {
566  sigcut+="&&"+precut;
567  bgcut +="&&"+precut;
568  }
569  cout <<sigcut <<" "<<bgcut<<endl;
570  t->Draw(">>els",sigcut);
571  t->Draw(">>elb",bgcut);
572 
573  Nsigev = countEvents(t, els);
574  Nbgev = countEvents(t, elb);
575 
576  cout <<"SIG:"<<Nsigev<<" ev "<<els.GetN()<<" cn BG:"<<Nbgev<<" ev "<<elb.GetN()<<" cn"<<endl;
577 
578  int cnt=1;
579 
580  for (i=0;i<MAX;++i) { idx[i]=i; signi[i]=0.;}
581 
582  for(i=0; i<=branches->GetLast(); ++i)
583  {
584 
585  TBranch* branch = (TBranch*)branches->UncheckedAt(i);
586  vars[i]=branch->GetName();
587  //t->SetBranchStatus(vars[i],1);
588 
589  if (vars[i]=="ev") continue;
590  if (vars[i]=="mode") continue;
591  if (vars[i]=="run") continue;
592  if (vars[i].Contains("pdg")) continue;
593  if (vars[i].Contains("beam")) continue;
594  if (vars[i].Contains("mct")) continue;
595  if (vars[i].BeginsWith("t") && vars[i]!="thr") continue;
596 
597  signi[i] = bestSuppressionEvt(t, vars[i], els, elb, cut[i], supr, i);
598 
599  }
600 
601  //t->SetBranchStatus("*",1);
602 
603  cout <<"\n\nBEST 12 vars:"<<endl<<endl;
604 
605  std::vector<int> myidx (idx, idx+MAX);
606 
607  std::sort(myidx.begin(), myidx.end(), mycompare);
608 
609  TLine l;
610  l.SetLineColor(6);
611  l.SetLineStyle(2);
612  l.SetLineWidth(2);
613  TLatex lt;
614  lt.SetTextSize(0.06);
615 
616  for (j=0;j<12;++j)
617  {
618  i=myidx[j];
619  cout << vars[i]<<" : qa = "<<signi[i]<<" cut = "<<cut[i]<<endl;
620  if (j<12)
621  {
622 /* findLimits(t,vars[i],sigcut, sigl, sigh);
623  findLimits(t,vars[i],bgcut, bgl, bgh);
624 
625  if (sigl>bgl) sigl=bgl;
626  if (sigh<bgh) sigh=bgh;
627  */
628  TH1F h1("h1",vars[i],200,minh[i],maxh[i]);
629  TH1F h2("h2",vars[i],200,minh[i],maxh[i]);
630  h2.SetLineColor(2);
631 
632  t->SetEventList(&els);
633  t->Project("h1",vars[i]);
634  t->SetEventList(&elb);
635  t->Project("h2",vars[i]);
636 
637  h1.Scale(1.0/h1.Integral());
638  h2.Scale(1.0/h2.Integral());
639  h1.SetTitleSize(0.05);
640 
641  c1->cd(j+1);
642  double maxi = h1.GetMaximum();
643  if (h2.GetMaximum()>maxi ) maxi = h2.GetMaximum();
644  maxi*=1.1;
645  h1.SetMaximum(maxi);
646  h2.SetMaximum(maxi);
647 
648  h1.DrawNormalized();
649  h2.DrawNormalized("same");
650 
651  double axmin = h1.GetXaxis()->GetXmin(),axmax = h1.GetXaxis()->GetXmax();
652 
653  l.DrawLine(cut[i],0, cut[i], maxi*0.9);
654 
655  lt.DrawLatex(axmin+(axmax-axmin)*0.6,1.01*maxi,TString::Format("sig = %6.4f",signi[i]));
656 
657  //c1->Update();
658  }
659 
660  }
661 
662  c1->Update();
663 
664  float nsig = countEvents(t,els);
665  float nbg = countEvents(t,elb);
666 
667  cout <<"SIG EVT : "<<nsig<<" "<<"BG EVT : "<<nbg<<endl;
668  cout <<"SIG EFF : "<<nsig/N0_sig<<" "<<"BG EFF : "<< nbg/N0_bg<<endl;
669  return 0;
670 }
els
TEventList * els
Definition: findcuts.C:52
MAX
#define MAX
Definition: cutfinder_fullev.C:20
i
Int_t i
Definition: run_full.C:25
minh
double minh[MAX]
Definition: cutfinder_fullev.C:31
vars
TString vars[MAX]
Definition: cutfinder_fullev.C:28
idx
int idx[MAX]
Definition: cutfinder_fullev.C:27
N0_bg
double N0_bg
Definition: cutfinder_fullev.C:33
h1
TH1F * h1
Definition: anal_hit_digi_cluster_fwendcap.C:121
maxh
double maxh[MAX]
Definition: cutfinder_fullev.C:31
cut
double cut[MAX]
Definition: cutfinder_fullev.C:29
N0_sig
double N0_sig
Definition: cutfinder_fullev.C:33
countEvents
int countEvents(TTree *t, TString ccut)
Definition: cutfinder_fullev.C:79
f
TFile * f
Definition: bump_analys.C:12
signi
double signi[MAX]
Definition: cutfinder_fullev.C:26
Nbgev
double Nbgev
Definition: cutfinder_fullev.C:34
c1
c1
Definition: plot_dirc.C:35
mycompare
bool mycompare(int i, int j)
Definition: cutfinder_fullev.C:36
nsig
int nsig
Definition: toy_core.C:46
h2
TH1F * h2
Definition: anal_hit_digi_cluster_fwendcap.C:122
elb
TEventList * elb
Definition: findcuts.C:52
cnt
Int_t cnt
Definition: hist-t7.C:106
Nsigev
double Nsigev
Definition: cutfinder_fullev.C:34
t
TTree * t
Definition: bump_analys.C:13
bestSuppressionEvt
double bestSuppressionEvt(TTree *t, TString precut, TString varname, TString sigcut, TString bgcut, double &bestcut, double eff, int idx)
Definition: cutfinder_fullev.C:235
double diffQA ( TH1F &  ht1,
TH1F &  ht2 
)

Definition at line 128 of file cutfinder_fullev.C.

References fabs(), h1, h2, i, and n.

129 {
130  int n=ht1.GetNbinsX();
131  TH1F h1=ht1;
132  TH1F h2=ht2;
133  h1.Scale(1.0/h1.Integral());
134  h2.Scale(1.0/h2.Integral());
135  h1.Add(&h2,-1);
136 
137  double sum=0.;
138 
139  for (int i=1;i<=n;++i)
140  {
141  sum+=fabs(h1.GetBinContent(i));
142  }
143 
144  return sum;
145 }
i
Int_t i
Definition: run_full.C:25
n
int n
Definition: HypStatDecay.cxx:21
h1
TH1F * h1
Definition: anal_hit_digi_cluster_fwendcap.C:121
fabs
friend F32vec4 fabs(const F32vec4 &a)
Definition: P4_F32vec4.h:47
h2
TH1F * h2
Definition: anal_hit_digi_cluster_fwendcap.C:122
void findLimits ( TTree *  t,
TString  var,
TString  ccut,
double &  low,
double &  high,
double  frac = 0.98 
)

Definition at line 500 of file cutfinder_fullev.C.

References htemp, and i.

501 {
502  low = t->GetMinimum(var);
503  high = t->GetMaximum(var);
504  return;
505 
506  double miss = (1.-frac)/2;
507 
508  t->SetEventList(0);
509  TEventList el("el");
510  t->Draw(">>el",ccut);
511  t->SetEventList(&el);
512 
513  double llow = t->GetMinimum(var);
514  double lhigh = t->GetMaximum(var);
515 
516  TH1F htemp("htemp","",500,llow, lhigh);
517  t->Project("htemp",var);
518 
519  double sum=0.;
520  double integ = htemp.Integral();
521 
522  int i=1;
523 
524  while (sum<miss) sum+=htemp.GetBinContent(i++)/integ;
525  low = htemp.GetBinCenter(i-2);
526  sum=0.; i=500;
527 
528  while (sum<miss) sum+=htemp.GetBinContent(i--)/integ;
529  high = htemp.GetBinCenter(i+2);
530 
531  t->SetEventList(0);
532 }
i
Int_t i
Definition: run_full.C:25
low
int low
Definition: anaMvdDigi.C:51
high
double high
Definition: full_core_ntp.C:35
t
TTree * t
Definition: bump_analys.C:13
htemp
TH1D * htemp[nsteps]
Definition: dedx_bands.C:61
bool mycompare ( int  i,
int  j 
)

Definition at line 36 of file cutfinder_fullev.C.

References i, and signi.

Referenced by cutfinder_fullev().

37 {
38  return signi[i]>signi[j];
39 }
i
Int_t i
Definition: run_full.C:25
signi
double signi[MAX]
Definition: cutfinder_fullev.C:26

Variable Documentation

double cut[MAX]

Definition at line 29 of file cutfinder_fullev.C.

Referenced by cutfinder_fullev().

std::map<int, int> evcnt

Definition at line 24 of file cutfinder_fullev.C.

Referenced by countEvents().

int idx[MAX]

Definition at line 27 of file cutfinder_fullev.C.

Referenced by bestSeparation(), bestSuppressionEvt(), and cutfinder_fullev().

double maxh[MAX]

Definition at line 31 of file cutfinder_fullev.C.

Referenced by bestSeparation(), bestSuppressionEvt(), and cutfinder_fullev().

double minh[MAX]

Definition at line 31 of file cutfinder_fullev.C.

Referenced by bestSeparation(), bestSuppressionEvt(), and cutfinder_fullev().

double N0_bg

Definition at line 33 of file cutfinder_fullev.C.

Referenced by cutfinder_fullev().

double N0_sig

Definition at line 33 of file cutfinder_fullev.C.

Referenced by cutfinder_fullev().

double Nbgev

Definition at line 34 of file cutfinder_fullev.C.

Referenced by bestSeparation(), bestSuppressionEvt(), and cutfinder_fullev().

double Nsigev

Definition at line 34 of file cutfinder_fullev.C.

Referenced by bestSeparation(), bestSuppressionEvt(), and cutfinder_fullev().

double signi[MAX]

Definition at line 26 of file cutfinder_fullev.C.

Referenced by cutfinder_fullev(), and mycompare().

TString vars[MAX]

Definition at line 28 of file cutfinder_fullev.C.

Referenced by cutfinder_fullev().

Generated on Sun Aug 18 2019 03:02:27 for FairRoot/PandaRoot by   doxygen 1.8.8