AliPhysics  958ad07 (958ad07)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
ComputeAcceptance.C
Go to the documentation of this file.
1 #if !defined(__CINT__) || defined(__MAKECINT__)
2 #include <TDatabasePDG.h>
3 #include <TRandom3.h>
4 #include <TMath.h>
5 #include <TLorentzVector.h>
6 #include <TParticle.h>
7 #include <TParticle.h>
8 #include <TClonesArray.h>
9 #include <TH1D.h>
10 #include <TH2D.h>
11 #include <TF1.h>
12 #include <TFile.h>
13 #include <TCanvas.h>
14 #include <TLatex.h>
15 #include <TSystem.h>
16 #include <TStyle.h>
17 #include <TPythia6Decayer.h>
18 #include <TPaveStats.h>
19 #endif
20 
21 enum EDDecay{kD0Kpi,kDplusKpipi,kDstarD0pi,kDsKKpi,kLcpKpi,kLcK0Sp};
22 enum EFidY{kFixedY,kPtDepY};
23 enum EPtShape{kFlat,kFONLL8TeV,kFONLL8TeVfeeddown,kFONLL7TeV,kPythia7TeV,kFONLL5TeV,kFONLL13TeVprompt,kFONLL13TeVfeeddown,kPythia13TeVprompt,kPythia13TeVfeeddown};
24 
25 // Configuration
26 Int_t fDDecay=kD0Kpi;
27 Double_t fPtMinDau=0.1;
28 Double_t fEtaMaxDau=0.9;
29 Int_t fOptionYFiducial=kFixedY;
30 Double_t fYMaxFidAccCut=0.8;
31 Int_t fPtShape=kFONLL7TeV;
32 TString fDecayTableFileName="$ALICE_PHYSICS/PWGHF/vertexingHF/macros/decaytable_acc.dat";
33 Int_t fDebugLevel=0;
34 Int_t totTrials=1000000;
35 
36 
37 Bool_t CountKpi(TClonesArray *array, Int_t nentries, Int_t &nPions, Int_t &nKaons, Int_t &nPionsInAcc, Int_t &nKaonsInAcc);
38 Bool_t IsInFiducialAcceptance(Double_t pt, Double_t y);
39 Bool_t CountPKpi(TClonesArray *array, Int_t nentries, Int_t &nPions, Int_t &nKaons, Int_t &nProtons, Int_t &nPionsInAcc, Int_t &nKaonsInAcc, Int_t &nProtonsInAcc);
40 
41 
42 // Pt-shape histograms
43 TH1D* LoadFONLL13TeV_promptD0();
44 TH1D* LoadFONLL13TeV_promptDplus();
45 TH1D* LoadFONLL13TeV_promptDstar();
46 TH1D* LoadFONLL13TeV_feeddownD();
47 TH1D* LoadFONLL13TeV_feeddownDstar();
48 TH1D* LoadPYTHIA13TeV_promptD0();
49 TH1D* LoadPYTHIA13TeV_promptDplus();
50 TH1D* LoadPYTHIA13TeV_promptDstar();
51 TH1D* LoadPYTHIA13TeV_promptDs();
52 TH1D* LoadPYTHIA13TeV_feeddownD0();
53 TH1D* LoadPYTHIA13TeV_feeddownDplus();
54 TH1D* LoadPYTHIA13TeV_feeddownDstar();
55 TH1D* LoadPYTHIA13TeV_feeddownDs();
56 
57 
58 
59 void ComputeAcceptance(){
60  // main function
61 
62  gSystem->Load("liblhapdf.so"); // Parton density functions
63  gSystem->Load("libEGPythia6.so"); // TGenerator interface
64  gSystem->Load("libpythia6.so"); // Pythia
65  gSystem->Load("libAliPythia6"); // ALICE specific implementations
66 
67  TPythia6Decayer* pdec=TPythia6Decayer::Instance();
68  if(fDecayTableFileName.CompareTo("")!=0){
69  if(fDecayTableFileName.Contains("ALICE_PHYSICS")){
70  gSystem->Exec(Form("cp %s .",fDecayTableFileName.Data()));
71  fDecayTableFileName.ReplaceAll("$ALICE_PHYSICS/PWGHF/vertexingHF/macros/","./");
72  }
73  pdec->SetDecayTableFile(fDecayTableFileName.Data());
74  pdec->ReadDecayTable();
75  }
76  pdec->Init();
77 
78  Int_t pdgCode=0;
79  Int_t nPionDau=-1;
80  Int_t nProtonDau=-1;
81  Int_t nKaonDau=-1;
82  TString outFileName="Acceptance_Toy_";
83  if(fDDecay==kD0Kpi){
84  pdgCode=421;
85  nPionDau=1;
86  nKaonDau=1;
87  nProtonDau=0;
88  outFileName.Append("D0Kpi_");
89  }else if(fDDecay==kDplusKpipi){
90  pdgCode=411;
91  nPionDau=2;
92  nKaonDau=1;
93  nProtonDau=0;
94  outFileName.Append("DplusKpipi_");
95  }else if(fDDecay==kDstarD0pi){
96  pdgCode=413;
97  nPionDau=2;
98  nKaonDau=1;
99  nProtonDau=0;
100  outFileName.Append("DStarD0pi_");
101  }else if(fDDecay==kDsKKpi){
102  pdgCode=431;
103  nPionDau=1;
104  nKaonDau=2;
105  nProtonDau=0;
106  outFileName.Append("DsKKpi_");
107  }else if(fDDecay==kLcpKpi){
108  pdgCode=4122;
109  nPionDau=1;
110  nKaonDau=1;
111  nProtonDau=1;
112  outFileName.Append("LcpKpi_");
113  }else if(fDDecay==kLcK0Sp){
114  pdgCode=4122;
115  nPionDau=2;
116  nKaonDau=0;
117  nProtonDau=1;
118  outFileName.Append("LcK0Sp_");
119  }else{
120  printf("ERROR: Wrong decay selected\n");
121  return;
122  }
123  if(fOptionYFiducial==kFixedY) outFileName.Append(Form("yfid%02d_",(Int_t)(fYMaxFidAccCut*10)));
124  else outFileName.Append("yfidPtDep_");
125  outFileName.Append(Form("etaDau%02d_",(Int_t)(fEtaMaxDau*10)));
126  outFileName.Append(Form("ptDau%d_",(Int_t)(fPtMinDau*1000)));
127  TDatabasePDG* db=TDatabasePDG::Instance();
128  Float_t massD=db->GetParticle(pdgCode)->Mass();
129  TClonesArray *array = new TClonesArray("TParticle",100);
130 
131  TH2D* hPtVsYGen=new TH2D("hPtVsYGen","",400,0.,40.,20.,-1.,1.);
132  hPtVsYGen->GetXaxis()->SetTitle("p_{T} (GeV/c)");
133  hPtVsYGen->GetYaxis()->SetTitle("y");
134  TH2D* hPtVsYGenLimAcc=new TH2D("hPtVsYGenLimAcc","",400,0.,40.,20.,-1.,1.);
135  hPtVsYGenLimAcc->GetXaxis()->SetTitle("p_{T} (GeV/c)");
136  hPtVsYGenLimAcc->GetYaxis()->SetTitle("y");
137  TH2D* hPtVsYGenAcc=new TH2D("hPtVsYGenAcc","",400,0.,40.,20.,-1.,1.);
138  hPtVsYGenAcc->GetXaxis()->SetTitle("p_{T} (GeV/c)");
139  hPtVsYGenAcc->GetYaxis()->SetTitle("y");
140 
141  TF1* funcPt=0x0;
142  TH1D* histPt=0x0;
143  if(fPtShape==kFONLL8TeV){
144  funcPt=new TF1("fFONLL","[0]*x/TMath::Power((1+TMath::Power(x/[1],[3])),[2])",0.,40.);
145  funcPt->SetParameters(0.518046,3.01138,3.38914,1.75899); // Prompt
146  outFileName.Append("FONLL8ptshape.root");
147  }else if(fPtShape==kFONLL8TeVfeeddown){
148  funcPt=new TF1("fFONLL","[0]*x/TMath::Power((1+TMath::Power(x/[1],[3])),[2])",0.,40.);
149  funcPt->SetParameters(0.398252, 3.9603, 3.915, 1.51853); // FeedDown
150  outFileName.Append("FONLL8ptshapeFeedDown.root");
151  }else if(fPtShape==kFONLL7TeV){
152  funcPt=new TF1("fFONLL","[0]*x/TMath::Power((1+TMath::Power(x/[1],[3])),[2])",0.,40.);
153  funcPt->SetParameters(0.322643,2.96275,2.30301,2.5);
154  outFileName.Append("FONLL7ptshape.root");
155  }else if(fPtShape==kFONLL5TeV){
156  funcPt=new TF1("fFONLL","[0]*x/TMath::Power((1+TMath::Power(x/[1],[3])),[2])",0.,40.);
157  funcPt->SetParameters(0.302879,2.9750,3.68139,1.68855);
158  outFileName.Append("FONLL5ptshape.root");
159  }else if(fPtShape==kPythia7TeV){
160  funcPt=new TF1("fFONLL","[0]*x/TMath::Power((1+TMath::Power(x/[1],[3])),[2])",0.,40.);
161  funcPt->SetParameters(0.322643,1.94635,1.40463,2.5);
162  outFileName.Append("PYTHIA7ptshape.root");
163  }else if(fPtShape==kFONLL13TeVprompt){
164  if(fDDecay==kDplusKpipi){
165  histPt = LoadFONLL13TeV_promptDplus();
166  outFileName.Append("promptDplus");
167  }else if (fDDecay==kDstarD0pi){
168  histPt = LoadFONLL13TeV_promptDstar();
169  outFileName.Append("promptDstar");
170  }else{
171  histPt = LoadFONLL13TeV_promptD0();
172  outFileName.Append("promptD0");
173  }
174  outFileName.Append("FONLL13ptshape.root");
175 
176  }else if (fPtShape==kFONLL13TeVfeeddown){
177  if (fDDecay==kDstarD0pi){
178  histPt = LoadFONLL13TeV_feeddownDstar();
179  outFileName.Append("feeddownDstar");
180  }else{
181  histPt = LoadFONLL13TeV_feeddownD();
182  outFileName.Append("feeddownD");
183  }
184  outFileName.Append("FONLL13ptshape.root");
185  }else if(fPtShape==kPythia13TeVprompt){
186  if(fDDecay==kDplusKpipi){
187  histPt = LoadPYTHIA13TeV_promptDplus();
188  outFileName.Append("promptDplus");
189  }else if (fDDecay==kDstarD0pi){
190  histPt = LoadPYTHIA13TeV_promptDstar();
191  outFileName.Append("promptDstar");
192  }else if (fDDecay==kDsKKpi){
193  histPt = LoadPYTHIA13TeV_promptDs();
194  outFileName.Append("promptDs");
195  }else{
196  histPt = LoadPYTHIA13TeV_promptD0();
197  outFileName.Append("promptD0");
198  }
199  outFileName.Append("PYTHIA13ptshape.root");
200  }else if (fPtShape==kPythia13TeVfeeddown){
201  if(fDDecay==kDplusKpipi){
202  histPt = LoadPYTHIA13TeV_feeddownDplus();
203  outFileName.Append("feeddownDplus");
204  }else if (fDDecay==kDstarD0pi){
205  histPt = LoadPYTHIA13TeV_feeddownDstar();
206  outFileName.Append("feeddownDstar");
207  }else if (fDDecay==kDsKKpi){
208  histPt = LoadPYTHIA13TeV_feeddownDs();
209  outFileName.Append("feeddownDs");
210  }else{
211  histPt = LoadPYTHIA13TeV_feeddownD0();
212  outFileName.Append("feeddownD0");
213  }
214  outFileName.Append("PYTHIA13ptshape.root");
215  }else{
216  funcPt=new TF1("fFlat","pol0",0.,40.);
217  funcPt->SetParameter(0,1.);
218  outFileName.Append("flatpt.root");
219  }
220 
221  if (funcPt) funcPt->SetNpx(10000);
222 
223  TRandom3* gener=new TRandom3(0);
224  TLorentzVector* vec=new TLorentzVector();
225 
226 
227  for(Int_t itry=0; itry<totTrials; itry++){
228  if(itry%10000==0) printf("Event %d\n",itry);
229  Float_t ptD = funcPt ? funcPt->GetRandom() : histPt->GetRandom();
230  Float_t phiD=gener->Rndm()*2*TMath::Pi();
231  Float_t yD=gener->Rndm()*2.-1.; // flat in -1<y<1
232  Float_t px=ptD*TMath::Cos(phiD);
233  Float_t py=ptD*TMath::Sin(phiD);
234  Float_t mt=TMath::Sqrt(massD*massD+ptD*ptD);
235  Float_t pz=mt*TMath::SinH(yD);
236  Float_t E=TMath::Sqrt(massD*massD+px*px+py*py+pz*pz);
237 
238  // TLorentzVector* vec=new TLorentzVector(px,py,pz,E);
239  vec->SetPxPyPzE(px,py,pz,E);
240  pdec->Decay(pdgCode,vec);
241  array->Clear();
242  Int_t nentries = pdec->ImportParticles(array);
243  TParticle* dmes=(TParticle*)array->At(0);
244  Int_t nDaughters=dmes->GetNDaughters();
245  if(fDDecay==kD0Kpi && nDaughters!=2) continue;
246  if(fDDecay==kLcK0Sp && nentries>6) continue;
247  Int_t nPionsInAcc=0;
248  Int_t nProtonsInAcc=0;
249  Int_t nKaonsInAcc=0;
250  Int_t nPions=0;
251  Int_t nProtons=0;
252  Int_t nKaons=0;
253  Bool_t isOk=CountPKpi(array,nentries,nPions,nKaons,nProtons,nPionsInAcc,nKaonsInAcc,nProtonsInAcc);
254 
255  if(isOk){
256  if(nPions==nPionDau && nKaons==nKaonDau && nProtons==nProtonDau){
257  hPtVsYGen->Fill(ptD,yD);
258  if(TMath::Abs(yD)<0.5){
259  hPtVsYGenLimAcc->Fill(ptD,yD);
260  }
261  if(IsInFiducialAcceptance(ptD,yD)){
262  if(nPionsInAcc==nPionDau && nKaonsInAcc==nKaonDau && nProtonsInAcc==nProtonDau){
263  hPtVsYGenAcc->Fill(ptD,yD);
264  }
265  }
266  }
267  }
268  // delete vec;
269  }
270 
271  TH1D* hPtGenAcc=(TH1D*)hPtVsYGenAcc->ProjectionX("hPtGenAcc");
272  hPtGenAcc->GetYaxis()->SetTitle("Entries");
273  TH1D* hPtGenLimAcc=(TH1D*)hPtVsYGenLimAcc->ProjectionX("hPtGenLimAcc");
274  hPtGenLimAcc->GetYaxis()->SetTitle("Entries");
275  hPtGenAcc->Sumw2();
276  hPtGenLimAcc->Sumw2();
277  TH1D* hAccVsPt=(TH1D*)hPtGenAcc->Clone("hAccVsPt");
278  hAccVsPt->Divide(hPtGenAcc,hPtGenLimAcc,1,1,"B");
279  hAccVsPt->GetYaxis()->SetTitle("Acceptance");
280  hAccVsPt->SetStats(0);
281 
282  TCanvas* c2d=new TCanvas("c2d","Pt vs y",1200,600);
283  c2d->Divide(3,1);
284  c2d->cd(1);
285  hPtVsYGen->Draw("colz");
286  c2d->cd(2);
287  hPtVsYGenLimAcc->Draw("colz");
288  c2d->cd(3);
289  hPtVsYGenAcc->Draw("colz");
290 
291  TCanvas* c1d=new TCanvas("c1d","Acceptance",1200,600);
292  c1d->Divide(2,1);
293  c1d->cd(1);
294  hPtGenLimAcc->Draw("");
295  Double_t ymax=1.2*TMath::Max(hPtGenLimAcc->GetMaximum(),hPtGenAcc->GetMaximum());
296  hPtGenLimAcc->SetMaximum(ymax);
297  gPad->Update();
298  TPaveStats *st1=(TPaveStats*)hPtGenLimAcc->GetListOfFunctions()->FindObject("stats");
299  st1->SetY1NDC(0.71);
300  st1->SetY2NDC(0.9);
301  hPtGenAcc->SetLineColor(kRed+1);
302  hPtGenAcc->Draw("sames");
303  gPad->Update();
304  TPaveStats *st2=(TPaveStats*)hPtGenAcc->GetListOfFunctions()->FindObject("stats");
305  st2->SetY1NDC(0.51);
306  st2->SetY2NDC(0.7);
307  st2->SetTextColor(hPtGenAcc->GetLineColor());
308  gPad->Modified();
309  c1d->cd(2);
310  hAccVsPt->Draw();
311 
312 
313  TFile* outfil=new TFile(outFileName.Data(),"recreate");
314  hPtVsYGen->Write();
315  hPtVsYGenLimAcc->Write();
316  hPtVsYGenAcc->Write();
317  hPtGenLimAcc->Write();
318  hPtGenAcc->Write();
319  hAccVsPt->Write();
320  outfil->Close();
321 
322 }
323 
324 //___________________________________________________
325 Bool_t IsInFiducialAcceptance(Double_t pt, Double_t y){
326  // check fiducial acceptance
327 
328  if(fOptionYFiducial==kFixedY){
329  if(TMath::Abs(y) > fYMaxFidAccCut) return kFALSE;
330  else return kTRUE;
331  }
332 
333  if(pt > 5.) {
334  if (TMath::Abs(y) > 0.8) return kFALSE;
335  } else {
336  Double_t maxFiducialY = -0.2/15*pt*pt+1.9/15*pt+0.5;
337  Double_t minFiducialY = 0.2/15*pt*pt-1.9/15*pt-0.5;
338  if (y < minFiducialY || y > maxFiducialY) return kFALSE;
339  }
340 
341  return kTRUE;
342 }
343 //___________________________________________________
344 Bool_t CountKpi(TClonesArray *array, Int_t nentries, Int_t &nPions, Int_t &nKaons, Int_t &nPionsInAcc, Int_t &nKaonsInAcc){
345  // count K and pi in Acc
346 
347  TParticle* dmes=(TParticle*)array->At(0);
348  Double_t sumPx=0;
349  Double_t sumPy=0;
350  Double_t sumPz=0;
351  for(int j=0; j<nentries; j++){
352  TParticle * o = (TParticle*)array->At(j);
353  Int_t pdgdau=TMath::Abs(o->GetPdgCode());
354  if(fDebugLevel>0) printf("%d ",pdgdau);
355  if(pdgdau==130) {
356  if(fDebugLevel>0) printf("K0 dacaying into K0L\n");
357  return kFALSE;
358  }
359  Float_t ptdau=TMath::Sqrt(o->Px()*o->Px()+o->Py()*o->Py());
360  Float_t etadau=o->Eta();
361  if(pdgdau==211){
362  nPions++;
363  sumPx+=o->Px();
364  sumPy+=o->Py();
365  sumPz+=o->Pz();
366  }
367  if(pdgdau==321){
368  nKaons++;
369  sumPx+=o->Px();
370  sumPy+=o->Py();
371  sumPz+=o->Pz();
372  }
373  if(TMath::Abs(etadau)<fEtaMaxDau && ptdau>fPtMinDau){
374  if(pdgdau==211) nPionsInAcc++;
375  if(pdgdau==321) nKaonsInAcc++;
376  }
377  }
378  if(fDebugLevel>0) printf("\n");
379  if(TMath::Abs(sumPx-dmes->Px())>0.001 ||
380  TMath::Abs(sumPy-dmes->Py())>0.001 ||
381  TMath::Abs(sumPz-dmes->Pz())>0.001){
382  printf("Momentum conservation violation\n");
383  return kFALSE;
384  }
385  return kTRUE;
386 }
387 
388 
389 //___________________________________________________
390 Bool_t CountPKpi(TClonesArray *array, Int_t nentries, Int_t &nPions, Int_t &nKaons, Int_t &nProtons, Int_t &nPionsInAcc, Int_t &nKaonsInAcc, Int_t &nProtonsInAcc){
391  // count K and pi in Acc
392 
393  TParticle* dmes=(TParticle*)array->At(0);
394  Double_t sumPx=0;
395  Double_t sumPy=0;
396  Double_t sumPz=0;
397 
398  for(int j=0; j<nentries; j++){
399  TParticle * o = (TParticle*)array->At(j);
400  Int_t pdgdau=TMath::Abs(o->GetPdgCode());
401  if(fDebugLevel>0) printf("%d ",pdgdau);
402  if(pdgdau==130) {
403  if(fDebugLevel>0) printf("K0 dacaying into K0L\n");
404  return kFALSE;
405  }
406  Float_t ptdau=TMath::Sqrt(o->Px()*o->Px()+o->Py()*o->Py());
407  Float_t etadau=o->Eta();
408  if(pdgdau==211){
409  nPions++;
410  sumPx+=o->Px();
411  sumPy+=o->Py();
412  sumPz+=o->Pz();
413  }
414  if(pdgdau==321){
415  nKaons++;
416  sumPx+=o->Px();
417  sumPy+=o->Py();
418  sumPz+=o->Pz();
419  }
420  if(pdgdau==2212){
421  nProtons++;
422  sumPx+=o->Px();
423  sumPy+=o->Py();
424  sumPz+=o->Pz();
425  }
426  if(TMath::Abs(etadau)<fEtaMaxDau && ptdau>fPtMinDau){
427  if(pdgdau==211) nPionsInAcc++;
428  if(pdgdau==321) nKaonsInAcc++;
429  if(pdgdau==2212) nProtonsInAcc++;
430  }
431  }
432  if(fDebugLevel>0) printf("\n");
433  if(TMath::Abs(sumPx-dmes->Px())>0.001 ||
434  TMath::Abs(sumPy-dmes->Py())>0.001 ||
435  TMath::Abs(sumPz-dmes->Pz())>0.001){
436  printf("Momentum conservation violation\n");
437  return kFALSE;
438  }
439  return kTRUE;
440 }
441 
442 
443 
444 //___________________________________________________
445 TH1D* LoadFONLL13TeV_promptD0()
446 {
447  TH1D *hFONLL13 = new TH1D("hFONLL13TeV_D0", "", 80, 0., 40.);
448  Float_t val[80] = {
449  1.4686e+08, 3.9542e+08, 5.4901e+08, 5.2166e+08, 4.1083e+08, 2.9968e+08, 2.1299e+08, 1.5057e+08, 1.0701e+08, 7.6919e+07,
450  5.6121e+07, 4.1546e+07, 3.1184e+07, 2.3715e+07, 1.8253e+07, 1.4206e+07, 1.1175e+07, 8.8774e+06, 7.1169e+06, 5.7544e+06,
451  4.6899e+06, 3.8509e+06, 3.1841e+06, 2.6499e+06, 2.2189e+06, 1.8687e+06, 1.5823e+06, 1.3467e+06, 1.1516e+06, 9.8933e+05,
452  8.5356e+05, 7.3942e+05, 6.4302e+05, 5.6122e+05, 4.9154e+05, 4.3193e+05, 3.8074e+05, 3.3662e+05, 2.9846e+05, 2.6535e+05,
453  2.3653e+05, 2.1136e+05, 1.8932e+05, 1.6997e+05, 1.5293e+05, 1.3789e+05, 1.2458e+05, 1.1277e+05, 1.0227e+05, 9.2913e+04,
454  8.4560e+04, 7.7087e+04, 7.0388e+04, 6.4371e+04, 5.8956e+04, 5.4074e+04, 4.9667e+04, 4.5680e+04, 4.2068e+04, 3.8791e+04,
455  3.5813e+04, 3.3103e+04, 3.0633e+04, 2.8380e+04, 2.6321e+04, 2.4436e+04, 2.2710e+04, 2.1127e+04, 1.9673e+04, 1.8336e+04,
456  1.7106e+04, 1.5972e+04, 1.4926e+04, 1.3961e+04, 1.3069e+04, 1.2243e+04, 1.1479e+04, 1.0771e+04, 1.0113e+04, 9.5031e+03
457  };
458  for (Int_t ibin=0; ibin<80; ++ibin) hFONLL13->SetBinContent(ibin+1, val[ibin]);
459 
460  return hFONLL13;
461 }
462 
463 
464 
465 //___________________________________________________
466 TH1D* LoadFONLL13TeV_promptDplus()
467 {
468  TH1D *hFONLL13 = new TH1D("hFONLL13TeV_Dplus", "", 80, 0., 40.);
469  Float_t val[80] = {
470  1.5242e+08, 3.9396e+08, 5.3767e+08, 5.1263e+08, 4.0706e+08, 2.9926e+08, 2.1414e+08, 1.5230e+08, 1.0879e+08, 7.8507e+07,
471  5.7460e+07, 4.2651e+07, 3.2093e+07, 2.4462e+07, 1.8865e+07, 1.4708e+07, 1.1588e+07, 9.2185e+06, 7.3998e+06, 5.9900e+06,
472  4.8871e+06, 4.0167e+06, 3.3240e+06, 2.7686e+06, 2.3200e+06, 1.9552e+06, 1.6566e+06, 1.4107e+06, 1.2070e+06, 1.0374e+06,
473  8.9550e+05, 7.7610e+05, 6.7519e+05, 5.8953e+05, 5.1652e+05, 4.5404e+05, 4.0037e+05, 3.5409e+05, 3.1405e+05, 2.7929e+05,
474  2.4902e+05, 2.2258e+05, 1.9943e+05, 1.7909e+05, 1.6117e+05, 1.4535e+05, 1.3135e+05, 1.1892e+05, 1.0787e+05, 9.8024e+04,
475  8.9230e+04, 8.1359e+04, 7.4302e+04, 6.7963e+04, 6.2256e+04, 5.7112e+04, 5.2465e+04, 4.8261e+04, 4.4452e+04, 4.0996e+04,
476  3.7854e+04, 3.4995e+04, 3.2390e+04, 3.0012e+04, 2.7838e+04, 2.5849e+04, 2.4027e+04, 2.2355e+04, 2.0820e+04, 1.9408e+04,
477  1.8108e+04, 1.6910e+04, 1.5805e+04, 1.4785e+04, 1.3842e+04, 1.2969e+04, 1.2161e+04, 1.1411e+04, 1.0716e+04, 1.0070e+04
478  };
479  for (Int_t ibin=0; ibin<80; ++ibin) hFONLL13->SetBinContent(ibin+1, val[ibin]);
480 
481  return hFONLL13;
482 }
483 
484 
485 
486 //___________________________________________________
487 TH1D* LoadFONLL13TeV_promptDstar()
488 {
489  TH1D *hFONLL13 = new TH1D("hFONLL13TeV_Dstar", "", 80, 0., 40.);
490  Float_t val[80] = {
491  1.2433e+08, 3.4512e+08, 5.0662e+08, 5.1020e+08, 4.2016e+08, 3.1661e+08, 2.3064e+08, 1.6632e+08, 1.2007e+08, 8.7340e+07,
492  6.4329e+07, 4.8000e+07, 3.6287e+07, 2.7772e+07, 2.1492e+07, 1.6808e+07, 1.3278e+07, 1.0589e+07, 8.5178e+06, 6.9084e+06,
493  5.6462e+06, 4.6479e+06, 3.8519e+06, 3.2124e+06, 2.6951e+06, 2.2738e+06, 1.9285e+06, 1.6437e+06, 1.4076e+06, 1.2108e+06,
494  1.0459e+06, 9.0711e+05, 7.8968e+05, 6.8993e+05, 6.0484e+05, 5.3197e+05, 4.6933e+05, 4.1529e+05, 3.6850e+05, 3.2786e+05,
495  2.9246e+05, 2.6152e+05, 2.3441e+05, 2.1058e+05, 1.8958e+05, 1.7104e+05, 1.5461e+05, 1.4003e+05, 1.2706e+05, 1.1550e+05,
496  1.0517e+05, 9.5920e+04, 8.7626e+04, 8.0172e+04, 7.3461e+04, 6.7409e+04, 6.1940e+04, 5.6992e+04, 5.2508e+04, 4.8437e+04,
497  4.4736e+04, 4.1367e+04, 3.8296e+04, 3.5492e+04, 3.2929e+04, 3.0583e+04, 2.8433e+04, 2.6460e+04, 2.4648e+04, 2.2981e+04,
498  2.1446e+04, 2.0032e+04, 1.8727e+04, 1.7522e+04, 1.6407e+04, 1.5376e+04, 1.4420e+04, 1.3535e+04, 1.2713e+04, 1.1950e+04
499  };
500  for (Int_t ibin=0; ibin<80; ++ibin) hFONLL13->SetBinContent(ibin+1, val[ibin]);
501 
502  return hFONLL13;
503 }
504 
505 
506 
507 //___________________________________________________
508 TH1D* LoadFONLL13TeV_feeddownD()
509 { // B->D with B.R=1
510  TH1D *hFONLL13 = new TH1D("FONLL13TeV_feeddownD", "", 80, 0., 40.);
511  Float_t val[80] = {
512  1.0310e+07, 2.6790e+07, 3.4480e+07, 3.4430e+07, 3.0200e+07, 2.4740e+07, 1.9600e+07, 1.5300e+07, 1.1880e+07, 9.2260e+06,
513  7.1950e+06, 5.6470e+06, 4.4650e+06, 3.5580e+06, 2.8570e+06, 2.3120e+06, 1.8860e+06, 1.5490e+06, 1.2810e+06, 1.0660e+06,
514  8.9260e+05, 7.5170e+05, 6.3650e+05, 5.4160e+05, 4.6310e+05, 3.9770e+05, 3.4300e+05, 2.9700e+05, 2.5820e+05, 2.2520e+05,
515  1.9710e+05, 1.7310e+05, 1.5240e+05, 1.3460e+05, 1.1920e+05, 1.0590e+05, 9.4290e+04, 8.4160e+04, 7.5290e+04, 6.7500e+04,
516  6.0650e+04, 5.4610e+04, 4.9270e+04, 4.4530e+04, 4.0320e+04, 3.6580e+04, 3.3240e+04, 3.0250e+04, 2.7570e+04, 2.5170e+04,
517  2.3020e+04, 2.1070e+04, 1.9320e+04, 1.7740e+04, 1.6310e+04, 1.5010e+04, 1.3830e+04, 1.2760e+04, 1.1790e+04, 1.0900e+04,
518  1.0090e+04, 9.3520e+03, 8.6760e+03, 8.0560e+03, 7.4890e+03, 6.9670e+03, 6.4880e+03, 6.0470e+03, 5.6410e+03, 5.2670e+03,
519  4.9220e+03, 4.6030e+03, 4.3080e+03, 4.0350e+03, 3.7830e+03, 3.5480e+03, 3.3310e+03, 3.1290e+03, 2.9410e+03, 2.7670e+03
520  };
521  for (Int_t ibin=0; ibin<80; ++ibin) hFONLL13->SetBinContent(ibin+1, val[ibin]);
522 
523  return hFONLL13;
524 }
525 
526 
527 
528 //___________________________________________________
529 TH1D* LoadFONLL13TeV_feeddownDstar()
530 { // B->D* with B.R=1
531  TH1D *hFONLL13 = new TH1D("FONLL13TeV_feeddownDstar", "", 80, 0., 40.);
532  Float_t val[80] = {
533  9.5260e+06, 2.5070e+07, 3.2890e+07, 3.3540e+07, 3.0010e+07, 2.5000e+07, 2.0090e+07, 1.5860e+07, 1.2430e+07, 9.7280e+06,
534  7.6360e+06, 6.0240e+06, 4.7820e+06, 3.8240e+06, 3.0800e+06, 2.5000e+06, 2.0430e+06, 1.6810e+06, 1.3930e+06, 1.1610e+06,
535  9.7390e+05, 8.2120e+05, 6.9610e+05, 5.9300e+05, 5.0750e+05, 4.3620e+05, 3.7650e+05, 3.2620e+05, 2.8370e+05, 2.4760e+05,
536  2.1680e+05, 1.9050e+05, 1.6780e+05, 1.4830e+05, 1.3140e+05, 1.1680e+05, 1.0400e+05, 9.2870e+04, 8.3100e+04, 7.4530e+04,
537  6.6990e+04, 6.0330e+04, 5.4440e+04, 4.9220e+04, 4.4580e+04, 4.0440e+04, 3.6760e+04, 3.3460e+04, 3.0510e+04, 2.7860e+04,
538  2.5470e+04, 2.3330e+04, 2.1390e+04, 1.9640e+04, 1.8060e+04, 1.6630e+04, 1.5320e+04, 1.4140e+04, 1.3060e+04, 1.2080e+04,
539  1.1190e+04, 1.0370e+04, 9.6190e+03, 8.9330e+03, 8.3050e+03, 7.7270e+03, 7.1970e+03, 6.7080e+03, 6.2580e+03, 5.8440e+03,
540  5.4610e+03, 5.1080e+03, 4.7810e+03, 4.4790e+03, 4.1980e+03, 3.9390e+03, 3.6970e+03, 3.4740e+03, 3.2650e+03, 3.0720e+03
541  };
542  for (Int_t ibin=0; ibin<80; ++ibin) hFONLL13->SetBinContent(ibin+1, val[ibin]);
543 
544  return hFONLL13;
545 }
546 
547 
548 
549 //___________________________________________________
550 TH1D* LoadPYTHIA13TeV_promptD0()
551 {
552  TH1D *hPYTHIA13 = new TH1D("hPYTHIA13TeV_D0", "", 40, 0., 40.);
553  Float_t val[40] = {
554  2617743, 3763836, 2235903, 1140807, 587707, 317881, 180756, 107762, 67299,42882,
555  28550, 19690, 13817, 9754, 7205, 5299, 4100, 3177, 2386,1822,
556  1484, 1248, 933, 798, 639, 534, 435, 405, 304,260,
557  210, 196, 161, 143, 139, 99, 89, 80, 66,63
558  };
559  for (Int_t ibin=0; ibin<40; ++ibin) hPYTHIA13->SetBinContent(ibin+1, val[ibin]);
560 
561  return hPYTHIA13;
562 }
563 
564 
565 
566 //___________________________________________________
567 TH1D* LoadPYTHIA13TeV_promptDplus()
568 {
569  TH1D *hPYTHIA13 = new TH1D("hPYTHIA13TeV_Dplus", "", 40, 0., 40.);
570  Float_t val[40] = {
571  1192379, 1736283, 1047485, 540499, 281523, 152153, 87150, 52200, 32709, 20991,
572  13986, 9611, 6832, 4705, 3484, 2650, 2001, 1518, 1197, 924,
573  735, 576, 497, 402, 292, 261, 214, 173, 169, 150,
574  106, 102, 95, 70, 67, 47, 53, 46, 35, 35
575  };
576  for (Int_t ibin=0; ibin<40; ++ibin) hPYTHIA13->SetBinContent(ibin+1, val[ibin]);
577 
578  return hPYTHIA13;
579 }
580 
581 
582 
583 //___________________________________________________
584 TH1D* LoadPYTHIA13TeV_promptDstar()
585 {
586  TH1D *hPYTHIA13 = new TH1D("hPYTHIA13TeV_Dstar", "", 40, 0., 40.);
587  Float_t val[40] = {
588  922985, 1419749, 904587, 484407, 257862, 142386, 82177, 49839, 31103, 20422,
589  13493, 9272, 6549, 4746, 3344, 2575, 1959, 1493, 1161, 935,
590  672, 585, 456, 415, 328, 260, 199, 180, 161, 134,
591  109, 97, 103, 72, 63, 56, 40, 46, 38, 27
592  };
593  for (Int_t ibin=0; ibin<40; ++ibin) hPYTHIA13->SetBinContent(ibin+1, val[ibin]);
594 
595  return hPYTHIA13;
596 }
597 
598 
599 
600 //___________________________________________________
601 TH1D* LoadPYTHIA13TeV_promptDs()
602 {
603  TH1D *hPYTHIA13 = new TH1D("hPYTHIA13TeV_Ds", "", 40, 0., 40.);
604  Float_t val[40] = {
605  346381, 519143, 320488, 167944, 87484, 47325, 26932, 16376, 10058, 6527,
606  4347, 3041, 2112, 1521, 1069, 849, 609, 468, 359, 275,
607  226, 197, 152, 124, 111, 90, 63, 60, 52, 46,
608  39, 35, 26, 29, 16, 20, 10, 13, 18, 9
609  };
610  for (Int_t ibin=0; ibin<40; ++ibin) hPYTHIA13->SetBinContent(ibin+1, val[ibin]);
611 
612  return hPYTHIA13;
613 }
614 
615 
616 
617 //___________________________________________________
618 TH1D* LoadPYTHIA13TeV_feeddownD0()
619 {
620  TH1D *hPYTHIA13 = new TH1D("hPYTHIA13_feeddownD0", "", 40, 0., 40.);
621  Float_t val[40] = {
622  1814030, 2939561, 2059721, 1197596, 684753, 397435, 237432, 146442, 93029, 60222,
623  40940, 28861, 19992, 14344, 10602, 7742, 5781, 4395, 3421, 2732,
624  2216, 1763, 1395, 1208, 888, 735, 599, 509, 418, 374,
625  302, 281, 242, 200, 174, 158, 133, 104, 98, 88
626  };
627  for (Int_t ibin=0; ibin<40; ++ibin) hPYTHIA13->SetBinContent(ibin+1, val[ibin]);
628 
629  return hPYTHIA13;
630 }
631 
632 
633 
634 //___________________________________________________
635 TH1D* LoadPYTHIA13TeV_feeddownDplus()
636 {
637  TH1D *hPYTHIA13 = new TH1D("hPYTHIA13_feeddownDplus", "", 40, 0., 40.);
638  Float_t val[40] = {
639  779389, 1268000, 899979, 527442, 302782, 177610, 105677, 66601, 42008, 27618,
640  18802, 13024, 8975, 6530, 4773, 3495, 2708, 2020, 1634, 1212,
641  1017, 809, 665, 507, 389, 355, 285, 250, 218, 180,
642  143, 137, 106, 85, 71, 64, 57, 43, 50, 47
643  };
644  for (Int_t ibin=0; ibin<40; ++ibin) hPYTHIA13->SetBinContent(ibin+1, val[ibin]);
645 
646  return hPYTHIA13;
647 }
648 
649 
650 
651 //___________________________________________________
652 TH1D* LoadPYTHIA13TeV_feeddownDstar()
653 {
654  TH1D *hPYTHIA13 = new TH1D("hPYTHIA13_feeddownDstar", "", 40, 0., 40.);
655  Float_t val[40] = {
656  665735, 1144214, 856647, 523704, 309061, 183569, 111994, 70187, 45266, 29601,
657  20079, 14042, 10338, 7184, 5234, 4002, 2935, 2301, 1788, 1349,
658  1088, 895, 748, 591, 480, 378, 302, 263, 216, 182,
659  172, 122, 133, 100, 92, 87, 65, 50, 42, 45
660  };
661  for (Int_t ibin=0; ibin<40; ++ibin) hPYTHIA13->SetBinContent(ibin+1, val[ibin]);
662 
663  return hPYTHIA13;
664 }
665 
666 
667 
668 //___________________________________________________
669 TH1D* LoadPYTHIA13TeV_feeddownDs()
670 {
671  TH1D *hPYTHIA13 = new TH1D("hPYTHIA13_feeddownDs", "", 40, 0., 40.);
672  Float_t val[40] = {
673  377925, 686435, 518802, 318819, 189396, 113651, 68963, 43048, 27776, 18177,
674  12356, 8643, 6100, 4512, 3170, 2479, 1835, 1424, 1086, 825,
675  659, 490, 465, 361, 308, 246, 181, 170, 117, 111,
676  112, 91, 69, 66, 54, 57, 39, 33, 31, 25
677  };
678  for (Int_t ibin=0; ibin<40; ++ibin) hPYTHIA13->SetBinContent(ibin+1, val[ibin]);
679 
680  return hPYTHIA13;
681 }
fPtShape
Int_t fPtShape
Definition: ComputeAcceptance.C:31
fEtaMaxDau
Double_t fEtaMaxDau
Definition: ComputeAcceptance.C:28
CountKpi
Bool_t CountKpi(TClonesArray *array, Int_t nentries, Int_t &nPions, Int_t &nKaons, Int_t &nPionsInAcc, Int_t &nKaonsInAcc)
Definition: ComputeAcceptance.C:344
kDstarD0pi
Definition: ComputeAcceptance.C:21
ymax
const Double_t ymax
Definition: AddTaskCFDStar.C:7
Double_t
double Double_t
Definition: External.C:58
kFONLL8TeV
Definition: ComputeAcceptance.C:23
LoadPYTHIA13TeV_promptD0
TH1D * LoadPYTHIA13TeV_promptD0()
Definition: ComputeAcceptance.C:550
gSystem
TSystem * gSystem
fOptionYFiducial
Int_t fOptionYFiducial
Definition: ComputeAcceptance.C:29
LoadFONLL13TeV_feeddownDstar
TH1D * LoadFONLL13TeV_feeddownDstar()
Definition: ComputeAcceptance.C:529
kPtDepY
Definition: ComputeAcceptance.C:22
ComputeAcceptance
void ComputeAcceptance()
Definition: ComputeAcceptance.C:59
kPythia13TeVfeeddown
Definition: ComputeAcceptance.C:23
fPtMinDau
Double_t fPtMinDau
Definition: ComputeAcceptance.C:27
EDDecay
EDDecay
Definition: ComputeAcceptance.C:21
kFONLL8TeVfeeddown
Definition: ComputeAcceptance.C:23
massD
Double_t massD
Definition: ProjectCombinHFAndFit.C:50
Int_t
int Int_t
Definition: External.C:63
kFONLL7TeV
Definition: ComputeAcceptance.C:23
Float_t
float Float_t
Definition: External.C:68
kLcK0Sp
Definition: ComputeAcceptance.C:21
kFixedY
Definition: ComputeAcceptance.C:22
kLcpKpi
Definition: ComputeAcceptance.C:21
TH2D
Definition: External.C:228
TH1D
Definition: External.C:212
totTrials
Int_t totTrials
Definition: ComputeAcceptance.C:34
LoadPYTHIA13TeV_feeddownDs
TH1D * LoadPYTHIA13TeV_feeddownDs()
Definition: ComputeAcceptance.C:669
LoadPYTHIA13TeV_feeddownDplus
TH1D * LoadPYTHIA13TeV_feeddownDplus()
Definition: ComputeAcceptance.C:635
LoadPYTHIA13TeV_feeddownD0
TH1D * LoadPYTHIA13TeV_feeddownD0()
Definition: ComputeAcceptance.C:618
fDebugLevel
Int_t fDebugLevel
Definition: ComputeAcceptance.C:33
LoadPYTHIA13TeV_promptDs
TH1D * LoadPYTHIA13TeV_promptDs()
Definition: ComputeAcceptance.C:601
kFONLL5TeV
Definition: ComputeAcceptance.C:23
EPtShape
EPtShape
Definition: ComputeAcceptance.C:23
kPythia13TeVprompt
Definition: ComputeAcceptance.C:23
fDDecay
Int_t fDDecay
Definition: ComputeAcceptance.C:26
LoadFONLL13TeV_promptD0
TH1D * LoadFONLL13TeV_promptD0()
Definition: ComputeAcceptance.C:445
LoadFONLL13TeV_promptDstar
TH1D * LoadFONLL13TeV_promptDstar()
Definition: ComputeAcceptance.C:487
fDecayTableFileName
TString fDecayTableFileName
Definition: ComputeAcceptance.C:32
fYMaxFidAccCut
Double_t fYMaxFidAccCut
Definition: ComputeAcceptance.C:30
kPythia7TeV
Definition: ComputeAcceptance.C:23
kFlat
Definition: ComputeAcceptance.C:23
kD0Kpi
Definition: ComputeAcceptance.C:21
kFONLL13TeVprompt
Definition: ComputeAcceptance.C:23
kDplusKpipi
Definition: ComputeAcceptance.C:21
kDsKKpi
Definition: ComputeAcceptance.C:21
LoadFONLL13TeV_promptDplus
TH1D * LoadFONLL13TeV_promptDplus()
Definition: ComputeAcceptance.C:466
LoadFONLL13TeV_feeddownD
TH1D * LoadFONLL13TeV_feeddownD()
Definition: ComputeAcceptance.C:508
kFONLL13TeVfeeddown
Definition: ComputeAcceptance.C:23
LoadPYTHIA13TeV_promptDstar
TH1D * LoadPYTHIA13TeV_promptDstar()
Definition: ComputeAcceptance.C:584
Bool_t
bool Bool_t
Definition: External.C:53
LoadPYTHIA13TeV_feeddownDstar
TH1D * LoadPYTHIA13TeV_feeddownDstar()
Definition: ComputeAcceptance.C:652
CountPKpi
Bool_t CountPKpi(TClonesArray *array, Int_t nentries, Int_t &nPions, Int_t &nKaons, Int_t &nProtons, Int_t &nPionsInAcc, Int_t &nKaonsInAcc, Int_t &nProtonsInAcc)
Definition: ComputeAcceptance.C:390
EFidY
EFidY
Definition: ComputeAcceptance.C:22
IsInFiducialAcceptance
Bool_t IsInFiducialAcceptance(Double_t pt, Double_t y)
Definition: ComputeAcceptance.C:325
LoadPYTHIA13TeV_promptDplus
TH1D * LoadPYTHIA13TeV_promptDplus()
Definition: ComputeAcceptance.C:567