FairRoot/PandaRoot
full_core_ntp.C
Go to the documentation of this file.
1 #include <iostream>
2 #include <fstream>
3 #include <vector>
4 #include <map>
5 
6 #include "TFile.h"
7 #include "TTree.h"
8 #include "TChain.h"
9 #include "TRandom3.h"
10 #include "TDatabasePDG.h"
11 #include "TParticlePDG.h"
12 #include "TCanvas.h"
13 #include "TH1F.h"
14 #include "TH2F.h"
15 #include "TLegend.h"
16 #include "TParticle.h"
17 #include "TClonesArray.h"
18 #include "TF1.h"
19 #include "TStyle.h"
20 #include "TNtuple.h"
21 #include "TLatex.h"
22 #include "TGaxis.h"
23 #include "TEventList.h"
24 
25 using std::cout;
26 using std::endl;
27 
28 TRandom3 fRand;
29 TLorentzVector fIni;
30 
31 // *** some program parameters
32 const int fMAX = 100; // *** max tracks per event
33 int nbins = 300; // *** bins of mass histograms
34 double low = 0.0; // *** lower edge of histo
35 double high = 6.; // *** higher edge of histo (will be modified acc. to resonance mass)
36 
37 std::map<int, int> evcnt;
38 std::map<int, int> treeidx;
39 std::map<TString, TString> cuts;
40 std::map<TString, TString> modes;
41 std::map<TString, TString> mctvar;
42 
43 TString ntpname[10] = {"nphi","nlam","n2e","njpsi1","njpsi2","netac2","nd01","ndpm1","nds1","nlamc"};
44 TString plotname[10] = {"phim","lamm","eem","jpsim","jpsim","etacm","d0m","dpmm","dsm","lamcm"};
45 
46 // *** setup style
47 void setStyle()
48 {
49  gStyle->SetPadTopMargin(0.13);
50  gStyle->SetPadBottomMargin(0.18);
51  gStyle->SetPadLeftMargin(0.2);
52  gStyle->SetPadRightMargin(0.016);
53 
54  gStyle->SetLabelSize(0.076,"X");
55  gStyle->SetLabelSize(0.07,"Y");
56 
57  gStyle->SetStatY(0.87);
58  gStyle->SetStatX(0.986);
59  gStyle->SetStatW(0.4);
60  gStyle->SetStatH(0.12);
61  gStyle->SetOptStat("e"); // only entries
62 
63  gStyle->SetTitleH(0.1);
64  gStyle->SetTitleX(0.08);
65  gStyle->SetTitleOffset(1.1,"x");
66  gStyle->SetTitleXSize(0.076);
67 
68  gStyle->SetNdivisions(505);
69 }
70 
71 double getPos(TH1* h)
72 {
73  return (h->GetBinCenter(h->GetNbinsX())-h->GetBinCenter(1))*0.05+h->GetBinCenter(1);
74 }
75 
76 void configHisto(TH1* h)
77 {
78  h->SetFillStyle(3001);
79  h->SetFillColor(602);
80  h->SetLineStyle(3);
81  h->SetLineColor(602);
82 /* TGaxis *ax = (TGaxis*)h->GetYaxis();
83  ax->SetMaxDigits(2);*/
84 }
85 
86 void init_modes()
87 {
88  modes["nphi"] = "000";
89  modes["nlam"] = "400";
90  modes["n2e"] = "300";
91  modes["njpsi1"] = "200";
92  modes["njpsi2"] = "201";
93  modes["netac2"] = "221";
94  modes["nd01"] = "100";
95  modes["ndpm1"] = "120";
96  modes["nds1"] = "140";
97  modes["nlamc"] = "420";
98 
99  modes["DPM"] = "900";
100 
101  mctvar["nphi"] = "phi";
102  mctvar["nlam"] = "lam";
103  mctvar["n2e"] = "ee";
104  mctvar["njpsi1"]= "jpsi";
105  mctvar["njpsi2"]= "jpsi";
106  mctvar["netac2"]= "etac";
107  mctvar["nd01"] = "d0";
108  mctvar["ndpm1"] = "dpm";
109  mctvar["nds1"] = "ds";
110  mctvar["nlamc"] = "lamc";
111 }
112 
113 void init_cuts_1000(TString energy)
114 {
115  if (energy=="24")
116  {
117  cuts["nphi"] = "eslnpidk>1&&phipcm>0.5&&phipcm<0.75";
118  cuts["nlam"] = "lampcm>0.38&&lampcm<0.48&&eslnpidp>0&&lamp>0.6";
119  cuts["n2e"] = "eslnpide>1";
120 
121  }
122  else if (energy=="38")
123  {
124  cuts["nphi"] = "eslnpidk>1&&phipcm>1.4";
125  cuts["nlam"] = "eslnpidp>0&&lampcm>1.25&&esthr>0.93&&lampt>0.75";
126  cuts["n2e"] = "eslnpide>1";
127  cuts["njpsi1"] = "eslnpide>1&&jpsipcm<0.9";
128  cuts["njpsi2"] = "eslnpidmu>1&&jpsipcm<0.9";
129  cuts["netac2"] = "eslnpidk>0&&eslnpidpi>2";
130  cuts["nd01"] = "eslnpidk>0&&d0pcm<0.4&&d0d0pt>0.4";
131  cuts["ndpm1"] = "eslnpidk>0&&dpmpcm<0.3";
132 
133  }
134  else if (energy=="45")
135  {
136  cuts["nphi"] = "eslnpidk>1&&phipcm>1.55&&esthr>0.85";
137  cuts["nlam"] = "esl1npidp>0&&esthr>0.93&&lampt>1";
138  cuts["n2e"] = "eslnpide>1&&mmiss>-1";
139  cuts["njpsi1"] = "eslnpide>1&&jpsipcm<1.8&&detemcmax>2.5";
140  cuts["njpsi2"] = "eslnpidmu>1&&jpsipcm<1.5";
141  cuts["netac2"] = "eslnpidk>0&&eslnpidpi>2&&etacd1pt>0.35&&etacd0pt>0.35";
142  cuts["nd01"] = "eslnpidk>0&&d0pt>1&&mmiss>1.5&&d0d0pidk>0.25";
143  cuts["ndpm1"] = "eslnpidk>0&&dpmpcm>1.1&&dpmpcm<1.5&&dpmpt>0.9&&dpmd0pidk>0.3";
144  cuts["nds1"] = "eslnpidk>1&&dspcm>1&&dspcm<1.17&&dspt>0.5";
145 
146  }
147  else if (energy=="55")
148  {
149  cuts["nphi"] = "eslnpidk>1&&phipcm>1.5&&esthr>0.89";
150  cuts["nlam"] = "eslnpidp>0&&esthr>0.97&&lampcm>2.1&&lampt>1";
151  cuts["n2e"] = "eslnpide>1&&mmiss>-1";
152  cuts["njpsi1"] = "eslnpide>1&&jpsipcm<2.2&&detemcmax>4&&jpsid1tht>0.07";
153  cuts["njpsi2"] = "eslnpidmu>1&&jpsipcm<2";
154  cuts["netac2"] = "eslnpidk>0&&eslnpidpi>2&&esptmax>1.1&&esnpart<15";
155  cuts["nd01"] = "eslnpidk>0&&d0pcm>1.9&&mmiss>1.5&&d0pt>1.6&&d0pcm<2.15";
156  cuts["ndpm1"] = "eslnpidk>0&&dpmpcm>1.95&&dpmpcm<2.15&&dpmpt>0.8&&dpmd0pt>0.4";
157  cuts["nds1"] = "eslnpidk>1&&dspcm>1.85&&dspcm<2.1&&dspt>1.3";
158  cuts["nlamc"] = "eslnpidk>0&&eslnpidp>0&&lamcpt>.9&&lamcd0pidp>0.21&&lamcd0thtcm<2";
159 
160  }
161  else
162  {
163  cuts["nphi"] = "1";
164  cuts["nlam"] = "1";
165  cuts["n2e"] = "1";
166  cuts["njpsi1"] = "1";
167  cuts["njpsi2"] = "1";
168  cuts["netac2"] = "1";
169  cuts["nd01"] = "1";
170  cuts["ndpm1"] = "1";
171  cuts["nds1"] = "1";
172  cuts["nlamc"] = "1";
173  }
174  // factor 1000
175 /* cuts["nphi"] = "fw1>0.0129387&&phipt>0.2";
176  cuts["nlam"] = "fw3>0.2&&fw4>0.3&&sumpc>0.9";
177  cuts["nee"] = "1";*/
178 
179  // factor 1000
180 /* cuts["nphi"] = "phipcm>1.31489&&phipt>0.698733";
181  cuts["nlam"] = "fw2>0.688083&&lampt>0.990166";
182  cuts["nee"] = "1";
183  cuts["njpsi1"] = "1";
184  cuts["njpsi2"] = "1";
185  cuts["netac"] = "ptmax>0.84&&sumptc>2.9";
186  cuts["nd0"] = "d0pcm<0.35&&ptmax>0.8";
187  cuts["ndpm"] = "dpcm<0.34&&ptmax>0.67&&fw1>0";*/
188 
189  // first attempt
190 /* cuts["netac"] = "ptmax>0.7&&sumptc>2.7";
191  cuts["nd0"] = "d0pcm<0.33&&ptmax>0.75";
192  cuts["ndpm"] = "dpcm<0.31&&ptmax>0.6&&pmin>0.1";*/
193  // factor 1000
194 // cuts["nphi"] = "phipcm>1.9&&phipt>1";
195 // cuts["nlam"] = "lampcm>1.8&&thr>0.95&&lampt>1";
196 // cuts["nee"] = "1";
197 // cuts["njpsi1"] = "sumpc>3.5";
198 // cuts["njpsi2"] = "sumpc>3.5";
199 // cuts["netac"] = "ptmax>1&&sumpc>4";
200 // cuts["nd0"] = "abs(d0pcm-1.24)<0.15&&ptmax>1.05&&sumptc>3.1";
201 // cuts["ndpm"] = "abs(dpcm-1.24)<0.12&&dpt>1&&sumptc>3";
202 // cuts["nds"] = "abs(dspcm-1.1)<0.09&&ptmax>0.93";
203 
204  // first attempt cuts["nphi"] = "phipcm>1.9&&phipt>1";
205 /* cuts["nd0"] = "abs(d0pcm-1.24)<0.15&&ptmax>1.2";
206  cuts["ndpm"] = "abs(dpcm-1.24)<0.15&&dpt>1&&pmin>0.15";
207  cuts["nds"] = "abs(dspcm-1.1)<0.1&&ptmax>0.8";*/
208  // factor 1000
209 /* cuts["nphi"] = "phipcm>2.4&&phipt>0.5";
210  cuts["nlam"] = "lampcm>0.23&&thr>0.98&&lampt>1.2";
211  cuts["nee"] = "1";
212  cuts["njpsi1"] = "sumpc>4";
213  cuts["njpsi2"] = "sumpc>4";
214  cuts["netac"] = "ptmax>1.1&&sumpc>4.8";
215  cuts["nd0"] = "ptmax>1&&d0pcm>1.8&&d0pt>1.5";
216  cuts["ndpm"] = "ptmax>1&&dpcm>1.8&&dpt>1.5";
217  cuts["nds"] = "dspcm>1.8&&ptmax>0.8&&dspt>1.2";
218  cuts["nlamc"] = "abs(lamcpcm-1.53)<0.12&&lamcpt>1.3&&fw1>0";*/
219  // first attempt
220 /* cuts["nlamc"] = "abs(lamcpcm-1.53)<0.15&&lamcpt>1.3";*/
221 }
222 void init_cuts_95(TString energy)
223 {
224  if (energy=="24")
225  {
226  cuts["nphi"] = "eslnpidk>1&&phipcm>0.45&&phipcm<0.78";
227  cuts["nlam"] = "eslnpidp>0&&lampcm<0.55&&lampcm>0.25&&mmiss<1.2";
228  cuts["n2e"] = "eslnpide>1&&mmiss>-1";
229 
230  }
231  else if (energy=="38")
232  {
233  cuts["nphi"] = "eslnpidk>1&&phipcm>1.2";
234  cuts["nlam"] = "esthr>0.85&&lampcm>1.15&&lampt>0.2&&eslnpidp>0";
235  cuts["n2e"] = "eslnpide>1&&mmiss>-1";
236  cuts["njpsi1"] = "eslnpide>1&&jpsipcm<1";
237  cuts["njpsi2"] = "eslnpidmu>1&&jpsipcm<1";
238  cuts["netac2"] = "etacpcm<1.25&&essumptc>1.75";
239  cuts["nd01"] = "eslnpidk>0&&d0pcm<0.4&&d0d0pt>0.25";
240  cuts["ndpm1"] = "eslnpidk>0";
241 
242  }
243  else if (energy=="45")
244  {
245  cuts["nphi"] = "eslnpidk>1&&phipcm>1.55";
246  cuts["nlam"] = "esthr>0.87&&eslnpidp>0&&lampcm>1.5";
247  cuts["n2e"] = "eslnpide>1&&mmiss>-1";
248  cuts["njpsi1"] = "eslnpide>1&&jpsipcm<1.8&&detemcmax>2.5";
249  cuts["njpsi2"] = "eslnpidmu>1&&jpsipcm<1.5";
250  cuts["netac2"] = "etacpcm<1.25&&etacd1p>0.25";
251  cuts["nd01"] = "eslnpidk>0&&d0pcm<1.6&&d0pcm>0.8&&d0d1p<6";
252  cuts["ndpm1"] = "eslnpidk>0";
253  cuts["nds1"] = "eslnpidk>0&&dspcm>0.4&&dsp<8&&dsthtcm>0.5";
254 
255  }
256  else if (energy=="55")
257  {
258  cuts["nphi"] = "eslnpidk>1&&phipcm>1.75&&phipt>0.5";
259  cuts["nlam"] = "esthr>0.9&&lampcm>1&&mmiss<2.4";
260  cuts["n2e"] = "eslnpide>0&&mmiss>-1";
261  cuts["njpsi1"] = "eslnpide>1&&jpsipcm<1.9&&detemcmax>3.5";
262  cuts["njpsi2"] = "eslnpidmu>1&&jpsipcm<2.5";
263  cuts["netac2"] = "etacpcm<2&&etacd1p>0.4&&etacecm>2.85&&eslnpide<1";
264  cuts["nd01"] = "eslnpidk>0&&d0pcm<2.3&&d0pcm>0.95";
265  cuts["ndpm1"] = "eslnpidk>0";
266  cuts["nds1"] = "eslnpidk>0&&dspcm>0.5&&dsthtcm>0.5";
267  cuts["nlamc"] = "mmiss>2";
268 
269  }
270  else
271  {
272  cuts["nphi"] = "1";
273  cuts["nlam"] = "1";
274  cuts["n2e"] = "1";
275  cuts["njpsi1"] = "1";
276  cuts["njpsi2"] = "1";
277  cuts["netac2"] = "1";
278  cuts["nd01"] = "1";
279  cuts["ndpm1"] = "1";
280  cuts["nds1"] = "1";
281  cuts["nlamc"] = "1";
282  }
283 }
284 
285 int cnt_events(TChain* t, TString cut, bool globcnt=false)
286 {
287  treeidx.clear();
288  Float_t ev;
289  t->Draw(">>el",cut);
290 
291  TEventList *el=(TEventList*)gDirectory->Get("el");
292 
293  t->SetBranchStatus("*",0);
294  t->SetBranchStatus("ev",1);
295  t->SetBranchAddress("ev",&ev);
296 
297  std::map<int, int> sigcnt;
298  for (int j=0;j<el->GetN();++j)
299  {
300  t->GetEntry(el->GetEntry(j));
301  int tidx = t->GetTreeNumber();
302  sigcnt[(int)ev+tidx*1000]+=1;
303  if (globcnt) evcnt[(int)ev+tidx*1000]+=1;
304  }
305  t->SetBranchStatus("*",1);
306  return sigcnt.size();
307 }
308 
309 void fillHistos(TChain *t, TString var, TString cut, TString all, TString sel, TString sig)
310 {
311  t->Project(all,var,cut);
312  t->Project(sel,var,cut+" && tag");
313  t->Project(sig,var,cut+" && "+mctvar[t->GetName()]+"mct");
314 
315 }
316 double DbMass(TString name)
317 {
318  TDatabasePDG *fPdg = TDatabasePDG::Instance();
319  if (fPdg->GetParticle(name)) return fPdg->GetParticle(name)->Mass();
320 }
321 
322 int full_core_ntp(TString fname, int nev=50, bool tagonly=false)
323 {
324  int i,j;
325  nev*=1000;
326 
327  setStyle();
328  init_modes();
329 
330  TString mode(fname(fname.Last('/')+4,3));
331 
332  TString energy(fname(fname.Last('/')+2,2));
333  double sqrts=energy.Atof();
334 
335  if (tagonly)
336  init_cuts_1000("");
337  else
338  init_cuts_1000(energy);
339 
340  cout <<energy<<endl;
341 
342 
343  //TFile *f=new TFile(fname);
344  TChain *t[10];
345 
346  for (i=0;i<10;++i)
347  {
348  t[i]=new TChain(ntpname[i]);
349  t[i]->Add(fname);
350 /* t[i]=(TTree*)f->Get(ntpname[i]);
351  if (t==0) cout <<"not found :"<<ntpname[i].Data()<<endl;*/
352  }
353 
354  TCanvas *c2=new TCanvas("c2","c2",20,20,1800,800);
355  c2->Divide(5,2,0.005,0.01);
356 
357 
358  double PhiMass = DbMass("phi");
359 
360  double D0Mass = DbMass("D0");
361  double DpmMass = DbMass("D+");
362  double DsMass = DbMass("D_s+");
363 
364  double JpsiMass = DbMass("J/psi");
365  double EtacMass = DbMass("eta_c");
366  double Chic0Mass = DbMass("chi_0c");
367 
368  double LamMass = DbMass("Lambda0");
369  double LamcMass = DbMass("Lambda_c+");
370 
371  // *** resonances
372  double fPhiQaMin = PhiMass - 0.1;
373  double fPhiQaMax = PhiMass + 0.1;
374 
375  // *** open charm
376  double fD0QaMin = D0Mass - 0.3;
377  double fD0QaMax = D0Mass + 0.3;
378 
379  double fDpmQaMin = DpmMass - 0.3;
380  double fDpmQaMax = DpmMass + 0.3;
381 
382  double fDsQaMin = DsMass - 0.25;
383  double fDsQaMax = DsMass + 0.25;
384 
385  // *** baryons
386  double fLamcQaMin = LamcMass - 0.2;
387  double fLamcQaMax = LamcMass + 0.2;
388 
389  double fLamQaMin = LamMass - 0.15;
390  double fLamQaMax = LamMass + 0.15;
391 
392  // *** charmonia
393  double fEtacQaMin = EtacMass - 0.65;
394  double fEtacQaMax = EtacMass + 0.65;
395 
396  double fJpsiQaMin = JpsiMass - 1.00;
397  double fJpsiQaMax = JpsiMass + 0.75;
398 
399 
400  // *** electro-magnetic
401  double fEcm = sqrts/10.;
402  double f2eQaMin = fEcm - 1.5;
403  double f2eQaMax = fEcm + 1.0;
404 
405 
406  double low2=0.5;
407  // ****** Jpsi -> ee
408  TH1F *h_jpsi1=new TH1F("h_jpsi1","J/#psi #rightarrow e^{+}e^{-}",nbins,fJpsiQaMin,fJpsiQaMax);
409  TH1F *h_jpsi1sel=new TH1F("h_jpsi1sel","J/psi",nbins,fJpsiQaMin,fJpsiQaMax);
410  TH1F *h_jpsi1sig=new TH1F("h_jpsi1sig","J/psi",nbins,fJpsiQaMin,fJpsiQaMax);
411  h_jpsi1sig->SetLineColor(2);
412  configHisto(h_jpsi1sel);
413 
414  // ****** Jpsi -> mu mu
415  TH1F *h_jpsi2=new TH1F("h_jpsi2","J/#psi #rightarrow #mu^{+}#mu^{-}",nbins,fJpsiQaMin,fJpsiQaMax);
416  TH1F *h_jpsi2sel=new TH1F("h_jpsi2sel","J/psi",nbins,fJpsiQaMin,fJpsiQaMax);
417  TH1F *h_jpsi2sig=new TH1F("h_jpsi2sig","J/psi",nbins,fJpsiQaMin,fJpsiQaMax);
418  h_jpsi2sig->SetLineColor(2);
419  configHisto(h_jpsi2sel);
420 
421  // ****** eta_c -> Ks K pi
422  TH1F *h_etac=new TH1F("h_etac","#eta_{c} #rightarrow K_{S} K^{+} #pi^{-}",nbins,fEtacQaMin,fEtacQaMax);
423  TH1F *h_etacsel=new TH1F("h_etacsel","etac",nbins,fEtacQaMin,fEtacQaMax);
424  TH1F *h_etacsig=new TH1F("h_etacsig","etac",nbins,fEtacQaMin,fEtacQaMax);
425  h_etacsig->SetLineColor(2);
426  configHisto(h_etacsel);
427 
428  // ****** D0 -> K- pi+
429  TH1F *h_d0=new TH1F("h_d0","D^{0} #rightarrow K^{-}#pi^{+}",nbins,fD0QaMin,fD0QaMax);
430  TH1F *h_d0sel=new TH1F("h_d0sel","D0",nbins,fD0QaMin,fD0QaMax);
431  TH1F *h_d0sig=new TH1F("h_d0sig","D0",nbins,fD0QaMin,fD0QaMax);;
432  h_d0sig->SetLineColor(2);
433  configHisto(h_d0sel);
434 
435  // ****** D+ -> K- pi+ pi+
436  TH1F *h_dpm=new TH1F("h_dpm","D^{+} #rightarrow K^{-}#pi^{+}#pi^{+}",nbins,fDpmQaMin,fDpmQaMax);
437  TH1F *h_dpmsel=new TH1F("h_dpmsel","D#pm",nbins,fDpmQaMin,fDpmQaMax);
438  TH1F *h_dpmsig=new TH1F("h_dpmsig","D#pm",nbins,fDpmQaMin,fDpmQaMax);
439  h_dpmsig->SetLineColor(2);
440  configHisto(h_dpmsel);
441 
442  // ****** Ds+ -> K+ K- pi+
443  TH1F *h_ds=new TH1F("h_ds","D_{s}^{+} #rightarrow K^{+}K^{-}#pi^{+}",nbins,fDsQaMin,fDsQaMax);
444  TH1F *h_dssel=new TH1F("h_dssel","Ds",nbins,fDsQaMin,fDsQaMax);
445  TH1F *h_dssig=new TH1F("h_dssig","Ds",nbins,fDsQaMin,fDsQaMax);
446  h_dssig->SetLineColor(2);
447  configHisto(h_dssel);
448 
449  // ****** phi -> K+ K-
450  TH1F *h_phi=new TH1F("h_phi","#phi #rightarrow K^{+}K^{-}",nbins,fPhiQaMin,fPhiQaMax);
451  TH1F *h_phisel=new TH1F("h_phisel","Phi",nbins,fPhiQaMin,fPhiQaMax);
452  TH1F *h_phisig=new TH1F("h_phisig","Phi",nbins,fPhiQaMin,fPhiQaMax);
453  h_phisig->SetLineColor(2);
454  configHisto(h_phisel);
455 
456  // ****** Lambda_c -> p K- pi+
457  TH1F *h_lamc=new TH1F("h_lamc","#Lambda_{c} #rightarrow pK^{-}#pi^{+}",nbins,fLamcQaMin,fLamcQaMax);
458  TH1F *h_lamcsel=new TH1F("h_lamcsel","Lam_c",nbins,fLamcQaMin,fLamcQaMax);
459  TH1F *h_lamcsig=new TH1F("h_lamcsig","Lam_c",nbins,fLamcQaMin,fLamcQaMax);
460  h_lamcsig->SetLineColor(2);
461  configHisto(h_lamcsel);
462 
463  // ****** Lambda -> p pi-
464  TH1F *h_lam=new TH1F("h_lam","#Lambda #rightarrow p#pi^{-}",nbins,fLamQaMin,fLamQaMax);
465  TH1F *h_lamsel=new TH1F("h_lamsel","Lam",nbins,fLamQaMin,fLamQaMax);
466  TH1F *h_lamsig=new TH1F("h_lamsig","Lam",nbins,fLamQaMin,fLamQaMax);
467  h_lamsig->SetLineColor(2);
468  configHisto(h_lamsel);
469 
470  // ****** pp -> e+ e-
471  TH1F *h_ee=new TH1F("h_ee","p#bar{p} #rightarrow e^{+}e^{-}",nbins,f2eQaMin,f2eQaMax);
472  TH1F *h_eesel=new TH1F("h_eesel","ee",nbins,f2eQaMin,f2eQaMax);
473  TH1F *h_eesig=new TH1F("h_eesig","ee",nbins,f2eQaMin,f2eQaMax);
474  h_eesig->SetLineColor(2);
475  configHisto(h_eesel);
476 
477  h_jpsi1->SetXTitle("m(e^{+}e^{-}) [GeV/c^{2}]");
478  h_jpsi2->SetXTitle("m(#mu^{+}#mu^{-}) [GeV/c^{2}]");
479  h_etac->SetXTitle("m(K_{S} K^{+}#pi^{-}) [GeV/c^{2}]");
480 
481  h_d0->SetXTitle("m(K^{-}#pi^{+}) [GeV/c^{2}]");
482 
483  h_dpm->SetXTitle("m(K^{-}#pi^{+}#pi^{+}) [GeV/c^{2}]");
484 
485  h_ds->SetXTitle("m(K^{+}K^{-}#pi^{+}) [GeV/c^{2}]");
486 
487  h_phi->SetXTitle("m(K^{+}K^{-}) [GeV/c^{2}]");
488  h_lamc->SetXTitle("m(pK^{-}#pi^{+}) [GeV/c^{2}]");
489  h_lam->SetXTitle("m(p#pi^{-}) [GeV/c^{2}]");
490  h_ee->SetXTitle("m(e^{+}e^{-}) [GeV/c^{2}]");
491 
492  // *** make histos
493  fillHistos(t[0], "phim", cuts["nphi"], "h_phi", "h_phisel", "h_phisig");
494  fillHistos(t[1], "lamm", cuts["nlam"], "h_lam", "h_lamsel", "h_lamsig");
495  fillHistos(t[2], "eem", cuts["n2e"], "h_ee", "h_eesel", "h_eesig");
496 
497  if (sqrts>24)
498  {
499  fillHistos(t[3], "jpsim", cuts["njpsi1"], "h_jpsi1", "h_jpsi1sel", "h_jpsi1sig");
500  fillHistos(t[4], "jpsim", cuts["njpsi2"], "h_jpsi2", "h_jpsi2sel", "h_jpsi2sig");
501  fillHistos(t[5], "etacm", cuts["netac2"], "h_etac", "h_etacsel", "h_etacsig");
502  fillHistos(t[6], "d0m", cuts["nd01"], "h_d0", "h_d0sel", "h_d0sig");
503  fillHistos(t[7], "dpmm", cuts["ndpm1"], "h_dpm", "h_dpmsel", "h_dpmsig");
504  }
505 
506  if (sqrts>38)
507  fillHistos(t[8], "dsm", cuts["nds1"], "h_ds", "h_dssel", "h_dssig");
508 
509  if (sqrts>45)
510  fillHistos(t[9], "lamcm", cuts["nlamc"], "h_lamc", "h_lamcsel", "h_lamcsig");
511 
512  // event counts
513 
514  double channelcount[10]={0,0,0,0,0,0,0,0,0,0};
515  double eff[10]={0,0,0,0,0,0,0,0,0,0};
516 
517  evcnt.clear();
518  for (i=0;i<3;++i)
519  {
520  channelcount[i]=cnt_events(t[i], cuts[ntpname[i]]+"&&tag", true);
521  cout <<i<< " "<<channelcount[i]<<endl;
522  }
523 
524  if (sqrts>24)
525  for (i=3;i<8;++i)
526  {
527  channelcount[i]=cnt_events(t[i], cuts[ntpname[i]]+"&&tag", true);
528  cout <<i<< " "<<channelcount[i]<<endl;
529  }
530  if (sqrts>38)
531  for (i=8;i<9;++i)
532  {
533  channelcount[i]=cnt_events(t[i], cuts[ntpname[i]]+"&&tag", true);
534  cout <<i<< " "<<channelcount[i]<<endl;
535  }
536  if (sqrts>45)
537  for (i=9;i<10;++i)
538  {
539  channelcount[i]=cnt_events(t[i], cuts[ntpname[i]]+"&&tag", true);
540  cout <<i<< " "<<channelcount[i]<<endl;
541  }
542 
543  double evCnt = evcnt.size();
544  for (i=0;i<10;++i) eff[i]=channelcount[i]/nev*100.;
545 
546  cout <<"all:"<<evCnt<<endl;
547 
548  // ** draw the plots
549  TLatex latex;
550  latex.SetTextSize(0.08);
551  char tmp[200];
552  double offset = 0;//0.1;
553  double posx = 0.24;
554  double posy = 0.79;
555 
556 // // *** plot all stuff and write out info
557 // c1->cd(1); h_pi0->Draw(); h_pi0sel->Draw("same"); h_pi0sig->Draw("same");
558 // c1->cd(2); h_ks->Draw(); h_kssel->Draw("same"); h_kssig->Draw("same");
559 
560  c2->cd(2); h_phi->Draw(); h_phisel->Draw("same"); h_phisig->Draw("same");
561  sprintf(tmp,"#epsilon = %4.1f%%",eff[0]);
562  latex.DrawLatexNDC(posx,posy,tmp);
563 
564  c2->cd(1); h_ee->Draw(); h_eesel->Draw("same"); h_eesig->Draw("same");
565  sprintf(tmp,"#epsilon = %4.1f%%",eff[2]);
566  latex.DrawLatexNDC(posx,posy,tmp);
567 
568  c2->cd(9); h_lam->Draw(); h_lamsel->Draw("same"); h_lamsig->Draw("same");
569  sprintf(tmp,"#epsilon = %4.1f%%",eff[1]);
570  latex.DrawLatexNDC(posx,posy,tmp);
571 
572  c2->cd(3); h_etac->Draw(); h_etacsel->Draw("same"); h_etacsig->Draw("same");
573  sprintf(tmp,"#epsilon = %4.1f%%",eff[5]);
574  latex.DrawLatexNDC(posx,posy,tmp);
575 
576  c2->cd(4); h_jpsi1->Draw(); h_jpsi1sel->Draw("same"); h_jpsi1sig->Draw("same");
577  sprintf(tmp,"#epsilon = %4.1f%%",eff[3]);
578  latex.DrawLatexNDC(posx,posy,tmp);
579 
580  c2->cd(5); h_jpsi2->Draw(); h_jpsi2sel->Draw("same"); h_jpsi2sig->Draw("same");
581  sprintf(tmp,"#epsilon = %4.1f%%",eff[4]);
582  latex.DrawLatexNDC(posx,posy,tmp);
583 
584  c2->cd(7); h_dpm->Draw(); h_dpmsel->Draw("same"); h_dpmsig->Draw("same");
585  sprintf(tmp,"#epsilon = %4.1f%%",eff[7]);
586  latex.DrawLatexNDC(posx,posy,tmp);
587 
588  c2->cd(6); h_d0->Draw(); h_d0sel->Draw("same"); h_d0sig->Draw("same");
589  sprintf(tmp,"#epsilon = %4.1f%%",eff[6]);
590  latex.DrawLatexNDC(posx,posy,tmp);
591 
592  c2->cd(8); h_ds->Draw(); h_dssel->Draw("same"); h_dssig->Draw("same");
593  sprintf(tmp,"#epsilon = %4.1f%%",eff[8]);
594  latex.DrawLatexNDC(posx,posy,tmp);
595 
596  c2->cd(10); h_lamc->Draw(); h_lamcsel->Draw("same"); h_lamcsig->Draw("same");
597  sprintf(tmp,"#epsilon = %4.1f%%",eff[9]);
598  latex.DrawLatexNDC(posx,posy,tmp);
599 
600  // total efficiency in first plot (phi)
601  c2->cd(5);
602  posy = h_jpsi2->GetMaximum()*1.12;
603  if (posy==0) posy=1.12;
604  sprintf(tmp,"(#epsilon_{t} = %4.1f%%)",double(evCnt)/double(nev)*100.);
605  latex.SetTextColor(kMagenta+2);
606  latex.DrawLatex(3.0,posy,tmp);
607 
608  c2->cd();
609 
610  //TString plot(fname(fname.Last('/')+1,fname.Last('.')-fname.Last('/')-1));
611 //fname(pos1+1,pos2-1));
612  TString plot(fname(fname.Last('/')+1,6));
613  if (!tagonly) plot+="_opt_supr";
614  //plot="figfull/"+plot+".pdf";
615 
616  //cout <<plot<<endl;
617 
618  c2->SaveAs("figfull_v2/full_"+plot+".gif");
619  c2->SaveAs("figfull_v2/full_"+plot+".pdf");
620 
621  std::ofstream ofs;
622 
623  ofs.open("full_supr1000.txt", std::ofstream::out | std::ofstream::app);
624  TString result=TString::Format("%s & %3.1f & %6.2f & %6.2f & %6.2f & %6.2f & %6.2f & %6.2f & %6.2f & %6.2f & %6.2f & %6.2f & %6.2f \\\\\n",
625  mode.Data(), sqrts/10., eff[2], eff[0], eff[5], eff[3], eff[4], eff[6], eff[7], eff[8], eff[1], eff[9], evCnt/nev*100);
626  //cout <<result;
627  ofs <<result;
628  ofs.close();
629  return 0;
630 }
631 
632 
i
Int_t i
Definition: run_full.C:25
plotname
TString plotname[10]
Definition: full_core_ntp.C:44
configHisto
void configHisto(TH1 *h)
Definition: full_core_ntp.C:76
fRand
TRandom3 fRand
Definition: full_core_ntp.C:28
init_cuts_1000
void init_cuts_1000(TString energy)
Definition: full_core_ntp.C:113
mctvar
std::map< TString, TString > mctvar
Definition: full_core_ntp.C:41
ev
int ev
Definition: GammaSpectraAnalysis.C:44
low
double low
Definition: full_core_ntp.C:34
fMAX
const int fMAX
Definition: full_core_ntp.C:32
c2
c2
Definition: plot_dirc.C:39
getPos
double getPos(TH1 *h)
Definition: full_core_ntp.C:71
ntpname
TString ntpname[10]
Definition: full_core_ntp.C:43
offset
TVector3 offset(2, 0, 0)
DbMass
double DbMass(TString name)
Definition: full_core_ntp.C:316
cut
double cut[MAX]
Definition: autocutx.C:36
mode
Int_t mode
Definition: autocutx.C:47
cuts
std::map< TString, TString > cuts
Definition: full_core_ntp.C:39
full_core_ntp
int full_core_ntp(TString fname, int nev=50, bool tagonly=false)
Definition: full_core_ntp.C:322
nbins
int nbins
Definition: full_core_ntp.C:33
cnt_events
int cnt_events(TChain *t, TString cut, bool globcnt=false)
Definition: full_core_ntp.C:285
out
TFile * out
Definition: reco_muo.C:20
name
TString name
Definition: createRootGeoFile.C:82
treeidx
std::map< int, int > treeidx
Definition: full_core_ntp.C:38
high
double high
Definition: full_core_ntp.C:35
evcnt
std::map< int, int > evcnt
Definition: full_core_ntp.C:37
modes
std::map< TString, TString > modes
Definition: full_core_ntp.C:40
fIni
TLorentzVector fIni
Definition: full_core_ntp.C:29
t
TTree * t
Definition: bump_analys.C:13
fillHistos
void fillHistos(TChain *t, TString var, TString cut, TString all, TString sel, TString sig)
Definition: full_core_ntp.C:309
setStyle
void setStyle()
Definition: full_core_ntp.C:47
init_modes
void init_modes()
Definition: full_core_ntp.C:86
init_cuts_95
void init_cuts_95(TString energy)
Definition: full_core_ntp.C:222
energy
Double_t energy
Definition: plot_dirc.C:15
Generated on Sun Aug 18 2019 03:01:12 for FairRoot/PandaRoot by   doxygen 1.8.8