AliPhysics  c923f52 (c923f52)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliAnalysisTaskSED0Mass.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Author: The ALICE Off-line Project. *
5  * Contributors are mentioned in the code where appropriate. *
6  * *
7  * Permission to use, copy, modify and distribute this software and its *
8  * documentation strictly for non-commercial purposes is hereby granted *
9  * without fee, provided that the above copyright notice appears in all *
10  * copies and that both the copyright notice and this permission notice *
11  * appear in the supporting documentation. The authors make no claims *
12  * about the suitability of this software for any purpose. It is *
13  * provided "as is" without express or implied warranty. *
14  **************************************************************************/
15 
16 /* $Id$ */
17 
19 //
20 // AliAnalysisTaskSE for D0 candidates invariant mass histogram
21 // and comparison with the MC truth and cut variables distributions.
22 //
23 // Authors: A.Dainese, andrea.dainese@lnl.infn.it
24 // Chiara Bianchin, chiara.bianchin@pd.infn.it (invariant mass)
25 // Carmelo Di Giglio, carmelo.digiglio@ba.infn.it (like sign)
26 // Jeremy Wilkinson, jwilkinson@physi.uni-heidelberg.de (weighted Bayesian
28 
29 #include <Riostream.h>
30 #include <TClonesArray.h>
31 #include <TCanvas.h>
32 #include <TNtuple.h>
33 #include <TTree.h>
34 #include <TList.h>
35 #include <TH1F.h>
36 #include <TH2F.h>
37 #include <TDatabasePDG.h>
38 #include <THnSparse.h>
39 #include "AliVertexingHFUtils.h"
40 #include <AliAnalysisDataSlot.h>
41 #include <AliAnalysisDataContainer.h>
42 #include "AliAnalysisManager.h"
43 #include "AliESDtrack.h"
44 #include "AliVertexerTracks.h"
45 #include "AliAODHandler.h"
46 #include "AliAODEvent.h"
47 #include "AliAODVertex.h"
48 #include "AliAODTrack.h"
49 #include "AliAODMCHeader.h"
50 #include "AliAODMCParticle.h"
52 #include "AliAODRecoCascadeHF.h"
53 #include "AliAnalysisVertexingHF.h"
54 #include "AliAnalysisTaskSE.h"
57 
58 using std::cout;
59 using std::endl;
60 
64 
65 //________________________________________________________________________
68  fOutputMass(0),
69  fOutputMassPt(0),
70  fOutputMassY(0),
71  fDistr(0),
72  fNentries(0),
73  fMCAccPrompt(0),
74  fMCAccBFeed(0),
75  fStepMCAcc(kTRUE),
76  fCuts(0),
77  fArray(0),
78  fReadMC(0),
79  fCutOnDistr(0),
80  fUsePid4Distr(0),
81  fCounter(0),
82  fNPtBins(1),
83  fLsNormalization(1.),
84  fFillOnlyD0D0bar(0),
85  fDaughterTracks(),
86  fIsSelectedCandidate(0),
87  fFillVarHists(kTRUE),
88  fSys(0),
89  fIsRejectSDDClusters(0),
90  fFillPtHist(kTRUE),
91  fFillYHist(kFALSE),
92  fFillImpParHist(kFALSE),
93  fUseSelectionBit(kTRUE),
94  fAODProtection(1),
95  fWriteVariableTree(kFALSE),
96  fVariablesTree(0),
97  fCandidateVariables(),
98  fPIDCheck(kFALSE),
99  fDrawDetSignal(kFALSE),
100  fUseQuarkTagInKine(kTRUE),
101  fhStudyImpParSingleTrackSign(0),
102  fhStudyImpParSingleTrackCand(0),
103  fhStudyImpParSingleTrackFd(0),
104  fDetSignal(0)
105 {
107  for(Int_t ih=0; ih<5; ih++) fHistMassPtImpParTC[ih]=0x0;
108 
109 }
110 
111 //________________________________________________________________________
113  AliAnalysisTaskSE(name),
114  fOutputMass(0),
115  fOutputMassPt(0),
116  fOutputMassY(0),
117  fDistr(0),
118  fNentries(0),
119  fMCAccPrompt(0),
120  fMCAccBFeed(0),
121  fStepMCAcc(kTRUE),
122  fCuts(0),
123  fArray(0),
124  fReadMC(0),
125  fCutOnDistr(0),
126  fUsePid4Distr(0),
127  fCounter(0),
128  fNPtBins(1),
129  fLsNormalization(1.),
130  fFillOnlyD0D0bar(0),
131  fDaughterTracks(),
132  fIsSelectedCandidate(0),
133  fFillVarHists(kTRUE),
134  fSys(0),
135  fIsRejectSDDClusters(0),
136  fFillPtHist(kTRUE),
137  fFillYHist(kFALSE),
138  fFillImpParHist(kFALSE),
139  fUseSelectionBit(kTRUE),
140  fAODProtection(1),
141  fWriteVariableTree(kFALSE),
142  fVariablesTree(0),
143  fCandidateVariables(),
144  fPIDCheck(kFALSE),
145  fDrawDetSignal(kFALSE),
146  fUseQuarkTagInKine(kTRUE),
147  fhStudyImpParSingleTrackSign(0),
148  fhStudyImpParSingleTrackCand(0),
149  fhStudyImpParSingleTrackFd(0),
150  fDetSignal(0)
151 {
153 
154  fNPtBins=cuts->GetNPtBins();
155 
156  fCuts=cuts;
157  for(Int_t ih=0; ih<5; ih++) fHistMassPtImpParTC[ih]=0x0;
158 
159  // Output slot #1 writes into a TList container (mass with cuts)
160  DefineOutput(1,TList::Class()); //My private output
161  // Output slot #2 writes into a TList container (distributions)
162  DefineOutput(2,TList::Class()); //My private output
163  // Output slot #3 writes into a TH1F container (number of events)
164  DefineOutput(3,TH1F::Class()); //My private output
165  // Output slot #4 writes into a TList container (cuts)
166  DefineOutput(4,AliRDHFCutsD0toKpi::Class()); //My private output
167  // Output slot #5 writes Normalization Counter
168  DefineOutput(5,AliNormalizationCounter::Class());
169  // Output slot #6 stores the mass vs pt and impact parameter distributions
170  DefineOutput(6,TList::Class()); //My private output
171  // Output slot #7 keeps a tree of the candidate variables after track selection
172  DefineOutput(7,TTree::Class()); //My private outpu
173  // Output slot #8 writes into a TList container (Detector signals)
174  DefineOutput(8, TList::Class()); //My private output
175  // Output slot #9 stores the mass vs rapidity (y) distributions
176  DefineOutput(9, TList::Class()); //My private output
177 }
178 
179 //________________________________________________________________________
181 {
182  if (fOutputMass) {
183  delete fOutputMass;
184  fOutputMass = 0;
185  }
186  if (fOutputMassPt) {
187  delete fOutputMassPt;
188  fOutputMassPt = 0;
189  }
190  if (fOutputMassY) {
191  delete fOutputMassY;
192  fOutputMassY = 0;
193  }
194  if (fDistr) {
195  delete fDistr;
196  fDistr = 0;
197  }
198  if (fCuts) {
199  delete fCuts;
200  fCuts = 0;
201  }
202  for(Int_t i=0; i<5; i++){
203  if(fHistMassPtImpParTC[i]) delete fHistMassPtImpParTC[i];
204  }
205  if (fNentries){
206  delete fNentries;
207  fNentries = 0;
208  }
209  if(fCounter){
210  delete fCounter;
211  fCounter=0;
212  }
213  if(fVariablesTree){
214  delete fVariablesTree;
215  fVariablesTree = 0;
216  }
217  if (fDetSignal) {
218  delete fDetSignal;
219  fDetSignal = 0;
220  }
221  delete fMCAccPrompt;
222  delete fMCAccBFeed;
226 
227 }
228 
229 //________________________________________________________________________
231 {
233 
234  if(fDebug > 1) printf("AnalysisTaskSED0Mass::Init() \n");
235 
236 
238  const char* nameoutput=GetOutputSlot(4)->GetContainer()->GetName();
239  copyfCuts->SetName(nameoutput);
240  // Post the data
241  PostData(4,copyfCuts);
242 
243 
244  return;
245 }
246 
247 //________________________________________________________________________
249 {
250 
252  //
253  if(fDebug > 1) printf("AnalysisTaskSED0Mass::UserCreateOutputObjects() \n");
254 
255  // Several histograms are more conveniently managed in a TList
256  fOutputMass = new TList();
257  fOutputMass->SetOwner();
258  fOutputMass->SetName("listMass");
259 
260  fOutputMassPt = new TList();
261  fOutputMassPt->SetOwner();
262  fOutputMassPt->SetName("listMassPt");
263 
264  fOutputMassY = new TList();
265  fOutputMassY->SetOwner();
266  fOutputMassY->SetName("listMassY");
267 
268  fDistr = new TList();
269  fDistr->SetOwner();
270  fDistr->SetName("distributionslist");
271 
272  fDetSignal = new TList();
273  fDetSignal->SetOwner();
274  fDetSignal->SetName("detectorsignals");
275 
276  TString nameMass=" ",nameSgn27=" ",nameSgn=" ", nameBkg=" ", nameRfl=" ",nameMassNocutsS =" ",nameMassNocutsB =" ", namedistr=" ";
277  TString nameMassPt="", nameSgnPt="", nameBkgPt="", nameRflPt="";
278  TString nameMassY="", nameSgnY="", nameBkgY="", nameRflY="";
279  Int_t nbins2dPt=60; Float_t binInPt=0., binFinPt=30.;
280  Int_t nbins2dY=60; Float_t binInY=-1.5, binFinY=1.5;
281 
282  for(Int_t i=0;i<fCuts->GetNPtBins();i++){
283 
284  nameMass="histMass_";
285  nameMass+=i;
286  nameSgn27="histSgn27_";
287  nameSgn27+=i;
288  nameSgn="histSgn_";
289  nameSgn+=i;
290  nameBkg="histBkg_";
291  nameBkg+=i;
292  nameRfl="histRfl_";
293  nameRfl+=i;
294  nameMassNocutsS="hMassS_";
295  nameMassNocutsS+=i;
296  nameMassNocutsB="hMassB_";
297  nameMassNocutsB+=i;
298 
299  //histograms of cut variable distributions
300  if(fFillVarHists){
301  if(fReadMC){
302 
303  namedistr="hNclsD0vsptS_";
304  namedistr+=i;
305  TH2F *hNclsD0vsptS = new TH2F(namedistr.Data(),"N cls distrubution [S];p_{T} [GeV/c];N cls",200,0.,20.,100,0.,200.);
306  namedistr="hNclsD0barvsptS_";
307  namedistr+=i;
308  TH2F *hNclsD0barvsptS = new TH2F(namedistr.Data(),"N cls distrubution [S];p_{T} [GeV/c];N cls",200,0.,20.,100,0.,200.);
309 
310  namedistr="hNITSpointsD0vsptS_";
311  namedistr+=i;
312  TH2F *hNITSpointsD0vsptS = new TH2F(namedistr.Data(),"N ITS points distrubution [S];p_{T} [GeV/c];N points",200,0.,20.,7,0.,7.);
313 
314  namedistr="hNSPDpointsD0S_";
315  namedistr+=i;
316  TH1I *hNSPDpointsD0S = new TH1I(namedistr.Data(),"N SPD points distrubution [S]; ;N tracks",4,0,4);
317  hNSPDpointsD0S->GetXaxis()->SetBinLabel(1, "no SPD");
318  hNSPDpointsD0S->GetXaxis()->SetBinLabel(2, "kOnlyFirst");
319  hNSPDpointsD0S->GetXaxis()->SetBinLabel(3, "kOnlySecond");
320  hNSPDpointsD0S->GetXaxis()->SetBinLabel(4, "kBoth");
321 
322  namedistr="hptD0S_";
323  namedistr+=i;
324  TH1F *hptD0S = new TH1F(namedistr.Data(), "p_{T} distribution [S];p_{T} [GeV/c]",200,0.,20.);
325  namedistr="hptD0barS_";
326  namedistr+=i;
327  TH1F *hptD0barS = new TH1F(namedistr.Data(), "p_{T} distribution [S];p_{T} [GeV/c]",200,0.,20.);
328 
329  namedistr="hphiD0S_";
330  namedistr+=i;
331  TH1F *hphiD0S = new TH1F(namedistr.Data(), "#phi distribution [S];#phi [rad]",100,0.,2*TMath::Pi());
332  namedistr="hphiD0barS_";
333  namedistr+=i;
334  TH1F *hphiD0barS = new TH1F(namedistr.Data(), "#phi distribution [S];#phi [rad]",100,0.,2*TMath::Pi());
335 
336 
337  namedistr="hetaphiD0candidateS_";
338  namedistr+=i;
339  TH2F *hetaphiD0candidateS = new TH2F(namedistr.Data(), "D^{0} candidates #eta #phi distribution [S];#eta;#phi [rad]",100, -1.5, 1.5, 100, 0.,2*TMath::Pi());
340  namedistr="hetaphiD0barcandidateS_";
341  namedistr+=i;
342  TH2F *hetaphiD0barcandidateS = new TH2F(namedistr.Data(), "anti-D^{0} candidates #eta #phi distribution [S];#eta;#phi [rad]",100, -1.5, 1.5, 100, 0.,2*TMath::Pi());
343 
344  namedistr="hetaphiD0candidatesignalregionS_";
345  namedistr+=i;
346  TH2F *hetaphiD0candidatesignalregionS = new TH2F(namedistr.Data(), "D^{0} candidates #eta #phi distribution [S] [mass cut];#eta;#phi [rad]",100, -1.5, 1.5, 100, 0.,2*TMath::Pi());
347  namedistr="hetaphiD0barcandidatesignalregionS_";
348  namedistr+=i;
349  TH2F *hetaphiD0barcandidatesignalregionS = new TH2F(namedistr.Data(), "anti-D^{0} candidates #eta #phi distribution [S] [mass cut];#eta;#phi [rad]",100, -1.5, 1.5, 100, 0.,2*TMath::Pi());
350 
351  // dca
352  namedistr="hdcaS_";
353  namedistr+=i;
354  TH1F *hdcaS = new TH1F(namedistr.Data(), "DCA distribution;dca [cm]",200,0.,0.1);
355  // impact parameter
356  namedistr="hd0piS_";
357  namedistr+=i;
358  TH1F *hd0piS = new TH1F(namedistr.Data(), "Impact parameter distribution (pions);d0(#pi) [cm]",200,-0.1,0.1);
359 
360  namedistr="hd0KS_";
361  namedistr+=i;
362  TH1F *hd0KS = new TH1F(namedistr.Data(), "Impact parameter distribution (kaons);d0(K) [cm]",200,-0.1,0.1);
363  namedistr="hd0d0S_";
364  namedistr+=i;
365  TH1F *hd0d0S = new TH1F(namedistr.Data(), "d_{0}#timesd_{0} distribution;d_{0}#timesd_{0} [cm^{2}]",200,-0.001,0.001);
366 
367  //decay lenght
368  namedistr="hdeclS_";
369  namedistr+=i;
370  TH1F *hdeclengthS = new TH1F(namedistr.Data(), "Decay Length^{2} distribution;Decay Length^{2} [cm]",200,0,0.015);
371 
372  namedistr="hnormdeclS_";
373  namedistr+=i;
374  TH1F *hnormdeclengthS = new TH1F(namedistr.Data(), "Normalized Decay Length^{2} distribution;(Decay Length/Err)^{2} ",200,0,12.);
375 
376  namedistr="hdeclxyS_";
377  namedistr+=i;
378  TH1F* hdeclxyS=new TH1F(namedistr.Data(),"Decay Length XY distribution;Decay Length XY [cm]",200,0,0.15);
379  namedistr="hnormdeclxyS_";
380  namedistr+=i;
381  TH1F* hnormdeclxyS=new TH1F(namedistr.Data(),"Normalized decay Length XY distribution;Decay Length XY/Err",200,0,6.);
382 
383  namedistr="hdeclxyd0d0S_";
384  namedistr+=i;
385  TH2F* hdeclxyd0d0S=new TH2F(namedistr.Data(),"Correlation decay Length XY - d_{0}#times d_{0};Decay Length XY [cm];d_{0}#times d_{0}[cm^{2}]",200,0,0.15,200,-0.001,0.001);
386 
387  namedistr="hnormdeclxyd0d0S_";
388  namedistr+=i;
389  TH2F* hnormdeclxyd0d0S=new TH2F(namedistr.Data(),"Correlation normalized decay Length XY - d_{0}#times d_{0};Decay Length XY/Err;d_{0}#times d_{0}[cm^{2}]",200,0,6,200,-0.001,0.001);
390 
391  // costhetapoint
392  namedistr="hcosthetapointS_";
393  namedistr+=i;
394  TH1F *hcosthetapointS = new TH1F(namedistr.Data(), "cos#theta_{Point} distribution;cos#theta_{Point}",200,0,1.);
395 
396  namedistr="hcosthetapointxyS_";
397  namedistr+=i;
398  TH1F *hcosthetapointxyS = new TH1F(namedistr.Data(), "cos#theta_{Point} XYdistribution;cos#theta_{Point}",300,0.,1.);
399 
400  TH1F* tmpMS = new TH1F(nameMassNocutsS.Data(),"D^{0} invariant mass; M [GeV]; Entries",600,1.6248,2.2248); //range (MD0-300MeV, mD0 + 300MeV)
401  tmpMS->Sumw2();
402 
403  fDistr->Add(hNclsD0vsptS);
404  fDistr->Add(hNclsD0barvsptS);
405  fDistr->Add(hNITSpointsD0vsptS);
406  fDistr->Add(hNSPDpointsD0S);
407  fDistr->Add(hptD0S);
408  fDistr->Add(hphiD0S);
409  fDistr->Add(hptD0barS);
410  fDistr->Add(hphiD0barS);
411  fDistr->Add(hetaphiD0candidateS);
412  fDistr->Add(hetaphiD0candidatesignalregionS);
413  fDistr->Add(hetaphiD0barcandidateS);
414  fDistr->Add(hetaphiD0barcandidatesignalregionS);
415 
416  fDistr->Add(hdcaS);
417 
418  fDistr->Add(hd0piS);
419  fDistr->Add(hd0KS);
420 
421  fDistr->Add(hd0d0S);
422 
423  fDistr->Add(hcosthetapointS);
424 
425  fDistr->Add(hcosthetapointxyS);
426 
427  fDistr->Add(hdeclengthS);
428 
429  fDistr->Add(hnormdeclengthS);
430 
431  fDistr->Add(hdeclxyS);
432 
433  fDistr->Add(hnormdeclxyS);
434 
435  fDistr->Add(hdeclxyd0d0S);
436  fDistr->Add(hnormdeclxyd0d0S);
437 
438  fDistr->Add(tmpMS);
439  }
440 
441 
442  //Ncls, phi, pt distrubutions
443 
444  namedistr="hNclsD0vsptB_";
445  namedistr+=i;
446  TH2F *hNclsD0vsptB = new TH2F(namedistr.Data(),"N cls distrubution [B];p_{T} [GeV/c];N cls",200,0.,20.,100,0.,200.);
447  namedistr="hNclsD0barvsptB_";
448  namedistr+=i;
449  TH2F *hNclsD0barvsptB = new TH2F(namedistr.Data(),"N cls distrubution [B];p_{T} [GeV/c];N cls",200,0.,20.,100,0.,200.);
450 
451  namedistr="hNITSpointsD0vsptB_";
452  namedistr+=i;
453  TH2F *hNITSpointsD0vsptB = new TH2F(namedistr.Data(),"N ITS points distrubution [B];p_{T} [GeV/c];N points",200,0.,20.,7,0.,7.);
454 
455  namedistr="hNSPDpointsD0B_";
456  namedistr+=i;
457  TH1I *hNSPDpointsD0B = new TH1I(namedistr.Data(),"N SPD points distrubution [B]; ;N tracks",4,0,4);
458  hNSPDpointsD0B->GetXaxis()->SetBinLabel(1, "no SPD");
459  hNSPDpointsD0B->GetXaxis()->SetBinLabel(2, "kOnlyFirst");
460  hNSPDpointsD0B->GetXaxis()->SetBinLabel(3, "kOnlySecond");
461  hNSPDpointsD0B->GetXaxis()->SetBinLabel(4, "kBoth");
462 
463  namedistr="hptD0B_";
464  namedistr+=i;
465  TH1F *hptD0B = new TH1F(namedistr.Data(), "p_{T} distribution [B];p_{T} [GeV/c]",200,0.,20.);
466  namedistr="hptD0barB_";
467  namedistr+=i;
468  TH1F *hptD0barB = new TH1F(namedistr.Data(), "p_{T} distribution [B];p_{T} [GeV/c]",200,0.,20.);
469 
470  namedistr="hphiD0B_";
471  namedistr+=i;
472  TH1F *hphiD0B = new TH1F(namedistr.Data(), "#phi distribution [B];#phi [rad]",100,0.,2*TMath::Pi());
473  namedistr="hphiD0barB_";
474  namedistr+=i;
475  TH1F *hphiD0barB = new TH1F(namedistr.Data(), "#phi distribution [B];#phi [rad]",100,0.,2*TMath::Pi());
476 
477  namedistr="hetaphiD0candidateB_";
478  namedistr+=i;
479  TH2F *hetaphiD0candidateB = new TH2F(namedistr.Data(), "D^{0} candidates #eta #phi distribution [B];#eta;#phi [rad]",100, -1.5, 1.5, 100, 0.,2*TMath::Pi());
480  namedistr="hetaphiD0barcandidateB_";
481  namedistr+=i;
482  TH2F *hetaphiD0barcandidateB = new TH2F(namedistr.Data(), "anti-D^{0} candidates #eta #phi distribution [B];#eta;#phi [rad]",100, -1.5, 1.5, 100, 0.,2*TMath::Pi());
483 
484  namedistr="hetaphiD0candidatesignalregionB_";
485  namedistr+=i;
486  TH2F *hetaphiD0candidatesignalregionB = new TH2F(namedistr.Data(), "D^{0} candidates #eta #phi distribution [B] [mass cut];#eta;#phi [rad]",100, -1.5, 1.5, 100, 0.,2*TMath::Pi());
487  namedistr="hetaphiD0barcandidatesignalregionB_";
488  namedistr+=i;
489  TH2F *hetaphiD0barcandidatesignalregionB = new TH2F(namedistr.Data(), "anti-D^{0} candidates #eta #phi distribution [B] [mass cut];#eta;#phi [rad]",100, -1.5, 1.5, 100, 0.,2*TMath::Pi());
490 
491  // dca
492  namedistr="hdcaB_";
493  namedistr+=i;
494  TH1F *hdcaB = new TH1F(namedistr.Data(), "DCA distribution;dca [cm]",200,0.,0.1);
495 
496  // impact parameter
497  namedistr="hd0B_";
498  namedistr+=i;
499  TH1F *hd0B = new TH1F(namedistr.Data(), "Impact parameter distribution (both);d0 [cm]",200,-0.1,0.1);
500 
501  namedistr="hd0d0B_";
502  namedistr+=i;
503  TH1F *hd0d0B = new TH1F(namedistr.Data(), "d_{0}#timesd_{0} distribution;d_{0}#timesd_{0} [cm^{2}]",200,-0.001,0.001);
504 
505  //decay lenght
506  namedistr="hdeclB_";
507  namedistr+=i;
508  TH1F *hdeclengthB = new TH1F(namedistr.Data(), "Decay Length^{2} distribution;Decay Length^{2} [cm^{2}]",200,0,0.015);
509 
510  namedistr="hnormdeclB_";
511  namedistr+=i;
512  TH1F *hnormdeclengthB = new TH1F(namedistr.Data(), "Normalized Decay Length distribution;(Decay Length/Err)^{2} ",200,0,12.);
513 
514  namedistr="hdeclxyB_";
515  namedistr+=i;
516  TH1F* hdeclxyB=new TH1F(namedistr.Data(),"Decay Length XY distribution;Decay Length XY [cm]",200,0,0.15);
517  namedistr="hnormdeclxyB_";
518  namedistr+=i;
519  TH1F* hnormdeclxyB=new TH1F(namedistr.Data(),"Normalized decay Length XY distribution;Decay Length XY/Err",200,0,6.);
520 
521  namedistr="hdeclxyd0d0B_";
522  namedistr+=i;
523  TH2F* hdeclxyd0d0B=new TH2F(namedistr.Data(),"Correlation decay Length XY - d_{0}#times d_{0};Decay Length XY [cm];d_{0}#times d_{0}[cm^{2}]",200,0,0.15,200,-0.001,0.001);
524 
525  namedistr="hnormdeclxyd0d0B_";
526  namedistr+=i;
527  TH2F* hnormdeclxyd0d0B=new TH2F(namedistr.Data(),"Correlation normalized decay Length XY - d_{0}#times d_{0};Decay Length XY/Err;d_{0}#times d_{0}[cm^{2}]",200,0,6,200,-0.001,0.001);
528 
529  // costhetapoint
530  namedistr="hcosthetapointB_";
531  namedistr+=i;
532  TH1F *hcosthetapointB = new TH1F(namedistr.Data(), "cos#theta_{Point} distribution;cos#theta_{Point}",200,0,1.);
533 
534  namedistr="hcosthetapointxyB_";
535  namedistr+=i;
536  TH1F *hcosthetapointxyB = new TH1F(namedistr.Data(), "cos#theta_{Point} XY distribution;cos#theta_{Point} XY",300,0.,1.);
537 
538  TH1F* tmpMB = new TH1F(nameMassNocutsB.Data(),"D^{0} invariant mass; M [GeV]; Entries",600,1.6248,2.2248); //range (MD0-300MeV, mD0 + 300MeV)
539  tmpMB->Sumw2();
540 
541 
542  fDistr->Add(hNclsD0vsptB);
543  fDistr->Add(hNclsD0barvsptB);
544  fDistr->Add(hNITSpointsD0vsptB);
545  fDistr->Add(hNSPDpointsD0B);
546  fDistr->Add(hptD0B);
547  fDistr->Add(hphiD0B);
548  fDistr->Add(hptD0barB);
549  fDistr->Add(hphiD0barB);
550  fDistr->Add(hetaphiD0candidateB);
551  fDistr->Add(hetaphiD0candidatesignalregionB);
552  fDistr->Add(hetaphiD0barcandidateB);
553  fDistr->Add(hetaphiD0barcandidatesignalregionB);
554 
555  fDistr->Add(hdcaB);
556 
557  fDistr->Add(hd0B);
558 
559  fDistr->Add(hd0d0B);
560 
561  fDistr->Add(hcosthetapointB);
562 
563  fDistr->Add(hcosthetapointxyB);
564 
565  fDistr->Add(hdeclengthB);
566 
567  fDistr->Add(hnormdeclengthB);
568 
569  fDistr->Add(hdeclxyB);
570 
571  fDistr->Add(hnormdeclxyB);
572 
573  fDistr->Add(hdeclxyd0d0B);
574  fDistr->Add(hnormdeclxyd0d0B);
575 
576  fDistr->Add(tmpMB);
577 
578  //histograms filled only when the secondary vertex is recalculated w/o the daughter tracks (as requested in the cut object)
579 
581  if(fReadMC){
582  namedistr="hd0vpiS_";
583  namedistr+=i;
584  TH1F *hd0vpiS = new TH1F(namedistr.Data(), "Impact parameter distribution (pions)(vtx w/o these tracks);d0(#pi) [cm]",200,-0.1,0.1);
585 
586  namedistr="hd0vKS_";
587  namedistr+=i;
588  TH1F *hd0vKS = new TH1F(namedistr.Data(), "Impact parameter distribution (kaons) (vtx w/o these tracks);d0(K) [cm]",200,-0.1,0.1);
589 
590  namedistr="hd0d0vS_";
591  namedistr+=i;
592  TH1F *hd0d0vS = new TH1F(namedistr.Data(), "d_{0}#timesd_{0} distribution (vtx w/o these tracks);d_{0}#timesd_{0} [cm^{2}]",200,-0.001,0.001);
593 
594  namedistr="hdeclvS_";
595  namedistr+=i;
596  TH1F *hdeclengthvS = new TH1F(namedistr.Data(), "Decay Length distribution (vtx w/o tracks);Decay Length [cm]",200,0,0.6);
597 
598  namedistr="hnormdeclvS_";
599  namedistr+=i;
600  TH1F *hnormdeclengthvS = new TH1F(namedistr.Data(), "Normalized Decay Length distribution (vtx w/o tracks);Decay Length/Err ",200,0,10.);
601 
602  fDistr->Add(hd0vpiS);
603  fDistr->Add(hd0vKS);
604  fDistr->Add(hd0d0vS);
605  fDistr->Add(hdeclengthvS);
606  fDistr->Add(hnormdeclengthvS);
607 
608  }
609 
610  namedistr="hd0vmoresB_";
611  namedistr+=i;
612  TH1F *hd0vmoresB = new TH1F(namedistr.Data(), "Impact parameter distribution (both);d0 [cm]",200,-0.1,0.1);
613 
614  namedistr="hd0d0vmoresB_";
615  namedistr+=i;
616  TH1F *hd0d0vmoresB = new TH1F(namedistr.Data(), "Impact parameter distribution (prong +);d0 [cm]",200,-0.001,0.001);
617 
618 
619  namedistr="hd0vB_";
620  namedistr+=i;
621  TH1F *hd0vB = new TH1F(namedistr.Data(), "Impact parameter distribution (vtx w/o these tracks);d0 [cm]",200,-0.1,0.1);
622 
623  namedistr="hd0vp0B_";
624  namedistr+=i;
625  TH1F *hd0vp0B = new TH1F(namedistr.Data(), "Impact parameter distribution (prong + ** vtx w/o these tracks);d0 [cm]",200,-0.1,0.1);
626 
627  namedistr="hd0vp1B_";
628  namedistr+=i;
629  TH1F *hd0vp1B = new TH1F(namedistr.Data(), "Impact parameter distribution (prong - ** vtx w/o these tracks);d0 [cm]",200,-0.1,0.1);
630 
631  namedistr="hd0d0vB_";
632  namedistr+=i;
633  TH1F *hd0d0vB = new TH1F(namedistr.Data(), "d_{0}#timesd_{0} distribution (vtx w/o these tracks);d_{0}#timesd_{0} [cm^{2}]",200,-0.001,0.001);
634 
635  namedistr="hdeclvB_";
636  namedistr+=i;
637  TH1F *hdeclengthvB = new TH1F(namedistr.Data(), "Decay Length distribution (vtx w/o tracks);Decay Length [cm]",200,0,0.6);
638 
639  namedistr="hnormdeclvB_";
640  namedistr+=i;
641  TH1F *hnormdeclengthvB = new TH1F(namedistr.Data(), "Normalized Decay Length distribution (vtx w/o tracks);Decay Length/Err ",200,0,10.);
642 
643  fDistr->Add(hd0vB);
644  fDistr->Add(hd0vp0B);
645  fDistr->Add(hd0vp1B);
646  fDistr->Add(hd0vmoresB);
647 
648  fDistr->Add(hd0d0vB);
649  fDistr->Add(hd0d0vmoresB);
650 
651  fDistr->Add(hdeclengthvB);
652 
653  fDistr->Add(hnormdeclengthvB);
654  }
655 
656 
657  }
658  //histograms of invariant mass distributions
659 
660 
661  //MC signal
662  if(fReadMC){
663  TH1F* tmpSt = new TH1F(nameSgn.Data(), "D^{0} invariant mass - MC; M [GeV]; Entries",600,1.6248,2.2248);
664 
665  TH1F *tmpSl=(TH1F*)tmpSt->Clone();
666  tmpSt->Sumw2();
667  tmpSl->Sumw2();
668 
669  //Reflection: histo filled with D0Mass which pass the cut (also) as D0bar and with D0bar which pass (also) the cut as D0
670  TH1F* tmpRt = new TH1F(nameRfl.Data(), "Reflected signal invariant mass - MC; M [GeV]; Entries",600,1.6248,2.2248);
671  //TH1F *tmpRl=(TH1F*)tmpRt->Clone();
672  TH1F* tmpBt = new TH1F(nameBkg.Data(), "Background invariant mass - MC; M [GeV]; Entries",600,1.6248,2.2248);
673  //TH1F *tmpBl=(TH1F*)tmpBt->Clone();
674  tmpBt->Sumw2();
675  //tmpBl->Sumw2();
676  tmpRt->Sumw2();
677  //tmpRl->Sumw2();
678 
679  fOutputMass->Add(tmpSt);
680  fOutputMass->Add(tmpRt);
681  fOutputMass->Add(tmpBt);
682 
683  }
684 
685  //mass
686  TH1F* tmpMt = new TH1F(nameMass.Data(),"D^{0} invariant mass; M [GeV]; Entries",600,1.6248,2.2248);
687  //TH1F *tmpMl=(TH1F*)tmpMt->Clone();
688  tmpMt->Sumw2();
689  //tmpMl->Sumw2();
690  //distribution w/o cuts large range
691  // TH1F* tmpMS = new TH1F(nameMassNocutsS.Data(),"D^{0} invariant mass; M [GeV]; Entries",300,0.7,3.);
692 
693  fOutputMass->Add(tmpMt);
694 
695  if(fSys==0){ //histograms filled only in pp to save time in PbPb
696  if(fFillVarHists){
697 
698  if(fReadMC){
699  // pT
700  namedistr="hptpiS_";
701  namedistr+=i;
702  TH1F *hptpiS = new TH1F(namedistr.Data(), "P_{T} distribution (pions);p_{T} [GeV/c]",200,0.,8.);
703 
704  namedistr="hptKS_";
705  namedistr+=i;
706  TH1F *hptKS = new TH1F(namedistr.Data(), "P_{T} distribution (kaons);p_{T} [GeV/c]",200,0.,8.);
707 
708  // costhetastar
709  namedistr="hcosthetastarS_";
710  namedistr+=i;
711  TH1F *hcosthetastarS = new TH1F(namedistr.Data(), "cos#theta* distribution;cos#theta*",200,-1.,1.);
712 
713  //pT no mass cut
714 
715  namedistr="hptpiSnoMcut_";
716  namedistr+=i;
717  TH1F *hptpiSnoMcut = new TH1F(namedistr.Data(), "P_{T} distribution (pions);p_{T} [GeV/c]",200,0.,8.);
718 
719  namedistr="hptKSnoMcut_";
720  namedistr+=i;
721  TH1F *hptKSnoMcut = new TH1F(namedistr.Data(), "P_{T} distribution (kaons);p_{T} [GeV/c]",200,0.,8.);
722 
723  fDistr->Add(hptpiS);
724  fDistr->Add(hptKS);
725  fDistr->Add(hcosthetastarS);
726 
727  fDistr->Add(hptpiSnoMcut);
728  fDistr->Add(hptKSnoMcut);
729 
730  // costhetapoint vs d0 or d0d0
731  namedistr="hcosthpointd0S_";
732  namedistr+=i;
733  TH2F *hcosthpointd0S= new TH2F(namedistr.Data(),"Correlation cos#theta_{Point}-d_{0};cos#theta_{Point};d_{0} [cm^{2}]",200,0,1.,200,-0.001,0.001);
734  namedistr="hcosthpointd0d0S_";
735  namedistr+=i;
736  TH2F *hcosthpointd0d0S= new TH2F(namedistr.Data(),"Correlation cos#theta_{Point}-d_{0}#timesd_{0};cos#theta_{Point};d_{0}#timesd_{0} [cm^{2}]",200,0,1.,200,-0.001,0.001);
737 
738  fDistr->Add(hcosthpointd0S);
739  fDistr->Add(hcosthpointd0d0S);
740 
741  //to compare with AliAnalysisTaskCharmFraction
742  TH1F* tmpS27t = new TH1F(nameSgn27.Data(),"D^{0} invariant mass in M(D^{0}) +/- 27 MeV - MC; M [GeV]; Entries",600,1.6248,2.2248);
743  TH1F *tmpS27l=(TH1F*)tmpS27t->Clone();
744  tmpS27t->Sumw2();
745  tmpS27l->Sumw2();
746 
747  fOutputMass->Add(tmpS27t);
748  fOutputMass->Add(tmpS27l);
749 
750  }
751 
752  // pT
753  namedistr="hptB_";
754  namedistr+=i;
755  TH1F *hptB = new TH1F(namedistr.Data(), "P_{T} distribution;p_{T} [GeV/c]",200,0.,8.);
756 
757  // costhetastar
758  namedistr="hcosthetastarB_";
759  namedistr+=i;
760  TH1F *hcosthetastarB = new TH1F(namedistr.Data(), "cos#theta* distribution;cos#theta*",200,-1.,1.);
761 
762  //pT no mass cut
763  namedistr="hptB1prongnoMcut_";
764  namedistr+=i;
765  TH1F *hptB1pnoMcut = new TH1F(namedistr.Data(), "P_{T} distribution;p_{T} [GeV/c]",200,0.,8.);
766 
767  namedistr="hptB2prongsnoMcut_";
768  namedistr+=i;
769  TH1F *hptB2pnoMcut = new TH1F(namedistr.Data(), "P_{T} distribution;p_{T} [GeV/c]",200,0.,8.);
770 
771  fDistr->Add(hptB);
772  fDistr->Add(hcosthetastarB);
773 
774  fDistr->Add(hptB1pnoMcut);
775  fDistr->Add(hptB2pnoMcut);
776 
777  //impact parameter of negative/positive track
778  namedistr="hd0p0B_";
779  namedistr+=i;
780  TH1F *hd0p0B = new TH1F(namedistr.Data(), "Impact parameter distribution (prong +);d0 [cm]",200,-0.1,0.1);
781 
782  namedistr="hd0p1B_";
783  namedistr+=i;
784  TH1F *hd0p1B = new TH1F(namedistr.Data(), "Impact parameter distribution (prong -);d0 [cm]",200,-0.1,0.1);
785 
786  //impact parameter corrected for strangeness
787  namedistr="hd0moresB_";
788  namedistr+=i;
789  TH1F *hd0moresB = new TH1F(namedistr.Data(), "Impact parameter distribution (both);d0 [cm]",200,-0.1,0.1);
790 
791  namedistr="hd0d0moresB_";
792  namedistr+=i;
793  TH1F *hd0d0moresB = new TH1F(namedistr.Data(), "Impact parameter distribution (prong +);d0 [cm]",200,-0.001,0.001);
794 
795 
796  namedistr="hcosthetapointmoresB_";
797  namedistr+=i;
798  TH1F *hcosthetapointmoresB = new TH1F(namedistr.Data(), "cos#theta_{Point} distribution;cos#theta_{Point}",200,0,1.);
799 
800  // costhetapoint vs d0 or d0d0
801  namedistr="hcosthpointd0B_";
802  namedistr+=i;
803  TH2F *hcosthpointd0B= new TH2F(namedistr.Data(),"Correlation cos#theta_{Point}-d_{0};cos#theta_{Point};d_{0} [cm^{2}]",200,0,1.,200,-0.001,0.001);
804 
805  namedistr="hcosthpointd0d0B_";
806  namedistr+=i;
807  TH2F *hcosthpointd0d0B= new TH2F(namedistr.Data(),"Correlation cos#theta_{Point}-d_{0}#timesd_{0};cos#theta_{Point};d_{0}#timesd_{0} [cm^{2}]",200,0,1.,200,-0.001,0.001);
808 
809  fDistr->Add(hd0p0B);
810  fDistr->Add(hd0p1B);
811 
812  fDistr->Add(hd0moresB);
813  fDistr->Add(hd0d0moresB);
814  fDistr->Add(hcosthetapointmoresB);
815 
816 
817  fDistr->Add(hcosthpointd0B);
818 
819 
820  fDistr->Add(hcosthpointd0d0B);
821  }
822 
823  } //end pp histo
824  }
825  //create THnSparse for impact parameter analysis in DATA sample.
826  //pt, normImpParTrk1, normImpParTrk2, decLXY, normDecLXY, massd0, cut, pid, D0D0bar
827  Int_t nbinsImpParStudy[9]= {50, 40, 40, 20, 15, 600, 3,4,2};
828  Double_t limitLowImpParStudy[9]={0, -5, -5, 0, 0,1.6248,1,0,0};
829  Double_t limitUpImpParStudy[9]= {50., 5, 5, 0.2,15,2.2248,4,4,2};
830  TString axTit[9]={"#it{p}_{T} (GeV/c)","normalized imp par residual, trk1","normalized imp par residual, trk2","#it{L}_{xy} (cm)","norm #it{L}_{xy}","MassD0_{k#pi} (GeV/#it{c}^{2})", "cutSel","PIDinfo","D0D0bar"};
831  fhStudyImpParSingleTrackCand=new THnSparseF("fhStudyImpParSingleTrackCand","fhStudyImpParSingleTrackCand",9,nbinsImpParStudy,limitLowImpParStudy,limitUpImpParStudy);
832  for(Int_t iax=0; iax<9; iax++) fhStudyImpParSingleTrackCand->GetAxis(iax)->SetTitle(axTit[iax].Data());
834  if(fReadMC){
835  //pt, ptB, normImpParTrk1, normImpParTrk2, decLXY, normDecLXY, iscut, ispid
836  Int_t nbinsImpParStudy[8]= {50,50,40, 40, 20, 15, 3, 4};
837  Double_t limitLowImpParStudy[8]={0, 0, -5,-5., 0., 0., 1.,0.};
838  Double_t limitUpImpParStudy[8]= {50.,50., 5, 5, 0.2, 15, 3.,4.};
839 
840  fhStudyImpParSingleTrackSign=new THnSparseF("fhStudyImpParSingleTrackSign","fhStudyImpParSingleTrackSign",8,nbinsImpParStudy,limitLowImpParStudy,limitUpImpParStudy);
841  TString axTitMC[8]={"#it{p}_{T} (GeV/c)","#it{p}_{T} (GeV/c)","normalized imp par residual, trk1","normalized imp par residual, trk2","#it{L}_{xy} (cm)","norm #it{L}_{xy}","cutSel","PIDinfo"};
842  for(Int_t iax=0; iax<8; iax++) fhStudyImpParSingleTrackSign->GetAxis(iax)->SetTitle(axTitMC[iax].Data());
844 
845 
846  fhStudyImpParSingleTrackFd=new THnSparseF("fhStudyImpParSingleTrackFd","fhStudyImpParSingleTrackFd",8,nbinsImpParStudy,limitLowImpParStudy,limitUpImpParStudy);
847  for(Int_t iax=0; iax<8; iax++) fhStudyImpParSingleTrackFd->GetAxis(iax)->SetTitle(axTitMC[iax].Data());
849 
851  }
852 
853 
854  //for Like sign analysis
855 
856  if(fArray==1){
857  namedistr="hpospair";
858  TH1F* hpospair=new TH1F(namedistr.Data(),"Number of positive pairs",fCuts->GetNPtBins(),-0.5,fCuts->GetNPtBins()-0.5);
859  namedistr="hnegpair";
860  TH1F* hnegpair=new TH1F(namedistr.Data(),"Number of negative pairs",fCuts->GetNPtBins(),-0.5,fCuts->GetNPtBins()-0.5);
861  fOutputMass->Add(hpospair);
862  fOutputMass->Add(hnegpair);
863  }
864 
865 
866  // 2D Pt distributions and impact parameter histograms
867  if(fFillPtHist) {
868 
869  nameMassPt="histMassPt";
870  nameSgnPt="histSgnPt";
871  nameBkgPt="histBkgPt";
872  nameRflPt="histRflPt";
873 
874  //MC signal
875  if(fReadMC){
876  TH2F* tmpStPt = new TH2F(nameSgnPt.Data(), "D^{0} invariant mass - MC; M [GeV]; Entries; Pt[GeV/c]",600,1.6248,2.2248,nbins2dPt,binInPt,binFinPt);
877  TH2F *tmpSlPt=(TH2F*)tmpStPt->Clone();
878  tmpStPt->Sumw2();
879  tmpSlPt->Sumw2();
880 
881  //Reflection: histo filled with D0MassV1 which pass the cut (also) as D0bar and with D0bar which pass (also) the cut as D0
882  TH2F* tmpRtPt = new TH2F(nameRflPt.Data(), "Reflected signal invariant mass - MC; M [GeV]; Entries; Pt[GeV/c]",600,1.6248,2.2248,nbins2dPt,binInPt,binFinPt);
883  TH2F* tmpBtPt = new TH2F(nameBkgPt.Data(), "Background invariant mass - MC; M [GeV]; Entries; Pt[GeV/c]",600,1.6248,2.2248,nbins2dPt,binInPt,binFinPt);
884  tmpBtPt->Sumw2();
885  tmpRtPt->Sumw2();
886 
887  fOutputMassPt->Add(tmpStPt);
888  fOutputMassPt->Add(tmpRtPt);
889  fOutputMassPt->Add(tmpBtPt);
890 
891  // cout<<endl<<endl<<"***************************************"<<endl;
892  // cout << " created and added to the list "<< nameSgnPt.Data() <<" "<< tmpStPt <<
893  // ", "<<nameRflPt.Data() <<" " << tmpRtPt<<", "<<nameBkgPt.Data()<< " " << tmpBtPt <<endl;
894  // cout<<"***************************************"<<endl<<endl;
895  }
896 
897  TH2F* tmpMtPt = new TH2F(nameMassPt.Data(),"D^{0} invariant mass; M [GeV]; Entries; Pt[GeV/c]",600,1.6248,2.2248,nbins2dPt,binInPt,binFinPt);
898  tmpMtPt->Sumw2();
899 
900  fOutputMassPt->Add(tmpMtPt);
901  }
902 
904 
905  // 2D Y distributions
906 
907  if(fFillYHist) {
908  for(Int_t i=0;i<fCuts->GetNPtBins();i++){
909  nameMassY="histMassY_";
910  nameMassY+=i;
911  nameSgnY="histSgnY_";
912  nameSgnY+=i;
913  nameBkgY="histBkgY_";
914  nameBkgY+=i;
915  nameRflY="histRflY_";
916  nameRflY+=i;
917  //MC signal
918  if(fReadMC){
919  TH2F* tmpStY = new TH2F(nameSgnY.Data(), "D^{0} invariant mass - MC; M [GeV]; Entries; y",600,1.6248,2.2248,nbins2dY,binInY,binFinY);
920  tmpStY->Sumw2();
921  //Reflection: histo filled with D0MassV1 which pass the cut (also) as D0bar and with D0bar which pass (also) the cut as D0
922  TH2F* tmpRtY = new TH2F(nameRflY.Data(), "Reflected signal invariant mass - MC; M [GeV]; Entries; y",600,1.6248,2.2248,nbins2dY,binInY,binFinY);
923  TH2F* tmpBtY = new TH2F(nameBkgY.Data(), "Background invariant mass - MC; M [GeV]; Entries; y",600,1.6248,2.2248,nbins2dY,binInY,binFinY);
924  tmpBtY->Sumw2();
925  tmpRtY->Sumw2();
926 
927  fOutputMassY->Add(tmpStY);
928  fOutputMassY->Add(tmpRtY);
929  fOutputMassY->Add(tmpBtY);
930  }
931  TH2F* tmpMtY = new TH2F(nameMassY.Data(),"D^{0} invariant mass; M [GeV]; Entries; y",600,1.6248,2.2248,nbins2dY,binInY,binFinY);
932  tmpMtY->Sumw2();
933  fOutputMassY->Add(tmpMtY);
934  }
935  }
936 
937 
938  const char* nameoutput=GetOutputSlot(3)->GetContainer()->GetName();
939 
940  fNentries=new TH1F(nameoutput, "Integral(1,2) = number of AODs *** Integral(2,3) = number of candidates selected with cuts *** Integral(3,4) = number of D0 selected with cuts *** Integral(4,5) = events with good vertex *** Integral(5,6) = pt out of bounds", 23,-0.5,22.5);
941 
942  fNentries->GetXaxis()->SetBinLabel(1,"nEventsAnal");
943  fNentries->GetXaxis()->SetBinLabel(2,"nCandSel(Cuts)");
944  if(fReadMC) fNentries->GetXaxis()->SetBinLabel(3,"nD0Selected");
945  else fNentries->GetXaxis()->SetBinLabel(3,"Dstar<-D0");
946  fNentries->GetXaxis()->SetBinLabel(4,"nEventsGoodVtxS");
947  fNentries->GetXaxis()->SetBinLabel(5,"ptbin = -1");
948  fNentries->GetXaxis()->SetBinLabel(6,"no daughter");
949  if(fSys==0) fNentries->GetXaxis()->SetBinLabel(7,"nCandSel(Tr)");
950  if(fFillVarHists || fPIDCheck){
951  fNentries->GetXaxis()->SetBinLabel(8,"PID=0");
952  fNentries->GetXaxis()->SetBinLabel(9,"PID=1");
953  fNentries->GetXaxis()->SetBinLabel(10,"PID=2");
954  fNentries->GetXaxis()->SetBinLabel(11,"PID=3");
955  }
956  if(fReadMC && fSys==0){
957  fNentries->GetXaxis()->SetBinLabel(12,"K");
958  fNentries->GetXaxis()->SetBinLabel(13,"Lambda");
959  }
960  fNentries->GetXaxis()->SetBinLabel(14,"Pile-up Rej");
961  fNentries->GetXaxis()->SetBinLabel(15,"N. of 0SMH");
962  if(fSys==1) fNentries->GetXaxis()->SetBinLabel(16,"Nev in centr");
963  if(fIsRejectSDDClusters) fNentries->GetXaxis()->SetBinLabel(17,"SDD-Cls Rej");
964  fNentries->GetXaxis()->SetBinLabel(18,"Phys.Sel.Rej");
965  fNentries->GetXaxis()->SetBinLabel(19,"D0 failed to be filled");
966  fNentries->GetXaxis()->SetBinLabel(20,"fisFilled is 0");
967  fNentries->GetXaxis()->SetBinLabel(21,"fisFilled is 1");
968  fNentries->GetXaxis()->SetBinLabel(22,"AOD/dAOD mismatch");
969  fNentries->GetXaxis()->SetBinLabel(23,"AOD/dAOD #events ok");
970  fNentries->GetXaxis()->SetNdivisions(1,kFALSE);
971 
972  fCounter = new AliNormalizationCounter(Form("%s",GetOutputSlot(5)->GetContainer()->GetName()));
973  fCounter->Init();
974 
975  //
976  // Output slot 7 : tree of the candidate variables after track selection
977  //
978  nameoutput = GetOutputSlot(7)->GetContainer()->GetName();
979  fVariablesTree = new TTree(nameoutput,"Candidates variables tree");
980  Int_t nVar = 15;
982  TString * fCandidateVariableNames = new TString[nVar];
983  fCandidateVariableNames[0] = "massD0";
984  fCandidateVariableNames[1] = "massD0bar";
985  fCandidateVariableNames[2] = "pt";
986  fCandidateVariableNames[3] = "dca";
987  fCandidateVariableNames[4] = "costhsD0";
988  fCandidateVariableNames[5] = "costhsD0bar";
989  fCandidateVariableNames[6] = "ptk";
990  fCandidateVariableNames[7] = "ptpi";
991  fCandidateVariableNames[8] = "d0k";
992  fCandidateVariableNames[9] = "d0pi";
993  fCandidateVariableNames[10] = "d0xd0";
994  fCandidateVariableNames[11] = "costhp";
995  fCandidateVariableNames[12] = "costhpxy";
996  fCandidateVariableNames[13] = "lxy";
997  fCandidateVariableNames[14] = "specialcuts";
998  for(Int_t ivar=0; ivar<nVar; ivar++){
999  fVariablesTree->Branch(fCandidateVariableNames[ivar].Data(),&fCandidateVariables[ivar],Form("%s/d",fCandidateVariableNames[ivar].Data()));
1000  }
1001 
1002 
1003  //
1004  // Output slot 8 : List for detector response histograms
1005  //
1006  if (fDrawDetSignal) {
1007  TH2F *TOFSigBefPID = new TH2F("TOFSigBefPID", "TOF signal of daughters before PID;p(daught)(GeV/c);Signal", 500, 0, 10, 5000, 0, 50e3);
1008  TH2F *TOFSigAftPID = new TH2F("TOFSigAftPID", "TOF signal after PID;p(daught)(GeV/c);Signal", 500, 0, 10, 5000, 0, 50e3);
1009 
1010  TH2F *TPCSigBefPID = new TH2F("TPCSigBefPID", "TPC dE/dx before PID;p(daught)(GeV/c);dE/dx (arb. units)", 1000, 0, 10, 1000, 0, 500);
1011  TH2F *TPCSigAftPID = new TH2F("TPCSigAftPID", "TPC dE/dx after PID;p(daught)(GeV/c);dE/dx (arb. units)", 1000, 0, 10, 1000, 0, 500);
1012 
1013  fDetSignal->Add(TOFSigBefPID);
1014  fDetSignal->Add(TOFSigAftPID);
1015  fDetSignal->Add(TPCSigBefPID);
1016  fDetSignal->Add(TPCSigAftPID);
1017  }
1018 
1019  // Post the data
1020  PostData(1,fOutputMass);
1021  PostData(2,fDistr);
1022  PostData(3,fNentries);
1023  PostData(5,fCounter);
1024  PostData(6,fOutputMassPt);
1025  PostData(7,fVariablesTree);
1026  PostData(8, fDetSignal);
1027  PostData(9,fOutputMassY);
1028 
1029  return;
1030 }
1031 
1032 //________________________________________________________________________
1034 
1035 {
1038 
1039  //cuts order
1040  // printf(" |M-MD0| [GeV] < %f\n",fD0toKpiCuts[0]);
1041  // printf(" dca [cm] < %f\n",fD0toKpiCuts[1]);
1042  // printf(" cosThetaStar < %f\n",fD0toKpiCuts[2]);
1043  // printf(" pTK [GeV/c] > %f\n",fD0toKpiCuts[3]);
1044  // printf(" pTpi [GeV/c] > %f\n",fD0toKpiCuts[4]);
1045  // printf(" |d0K| [cm] < %f\n",fD0toKpiCuts[5]);
1046  // printf(" |d0pi| [cm] < %f\n",fD0toKpiCuts[6]);
1047  // printf(" d0d0 [cm^2] < %f\n",fD0toKpiCuts[7]);
1048  // printf(" cosThetaPoint > %f\n",fD0toKpiCuts[8]);
1049 
1050  AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
1051 
1052  TString bname;
1053  if(fArray==0){ //D0 candidates
1054  // load D0->Kpi candidates
1055  //cout<<"D0 candidates"<<endl;
1056  bname="D0toKpi";
1057  } else { //LikeSign candidates
1058  // load like sign candidates
1059  //cout<<"LS candidates"<<endl;
1060  bname="LikeSign2Prong";
1061  }
1062  TClonesArray *inputArray=0;
1063  if(!aod && AODEvent() && IsStandardAOD()) {
1064  // In case there is an AOD handler writing a standard AOD, use the AOD
1065  // event in memory rather than the input (ESD) event.
1066  aod = dynamic_cast<AliAODEvent*> (AODEvent());
1067  // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
1068  // have to taken from the AOD event hold by the AliAODExtension
1069  AliAODHandler* aodHandler = (AliAODHandler*)
1070  ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
1071 
1072  if(aodHandler->GetExtensions()) {
1073  AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
1074  AliAODEvent* aodFromExt = ext->GetAOD();
1075  inputArray=(TClonesArray*)aodFromExt->GetList()->FindObject(bname.Data());
1076  }
1077  } else if(aod) {
1078  inputArray=(TClonesArray*)aod->GetList()->FindObject(bname.Data());
1079  }
1080 
1081  if(!inputArray || !aod) {
1082  printf("AliAnalysisTaskSED0Mass::UserExec: input branch not found!\n");
1083  return;
1084  }
1085  // fix for temporary bug in ESDfilter
1086  // the AODs with null vertex pointer didn't pass the PhysSel
1087  if(!aod->GetPrimaryVertex() || TMath::Abs(aod->GetMagneticField())<0.001) return;
1088 
1089  TClonesArray *mcArray = 0;
1090  AliAODMCHeader *mcHeader = 0;
1091 
1092  if(fReadMC) {
1093  // load MC particles
1094  mcArray = (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName());
1095  if(!mcArray) {
1096  printf("AliAnalysisTaskSED0Mass::UserExec: MC particles branch not found!\n");
1097  return;
1098  }
1099 
1100  // load MC header
1101  mcHeader = (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
1102  if(!mcHeader) {
1103  printf("AliAnalysisTaskSED0Mass::UserExec: MC header branch not found!\n");
1104  return;
1105  }
1106  }
1107  //printf("VERTEX Z %f %f\n",vtx1->GetZ(),mcHeader->GetVtxZ());
1108 
1109  //histogram filled with 1 for every AOD
1110  fNentries->Fill(0);
1111  if(fAODProtection>=0){
1112  // Protection against different number of events in the AOD and deltaAOD
1113  // In case of discrepancy the event is rejected.
1114  Int_t matchingAODdeltaAODlevel = AliRDHFCuts::CheckMatchingAODdeltaAODevents();
1115  if (matchingAODdeltaAODlevel<0 || (matchingAODdeltaAODlevel==0 && fAODProtection==1)) {
1116  // AOD/deltaAOD trees have different number of entries || TProcessID do not match while it was required
1117  fNentries->Fill(21);
1118  return;
1119  }
1120  fNentries->Fill(22);
1121  }
1122 
1124  //fCounter->StoreEvent(aod,fReadMC);
1125  // trigger class for PbPb C0SMH-B-NOPF-ALLNOTRD, C0SMH-B-NOPF-ALL
1126  TString trigclass=aod->GetFiredTriggerClasses();
1127  if(trigclass.Contains("C0SMH-B-NOPF-ALLNOTRD") || trigclass.Contains("C0SMH-B-NOPF-ALL")) fNentries->Fill(14);
1128  if(fReadMC && fStepMCAcc){
1129  FillMCAcceptanceHistos(mcArray, mcHeader);
1130  }
1131 
1132  if(!fCuts->IsEventSelected(aod)) {
1133  if(fCuts->GetWhyRejection()==1) // rejected for pileup
1134  fNentries->Fill(13);
1135  if(fSys==1 && (fCuts->GetWhyRejection()==2 || fCuts->GetWhyRejection()==3)) fNentries->Fill(15);
1136  if(fCuts->GetWhyRejection()==7) fNentries->Fill(17);
1137  return;
1138  }
1139  // Check the Nb of SDD clusters
1140  if (fIsRejectSDDClusters) {
1141  Bool_t skipEvent = kFALSE;
1142  Int_t ntracks = 0;
1143  if (aod) ntracks = aod->GetNumberOfTracks();
1144  for(Int_t itrack=0; itrack<ntracks; itrack++) { // loop on tacks
1145  // ... get the track
1146  AliAODTrack * track = dynamic_cast<AliAODTrack*>(aod->GetTrack(itrack));
1147  if(!track) AliFatal("Not a standard AOD");
1148  if(TESTBIT(track->GetITSClusterMap(),2) || TESTBIT(track->GetITSClusterMap(),3) ){
1149  skipEvent=kTRUE;
1150  fNentries->Fill(16);
1151  break;
1152  }
1153  }
1154  if (skipEvent) return;
1155  }
1156 
1157  // AOD primary vertex
1158  AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex();
1159 
1160  Bool_t isGoodVtx=kFALSE;
1161 
1162  //vtx1->Print();
1163  TString primTitle = vtx1->GetTitle();
1164  if(primTitle.Contains("VertexerTracks") && vtx1->GetNContributors()>0) {
1165  isGoodVtx=kTRUE;
1166  fNentries->Fill(3);
1167  }
1168 
1169  // loop over candidates
1170  Int_t nInD0toKpi = inputArray->GetEntriesFast();
1171  if(fDebug>2) printf("Number of D0->Kpi: %d\n",nInD0toKpi);
1172 
1173  // FILE *f=fopen("4display.txt","a");
1174  // printf("Number of D0->Kpi: %d\n",nInD0toKpi);
1175  Int_t nSelectedloose=0,nSelectedtight=0;
1176 
1177  // vHF object is needed to call the method that refills the missing info of the candidates
1178  // if they have been deleted in dAOD reconstruction phase
1179  // in order to reduce the size of the file
1181 
1182  for (Int_t iD0toKpi = 0; iD0toKpi < nInD0toKpi; iD0toKpi++) {
1183  AliAODRecoDecayHF2Prong *d = (AliAODRecoDecayHF2Prong*)inputArray->UncheckedAt(iD0toKpi);
1184 
1186  fNentries->Fill(2);
1187  continue; //skip the D0 from Dstar
1188  }
1189  if(d->GetIsFilled()==0)fNentries->Fill(19);//tmp check
1190  if(d->GetIsFilled()==1)fNentries->Fill(20);//tmp check
1191  if(!(vHF->FillRecoCand(aod,d))) {//Fill the data members of the candidate only if they are empty.
1192  fNentries->Fill(18); //monitor how often this fails
1193  continue;
1194  }
1195 
1196  if ( fCuts->IsInFiducialAcceptance(d->Pt(),d->Y(421)) ) {
1197  nSelectedloose++;
1198  nSelectedtight++;
1199  if(fSys==0){
1200  if(fCuts->IsSelected(d,AliRDHFCuts::kTracks,aod))fNentries->Fill(6);
1201  }
1202  Int_t ptbin=fCuts->PtBin(d->Pt());
1203  if(ptbin==-1) {fNentries->Fill(4); continue;} //out of bounds
1205  if(fFillVarHists) {
1206  //if(!fCutOnDistr || (fCutOnDistr && fIsSelectedCandidate)) {
1207  fDaughterTracks.AddAt((AliAODTrack*)d->GetDaughter(0),0);
1208  fDaughterTracks.AddAt((AliAODTrack*)d->GetDaughter(1),1);
1209  //check daughters
1210  if(!fDaughterTracks.UncheckedAt(0) || !fDaughterTracks.UncheckedAt(1)) {
1211  AliDebug(1,"at least one daughter not found!");
1212  fNentries->Fill(5);
1213  fDaughterTracks.Clear();
1214  continue;
1215  }
1216  //}
1217  FillVarHists(aod,d,mcArray,fCuts,fDistr);
1218  }
1219 
1220  if (fDrawDetSignal) {
1222  }
1223  FillMassHists(d,mcArray,mcHeader,fCuts,fOutputMass);
1224  NormIPvar(aod, d,mcArray);
1225  if (fPIDCheck) {
1226  Int_t isSelectedPIDfill = 3;
1227  if (!fReadMC || (fReadMC && fUsePid4Distr)) isSelectedPIDfill = fCuts->IsSelectedPID(d); //0 rejected,1 D0,2 Dobar, 3 both
1228 
1229  if (isSelectedPIDfill == 0)fNentries->Fill(7);
1230  if (isSelectedPIDfill == 1)fNentries->Fill(8);
1231  if (isSelectedPIDfill == 2)fNentries->Fill(9);
1232  if (isSelectedPIDfill == 3)fNentries->Fill(10);
1233  }
1234 
1235  }
1236 
1237  fDaughterTracks.Clear();
1238  //if(unsetvtx) d->UnsetOwnPrimaryVtx();
1239  } //end for prongs
1240  fCounter->StoreCandidates(aod,nSelectedloose,kTRUE);
1241  fCounter->StoreCandidates(aod,nSelectedtight,kFALSE);
1242  delete vHF;
1243  // Post the data
1244  PostData(1,fOutputMass);
1245  PostData(2,fDistr);
1246  PostData(3,fNentries);
1247  PostData(5,fCounter);
1248  PostData(6,fOutputMassPt);
1249  PostData(7,fVariablesTree);
1250  PostData(8, fDetSignal);
1251 
1252  return;
1253 }
1254 
1255 //____________________________________________________________________________
1257 {
1258  //
1260  //
1261  fDaughterTracks.AddAt((AliAODTrack*)part->GetDaughter(0), 0);
1262  fDaughterTracks.AddAt((AliAODTrack*)part->GetDaughter(1), 1);
1263 
1264  AliESDtrack *esdtrack1 = new AliESDtrack((AliVTrack*)fDaughterTracks.UncheckedAt(0));
1265  AliESDtrack *esdtrack2 = new AliESDtrack((AliVTrack*)fDaughterTracks.UncheckedAt(1));
1266 
1267 
1268  //For filling detector histograms
1269  Int_t isSelectedPIDfill = 3;
1270  if (!fReadMC || (fReadMC && fUsePid4Distr)) isSelectedPIDfill = fCuts->IsSelectedPID(part); //0 rejected,1 D0,2 Dobar, 3 both
1271 
1272 
1273  //fill "before PID" histos with every daughter
1274  ((TH2F*)ListDetSignal->FindObject("TOFSigBefPID"))->Fill(esdtrack1->P(), esdtrack1->GetTOFsignal());
1275  ((TH2F*)ListDetSignal->FindObject("TOFSigBefPID"))->Fill(esdtrack2->P(), esdtrack2->GetTOFsignal());
1276  ((TH2F*)ListDetSignal->FindObject("TPCSigBefPID"))->Fill(esdtrack1->P(), esdtrack1->GetTPCsignal());
1277  ((TH2F*)ListDetSignal->FindObject("TPCSigBefPID"))->Fill(esdtrack2->P(), esdtrack2->GetTPCsignal());
1278 
1279  if (isSelectedPIDfill != 0) { //fill "After PID" with everything that's not rejected
1280  ((TH2F*)ListDetSignal->FindObject("TOFSigAftPID"))->Fill(esdtrack1->P(), esdtrack1->GetTOFsignal());
1281  ((TH2F*)ListDetSignal->FindObject("TOFSigAftPID"))->Fill(esdtrack2->P(), esdtrack2->GetTOFsignal());
1282  ((TH2F*)ListDetSignal->FindObject("TPCSigAftPID"))->Fill(esdtrack1->P(), esdtrack1->GetTPCsignal());
1283  ((TH2F*)ListDetSignal->FindObject("TPCSigAftPID"))->Fill(esdtrack2->P(), esdtrack2->GetTPCsignal());
1284 
1285  }
1286 
1287  delete esdtrack1;
1288  delete esdtrack2;
1289 
1290  return;
1291 }
1292 
1293 //____________________________________________________________________________
1295  //
1297  //
1298 
1299 
1300  Int_t pdgDgD0toKpi[2]={321,211};
1301  Int_t lab=-9999;
1302  if(fReadMC) lab=part->MatchToMC(421,arrMC,2,pdgDgD0toKpi); //return MC particle label if the array corresponds to a D0, -1 if not (cf. AliAODRecoDecay.cxx)
1303  //Double_t pt = d->Pt(); //mother pt
1304  Int_t isSelectedPID=3;
1305  if(!fReadMC || (fReadMC && fUsePid4Distr)) isSelectedPID=cuts->IsSelectedPID(part); //0 rejected,1 D0,2 Dobar, 3 both
1306  if (!fPIDCheck) { //if fPIDCheck, this is already filled elsewhere
1307  if (isSelectedPID==0)fNentries->Fill(7);
1308  if (isSelectedPID==1)fNentries->Fill(8);
1309  if (isSelectedPID==2)fNentries->Fill(9);
1310  if (isSelectedPID==3)fNentries->Fill(10);
1311  //fNentries->Fill(8+isSelectedPID);
1312  }
1313 
1314  if(fCutOnDistr && !fIsSelectedCandidate) return;
1315  //printf("\nif no cuts or cuts passed\n");
1316 
1317 
1318  //add distr here
1319  UInt_t pdgs[2];
1320 
1321  Double_t mPDG=TDatabasePDG::Instance()->GetParticle(421)->Mass();
1322  pdgs[0]=211;
1323  pdgs[1]=321;
1324  Double_t minvD0 = part->InvMassD0();
1325  pdgs[1]=211;
1326  pdgs[0]=321;
1327  Double_t minvD0bar = part->InvMassD0bar();
1328  //cout<<"inside mass cut"<<endl;
1329 
1330  Double_t invmasscut=0.03;
1331 
1332  TString fillthispi="",fillthisK="",fillthis="", fillthispt="", fillthisetaphi="";
1333 
1334  Int_t ptbin=cuts->PtBin(part->Pt());
1335  Double_t pt = part->Pt();
1336 
1337  Double_t dz1[2],dz2[2],covar1[3],covar2[3];//,d0xd0proper,errd0xd0proper;
1338  dz1[0]=-99; dz2[0]=-99;
1339  Double_t d0[2];
1340  Double_t decl[2]={-99,-99};
1341  Bool_t recalcvtx=kFALSE;
1342 
1343 
1344 
1346  recalcvtx=kTRUE;
1347  //recalculate vertex
1348  AliAODVertex *vtxProp=0x0;
1349  vtxProp=GetPrimaryVtxSkipped(aod);
1350  if(vtxProp) {
1351  part->SetOwnPrimaryVtx(vtxProp);
1352  //Bool_t unsetvtx=kTRUE;
1353  //Calculate d0 for daughter tracks with Vtx Skipped
1354  AliESDtrack *esdtrack1=new AliESDtrack((AliVTrack*)fDaughterTracks.UncheckedAt(0));
1355  AliESDtrack *esdtrack2=new AliESDtrack((AliVTrack*)fDaughterTracks.UncheckedAt(1));
1356  esdtrack1->PropagateToDCA(vtxProp,aod->GetMagneticField(),1.,dz1,covar1);
1357  esdtrack2->PropagateToDCA(vtxProp,aod->GetMagneticField(),1.,dz2,covar2);
1358  delete vtxProp; vtxProp=NULL;
1359  delete esdtrack1;
1360  delete esdtrack2;
1361  }
1362 
1363  d0[0]=dz1[0];
1364  d0[1]=dz2[0];
1365 
1366  decl[0]=part->DecayLength2();
1367  decl[1]=part->NormalizedDecayLength2();
1368  part->UnsetOwnPrimaryVtx();
1369 
1370  }
1371 
1372  Double_t cosThetaStarD0 = 99;
1373  Double_t cosThetaStarD0bar = 99;
1374  Double_t cosPointingAngle = 99;
1375  Double_t normalizedDecayLength2 = -1, normalizedDecayLengthxy=-1;
1376  Double_t decayLength2 = -1, decayLengthxy=-1;
1377  Double_t ptProng[2]={-99,-99};
1378  Double_t d0Prong[2]={-99,-99};
1379  Double_t etaD = 99.;
1380  Double_t phiD = 99.;
1381 
1382 
1383  //disable the PID
1384  if(!fUsePid4Distr) isSelectedPID=0;
1385  if((lab>=0 && fReadMC) || (!fReadMC && isSelectedPID)){ //signal (from MC or PID)
1386 
1387  //check pdg of the prongs
1388  AliAODTrack *prong0=(AliAODTrack*)fDaughterTracks.UncheckedAt(0);
1389  AliAODTrack *prong1=(AliAODTrack*)fDaughterTracks.UncheckedAt(1);
1390 
1391  if(!prong0 || !prong1) {
1392  return;
1393  }
1394  Int_t labprong[2];
1395  if(fReadMC){
1396  labprong[0]=prong0->GetLabel();
1397  labprong[1]=prong1->GetLabel();
1398  }
1399  AliAODMCParticle *mcprong=0;
1400  Int_t pdgProng[2]={0,0};
1401 
1402  for (Int_t iprong=0;iprong<2;iprong++){
1403  if(fReadMC && labprong[iprong]>=0) {
1404  mcprong= (AliAODMCParticle*)arrMC->At(labprong[iprong]);
1405  pdgProng[iprong]=mcprong->GetPdgCode();
1406  }
1407  }
1408 
1409  if(fSys==0){
1410  //no mass cut ditributions: ptbis
1411 
1412  fillthispi="hptpiSnoMcut_";
1413  fillthispi+=ptbin;
1414 
1415  fillthisK="hptKSnoMcut_";
1416  fillthisK+=ptbin;
1417 
1418  if ((TMath::Abs(pdgProng[0]) == 211 && TMath::Abs(pdgProng[1]) == 321)
1419  || (isSelectedPID==1 || isSelectedPID==3)){
1420  ((TH1F*)listout->FindObject(fillthispi))->Fill(prong0->Pt());
1421  ((TH1F*)listout->FindObject(fillthisK))->Fill(prong1->Pt());
1422  }
1423 
1424  if ((TMath::Abs(pdgProng[0]) == 321 && TMath::Abs(pdgProng[1]) == 211)
1425  || (isSelectedPID==2 || isSelectedPID==3)){
1426  ((TH1F*)listout->FindObject(fillthisK))->Fill(prong0->Pt());
1427  ((TH1F*)listout->FindObject(fillthispi))->Fill(prong1->Pt());
1428  }
1429  }
1430 
1431  //no mass cut ditributions: mass
1432 
1433  etaD = part->Eta();
1434  phiD = part->Phi();
1435 
1436 
1437  fillthis="hMassS_";
1438  fillthis+=ptbin;
1439  fillthispt="histSgnPt";
1440 
1441  if ((fReadMC && ((AliAODMCParticle*)arrMC->At(lab))->GetPdgCode() == 421)
1442  || (!fReadMC && (isSelectedPID==1 || isSelectedPID==3))){//D0
1443  ((TH1F*)listout->FindObject(fillthis))->Fill(minvD0);
1444  if(fFillPtHist && fReadMC) ((TH2F*)fOutputMassPt->FindObject(fillthispt))->Fill(minvD0,pt);
1445 
1446  fillthisetaphi="hetaphiD0candidateS_";
1447  fillthisetaphi+=ptbin;
1448  ((TH2F*)listout->FindObject(fillthisetaphi))->Fill(etaD, phiD);
1449 
1450  if(TMath::Abs(minvD0-mPDG)<0.05){
1451  fillthisetaphi="hetaphiD0candidatesignalregionS_";
1452  fillthisetaphi+=ptbin;
1453  ((TH2F*)listout->FindObject(fillthisetaphi))->Fill(etaD, phiD);
1454  }
1455 
1456  }
1457  else { //D0bar
1458  if(fReadMC || (!fReadMC && isSelectedPID > 1)){
1459  ((TH1F*)listout->FindObject(fillthis))->Fill(minvD0bar);
1460  if(fFillPtHist && fReadMC) ((TH2F*)fOutputMassPt->FindObject(fillthispt))->Fill(minvD0bar,pt);
1461 
1462  fillthisetaphi="hetaphiD0barcandidateS_";
1463  fillthisetaphi+=ptbin;
1464  ((TH2F*)listout->FindObject(fillthisetaphi))->Fill(etaD, phiD);
1465 
1466  if(TMath::Abs(minvD0bar-mPDG)<0.05){
1467  fillthisetaphi="hetaphiD0barcandidatesignalregionS_";
1468  fillthisetaphi+=ptbin;
1469  ((TH2F*)listout->FindObject(fillthisetaphi))->Fill(etaD, phiD);
1470  }
1471 
1472  }
1473  }
1474 
1475  //apply cut on invariant mass on the pair
1476  if(TMath::Abs(minvD0-mPDG)<invmasscut || TMath::Abs(minvD0bar-mPDG)<invmasscut){
1477 
1478  cosThetaStarD0 = part->CosThetaStarD0();
1479  cosThetaStarD0bar = part->CosThetaStarD0bar();
1480  cosPointingAngle = part->CosPointingAngle();
1481  normalizedDecayLength2 = part->NormalizedDecayLength2();
1482  decayLength2 = part->DecayLength2();
1483  decayLengthxy = part->DecayLengthXY();
1484  normalizedDecayLengthxy=decayLengthxy/part->DecayLengthXYError();
1485 
1486  ptProng[0]=prong0->Pt(); ptProng[1]=prong1->Pt();
1487  d0Prong[0]=part->Getd0Prong(0); d0Prong[1]=part->Getd0Prong(1);
1488 
1489  if(fArray==1) cout<<"LS signal: ERROR"<<endl;
1490  for (Int_t iprong=0; iprong<2; iprong++){
1491  AliAODTrack *prong=(AliAODTrack*)fDaughterTracks.UncheckedAt(iprong);
1492  if (fReadMC) labprong[iprong]=prong->GetLabel();
1493 
1494  //cout<<"prong name = "<<prong->GetName()<<" label = "<<prong->GetLabel()<<endl;
1495  Int_t pdgprong=0;
1496  if(fReadMC && labprong[iprong]>=0) {
1497  mcprong= (AliAODMCParticle*)arrMC->At(labprong[iprong]);
1498  pdgprong=mcprong->GetPdgCode();
1499  }
1500 
1501  Bool_t isPionHere[2]={(isSelectedPID==1 || isSelectedPID==3) ? kTRUE : kFALSE,(isSelectedPID==2 || isSelectedPID==3) ? kTRUE : kFALSE};
1502 
1503  if(TMath::Abs(pdgprong)==211 || isPionHere[iprong]) {
1504  //cout<<"pi"<<endl;
1505 
1506  if(fSys==0){
1507  fillthispi="hptpiS_";
1508  fillthispi+=ptbin;
1509  ((TH1F*)listout->FindObject(fillthispi))->Fill(ptProng[iprong]);
1510  }
1511 
1512  fillthispi="hd0piS_";
1513  fillthispi+=ptbin;
1514  ((TH1F*)listout->FindObject(fillthispi))->Fill(d0Prong[iprong]);
1515  if(recalcvtx) {
1516 
1517  fillthispi="hd0vpiS_";
1518  fillthispi+=ptbin;
1519  ((TH1F*)listout->FindObject(fillthispi))->Fill(d0[iprong]);
1520  }
1521 
1522  }
1523 
1524  if(TMath::Abs(pdgprong)==321 || !isPionHere[iprong]) {
1525  //cout<<"kappa"<<endl;
1526  if(fSys==0){
1527  fillthisK="hptKS_";
1528  fillthisK+=ptbin;
1529  ((TH1F*)listout->FindObject(fillthisK))->Fill(ptProng[iprong]);
1530  }
1531 
1532 
1533  fillthisK="hd0KS_";
1534  fillthisK+=ptbin;
1535  ((TH1F*)listout->FindObject(fillthisK))->Fill(d0Prong[iprong]);
1536  if (recalcvtx){
1537  fillthisK="hd0vKS_";
1538  fillthisK+=ptbin;
1539  ((TH1F*)listout->FindObject(fillthisK))->Fill(d0[iprong]);
1540  }
1541  }
1542 
1543  if(fSys==0){
1544  fillthis="hcosthpointd0S_";
1545  fillthis+=ptbin;
1546  ((TH1F*)listout->FindObject(fillthis))->Fill(cosPointingAngle,d0Prong[iprong]);
1547  }
1548  } //end loop on prongs
1549 
1550  fillthis="hdcaS_";
1551  fillthis+=ptbin;
1552  ((TH1F*)listout->FindObject(fillthis))->Fill(part->GetDCA());
1553 
1554  fillthis="hcosthetapointS_";
1555  fillthis+=ptbin;
1556  ((TH1F*)listout->FindObject(fillthis))->Fill(cosPointingAngle);
1557 
1558  fillthis="hcosthetapointxyS_";
1559  fillthis+=ptbin;
1560  ((TH1F*)listout->FindObject(fillthis))->Fill(part->CosPointingAngleXY());
1561 
1562 
1563  fillthis="hd0d0S_";
1564  fillthis+=ptbin;
1565  ((TH1F*)listout->FindObject(fillthis))->Fill(part->Prodd0d0());
1566 
1567  fillthis="hdeclS_";
1568  fillthis+=ptbin;
1569  ((TH1F*)listout->FindObject(fillthis))->Fill(decayLength2);
1570 
1571  fillthis="hnormdeclS_";
1572  fillthis+=ptbin;
1573  ((TH1F*)listout->FindObject(fillthis))->Fill(normalizedDecayLength2);
1574 
1575  fillthis="hdeclxyS_";
1576  fillthis+=ptbin;
1577  ((TH1F*)listout->FindObject(fillthis))->Fill(decayLengthxy);
1578 
1579  fillthis="hnormdeclxyS_";
1580  fillthis+=ptbin;
1581  ((TH1F*)listout->FindObject(fillthis))->Fill(normalizedDecayLengthxy);
1582 
1583  fillthis="hdeclxyd0d0S_";
1584  fillthis+=ptbin;
1585  ((TH2F*)listout->FindObject(fillthis))->Fill(decayLengthxy,d0Prong[0]*d0Prong[1]);
1586 
1587  fillthis="hnormdeclxyd0d0S_";
1588  fillthis+=ptbin;
1589  ((TH2F*)listout->FindObject(fillthis))->Fill(normalizedDecayLengthxy,d0Prong[0]*d0Prong[1]);
1590 
1591  if(recalcvtx) {
1592  fillthis="hdeclvS_";
1593  fillthis+=ptbin;
1594  ((TH1F*)listout->FindObject(fillthis))->Fill(decl[0]);
1595 
1596  fillthis="hnormdeclvS_";
1597  fillthis+=ptbin;
1598  ((TH1F*)listout->FindObject(fillthis))->Fill(decl[1]);
1599 
1600  fillthis="hd0d0vS_";
1601  fillthis+=ptbin;
1602  ((TH1F*)listout->FindObject(fillthis))->Fill(d0[0]*d0[1]);
1603  }
1604 
1605  if(fSys==0){
1606  fillthis="hcosthetastarS_";
1607  fillthis+=ptbin;
1608  if ((fReadMC && ((AliAODMCParticle*)arrMC->At(lab))->GetPdgCode() == 421)) ((TH1F*)listout->FindObject(fillthis))->Fill(cosThetaStarD0);
1609  else {
1610  if (fReadMC || isSelectedPID>1)((TH1F*)listout->FindObject(fillthis))->Fill(cosThetaStarD0bar);
1611  if(isSelectedPID==1 || isSelectedPID==3)((TH1F*)listout->FindObject(fillthis))->Fill(cosThetaStarD0);
1612  }
1613 
1614  fillthis="hcosthpointd0d0S_";
1615  fillthis+=ptbin;
1616  ((TH2F*)listout->FindObject(fillthis))->Fill(cosPointingAngle,part->Prodd0d0());
1617  }
1618 
1619  if ((fReadMC && ((AliAODMCParticle*)arrMC->At(lab))->GetPdgCode() == 421)){
1620  for(Int_t it=0; it<2; it++){
1621  fillthis="hptD0S_";
1622  fillthis+=ptbin;
1623  ((TH1F*)listout->FindObject(fillthis))->Fill(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->Pt());
1624  fillthis="hphiD0S_";
1625  fillthis+=ptbin;
1626  ((TH1F*)listout->FindObject(fillthis))->Fill(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->Phi());
1627  Int_t nPointsITS = 0;
1628  for (Int_t il=0; il<6; il++){
1629  if(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->HasPointOnITSLayer(il)) nPointsITS++;
1630  }
1631  fillthis="hNITSpointsD0vsptS_";
1632  fillthis+=ptbin;
1633  ((TH2F*)listout->FindObject(fillthis))->Fill(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->Pt(),nPointsITS);
1634  fillthis="hNSPDpointsD0S_";
1635  fillthis+=ptbin;
1636  if(!(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->HasPointOnITSLayer(0)) && !(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->HasPointOnITSLayer(1))){ //no SPD points
1637  ((TH1I*)listout->FindObject(fillthis))->Fill(0);
1638  }
1639  if(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->HasPointOnITSLayer(0) && !(((AliAODTrack*)(fDaughterTracks.UncheckedAt(it)))->HasPointOnITSLayer(1))){ //kOnlyFirst
1640  ((TH1I*)listout->FindObject(fillthis))->Fill(1);
1641  }
1642  if(!(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->HasPointOnITSLayer(0)) && ((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->HasPointOnITSLayer(1)){ //kOnlySecond
1643  ((TH1I*)listout->FindObject(fillthis))->Fill(2);
1644  }
1645  if(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->HasPointOnITSLayer(0) && ((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->HasPointOnITSLayer(1)){ //kboth
1646  ((TH1I*)listout->FindObject(fillthis))->Fill(3);
1647  }
1648  fillthis="hNclsD0vsptS_";
1649  fillthis+=ptbin;
1650  Float_t mom = ((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->Pt();
1651  Float_t ncls = (Float_t)((AliAODTrack*)fDaughterTracks.UncheckedAt(0))->GetTPCNcls();
1652  ((TH2F*)listout->FindObject(fillthis))->Fill(mom, ncls);
1653  }
1654  }
1655  else {
1656  if (fReadMC || isSelectedPID>1){
1657  for(Int_t it=0; it<2; it++){
1658  fillthis="hptD0barS_";
1659  fillthis+=ptbin;
1660  ((TH1F*)listout->FindObject(fillthis))->Fill(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->Pt());
1661  fillthis="hphiD0barS_";
1662  fillthis+=ptbin;
1663  ((TH1F*)listout->FindObject(fillthis))->Fill(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->Phi());
1664  fillthis="hNclsD0barvsptS_";
1665  fillthis+=ptbin;
1666  Float_t mom = ((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->Pt();
1667  Float_t ncls = (Float_t)((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->GetTPCNcls();
1668  ((TH2F*)listout->FindObject(fillthis))->Fill(mom, ncls);
1669  }
1670  }
1671  if(isSelectedPID==1 || isSelectedPID==3){
1672  for(Int_t it=0; it<2; it++){
1673  fillthis="hptD0S_";
1674  fillthis+=ptbin;
1675  ((TH1F*)listout->FindObject(fillthis))->Fill(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->Pt());
1676  fillthis="hphiD0S_";
1677  fillthis+=ptbin;
1678  ((TH1F*)listout->FindObject(fillthis))->Fill(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->Phi());
1679  Int_t nPointsITS = 0;
1680  for (Int_t il=0; il<6; il++){
1681  if(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->HasPointOnITSLayer(il)) nPointsITS++;
1682  }
1683  fillthis="hNITSpointsD0vsptS_";
1684  fillthis+=ptbin;
1685  ((TH2F*)listout->FindObject(fillthis))->Fill(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->Pt(), nPointsITS);
1686  fillthis="hNSPDpointsD0S_";
1687  fillthis+=ptbin;
1688  if(!(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->HasPointOnITSLayer(0)) && !(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->HasPointOnITSLayer(1))){ //no SPD points
1689  ((TH1I*)listout->FindObject(fillthis))->Fill(0);
1690  }
1691  if(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->HasPointOnITSLayer(0) && !(((AliAODTrack*)(fDaughterTracks.UncheckedAt(it)))->HasPointOnITSLayer(1))){ //kOnlyFirst
1692  ((TH1I*)listout->FindObject(fillthis))->Fill(1);
1693  }
1694  if(!(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->HasPointOnITSLayer(0)) && ((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->HasPointOnITSLayer(1)){ //kOnlySecond
1695  ((TH1I*)listout->FindObject(fillthis))->Fill(2);
1696  }
1697  if(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->HasPointOnITSLayer(0) && ((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->HasPointOnITSLayer(1)){ //kboth
1698  ((TH1I*)listout->FindObject(fillthis))->Fill(3);
1699  }
1700  fillthis="hNclsD0vsptS_";
1701  fillthis+=ptbin;
1702  Float_t mom = ((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->Pt();
1703  Float_t ncls = (Float_t)((AliAODTrack*)fDaughterTracks.UncheckedAt(0))->GetTPCNcls();
1704  ((TH2F*)listout->FindObject(fillthis))->Fill(mom, ncls);
1705  }
1706  }
1707  }
1708 
1709 
1710  } //end mass cut
1711 
1712  } else{ //Background or LS
1713  //if(!fReadMC){
1714  //cout<<"is background"<<endl;
1715 
1716  etaD = part->Eta();
1717  phiD = part->Phi();
1718 
1719  //no mass cut distributions: mass, ptbis
1720  fillthis="hMassB_";
1721  fillthis+=ptbin;
1722  fillthispt="histBkgPt";
1723 
1725  ((TH1F*)listout->FindObject(fillthis))->Fill(minvD0);
1726  if(fFillPtHist && fReadMC) ((TH2F*)fOutputMassPt->FindObject(fillthispt))->Fill(minvD0,pt);
1727 
1728  fillthisetaphi="hetaphiD0candidateB_";
1729  fillthisetaphi+=ptbin;
1730  ((TH2F*)listout->FindObject(fillthisetaphi))->Fill(etaD, phiD);
1731 
1732  if(TMath::Abs(minvD0-mPDG)<0.05){
1733  fillthisetaphi="hetaphiD0candidatesignalregionB_";
1734  fillthisetaphi+=ptbin;
1735  ((TH2F*)listout->FindObject(fillthisetaphi))->Fill(etaD, phiD);
1736  }
1737  }
1738  if (!fCutOnDistr || (fCutOnDistr && fIsSelectedCandidate>1)) {
1739  ((TH1F*)listout->FindObject(fillthis))->Fill(minvD0bar);
1740  if(fFillPtHist && fReadMC) ((TH2F*)fOutputMassPt->FindObject(fillthispt))->Fill(minvD0bar,pt);
1741 
1742  fillthisetaphi="hetaphiD0barcandidateB_";
1743  fillthisetaphi+=ptbin;
1744  ((TH2F*)listout->FindObject(fillthisetaphi))->Fill(etaD, phiD);
1745 
1746  if(TMath::Abs(minvD0bar-mPDG)<0.05){
1747  fillthisetaphi="hetaphiD0barcandidatesignalregionB_";
1748  fillthisetaphi+=ptbin;
1749  ((TH2F*)listout->FindObject(fillthisetaphi))->Fill(etaD, phiD);
1750  }
1751 
1752  }
1753  if(fSys==0){
1754  fillthis="hptB1prongnoMcut_";
1755  fillthis+=ptbin;
1756 
1757  ((TH1F*)listout->FindObject(fillthis))->Fill(((AliAODTrack*)fDaughterTracks.UncheckedAt(0))->Pt());
1758 
1759  fillthis="hptB2prongsnoMcut_";
1760  fillthis+=ptbin;
1761  ((TH1F*)listout->FindObject(fillthis))->Fill(((AliAODTrack*)fDaughterTracks.UncheckedAt(0))->Pt());
1762  ((TH1F*)listout->FindObject(fillthis))->Fill(((AliAODTrack*)fDaughterTracks.UncheckedAt(0))->Pt());
1763  }
1764 
1765 
1766  //apply cut on invariant mass on the pair
1767  if(TMath::Abs(minvD0-mPDG)<invmasscut || TMath::Abs(minvD0bar-mPDG)<invmasscut){
1768  if(fSys==0){
1769  ptProng[0]=((AliAODTrack*)fDaughterTracks.UncheckedAt(0))->Pt(); ptProng[1]=((AliAODTrack*)fDaughterTracks.UncheckedAt(0))->Pt();
1770  cosThetaStarD0 = part->CosThetaStarD0();
1771  cosThetaStarD0bar = part->CosThetaStarD0bar();
1772  }
1773 
1774  cosPointingAngle = part->CosPointingAngle();
1775  normalizedDecayLength2 = part->NormalizedDecayLength2();
1776  decayLength2 = part->DecayLength2();
1777  decayLengthxy = part->DecayLengthXY();
1778  normalizedDecayLengthxy=decayLengthxy/part->DecayLengthXYError();
1779  d0Prong[0]=part->Getd0Prong(0); d0Prong[1]=part->Getd0Prong(1);
1780 
1781 
1782  AliAODTrack *prongg=(AliAODTrack*)fDaughterTracks.UncheckedAt(0);
1783  if(!prongg) {
1784  if(fDebug>2) cout<<"No daughter found";
1785  return;
1786  }
1787  else{
1788  if(fArray==1){
1789  if(prongg->Charge()==1) {
1790  //fTotPosPairs[ptbin]++;
1791  ((TH1F*)fOutputMass->FindObject("hpospair"))->Fill(ptbin);
1792  } else {
1793  //fTotNegPairs[ptbin]++;
1794  ((TH1F*)fOutputMass->FindObject("hnegpair"))->Fill(ptbin);
1795  }
1796  }
1797  }
1798 
1799  //fill pt and phi distrib for prongs with M cut
1800 
1802  for(Int_t it=0; it<2; it++){
1803  fillthis="hptD0B_";
1804  fillthis+=ptbin;
1805  ((TH1F*)listout->FindObject(fillthis))->Fill(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->Pt());
1806  fillthis="hphiD0B_";
1807  fillthis+=ptbin;
1808  ((TH1F*)listout->FindObject(fillthis))->Fill(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->Phi());
1809 
1810  Int_t nPointsITS = 0;
1811  for (Int_t il=0; il<6; il++){
1812  if(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->HasPointOnITSLayer(il)) nPointsITS++;
1813  }
1814  fillthis="hNITSpointsD0vsptB_";
1815  fillthis+=ptbin;
1816  ((TH2F*)listout->FindObject(fillthis))->Fill(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->Pt(), nPointsITS);
1817  fillthis="hNSPDpointsD0B_";
1818  fillthis+=ptbin;
1819  if(!(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->HasPointOnITSLayer(0)) && !(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->HasPointOnITSLayer(1))){ //no SPD points
1820  ((TH1I*)listout->FindObject(fillthis))->Fill(0);
1821  }
1822  if(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->HasPointOnITSLayer(0) && !(((AliAODTrack*)(fDaughterTracks.UncheckedAt(it)))->HasPointOnITSLayer(1))){ //kOnlyFirst
1823  ((TH1I*)listout->FindObject(fillthis))->Fill(1);
1824  }
1825  if(!(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->HasPointOnITSLayer(0)) && ((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->HasPointOnITSLayer(1)){ //kOnlySecond
1826  ((TH1I*)listout->FindObject(fillthis))->Fill(2);
1827  }
1828  if(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->HasPointOnITSLayer(0) && ((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->HasPointOnITSLayer(1)){ //kboth
1829  ((TH1I*)listout->FindObject(fillthis))->Fill(3);
1830  }
1831  fillthis="hNclsD0vsptB_";
1832  fillthis+=ptbin;
1833  Float_t mom = ((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->Pt();
1834  Float_t ncls = (Float_t)((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->GetTPCNcls();
1835  ((TH2F*)listout->FindObject(fillthis))->Fill(mom, ncls);
1836  }
1837 
1838 
1839  }
1840 
1841  if (!fCutOnDistr || (fCutOnDistr && fIsSelectedCandidate>1)) {
1842  for(Int_t it=0; it<2; it++){
1843  fillthis="hptD0barB_";
1844  fillthis+=ptbin;
1845  ((TH1F*)listout->FindObject(fillthis))->Fill(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->Pt());
1846  fillthis="hphiD0barB_";
1847  fillthis+=ptbin;
1848  ((TH1F*)listout->FindObject(fillthis))->Fill(((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->Phi());
1849  fillthis="hNclsD0barvsptB_";
1850  fillthis+=ptbin;
1851  Float_t mom = ((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->Pt();
1852  Float_t ncls = (Float_t)((AliAODTrack*)fDaughterTracks.UncheckedAt(it))->GetTPCNcls();
1853  ((TH2F*)listout->FindObject(fillthis))->Fill(mom, ncls);
1854  }
1855  }
1856 
1857  fillthis="hd0B_";
1858  fillthis+=ptbin;
1859  ((TH1F*)listout->FindObject(fillthis))->Fill(d0Prong[0]);
1860  ((TH1F*)listout->FindObject(fillthis))->Fill(d0Prong[1]);
1861 
1862  if(fReadMC){
1863  Int_t pdgMother[2]={0,0};
1864  Double_t factor[2]={1,1};
1865 
1866  for(Int_t iprong=0;iprong<2;iprong++){
1867  AliAODTrack *prong=(AliAODTrack*)fDaughterTracks.UncheckedAt(iprong);
1868  lab=prong->GetLabel();
1869  if(lab>=0){
1870  AliAODMCParticle* mcprong=(AliAODMCParticle*)arrMC->At(lab);
1871  if(mcprong){
1872  Int_t labmom=mcprong->GetMother();
1873  if(labmom>=0){
1874  AliAODMCParticle* mcmother=(AliAODMCParticle*)arrMC->At(labmom);
1875  if(mcmother) pdgMother[iprong]=mcmother->GetPdgCode();
1876  }
1877  }
1878  }
1879 
1880  if(fSys==0){
1881 
1882  fillthis="hd0moresB_";
1883  fillthis+=ptbin;
1884 
1885  if(TMath::Abs(pdgMother[iprong])==310 || TMath::Abs(pdgMother[iprong])==130 || TMath::Abs(pdgMother[iprong])==321){ //K^0_S, K^0_L, K^+-
1886  if(ptProng[iprong]<=1)factor[iprong]=1./.7;
1887  else factor[iprong]=1./.6;
1888  fNentries->Fill(11);
1889  }
1890 
1891  if(TMath::Abs(pdgMother[iprong])==3122) { //Lambda
1892  factor[iprong]=1./0.25;
1893  fNentries->Fill(12);
1894  }
1895  fillthis="hd0moresB_";
1896  fillthis+=ptbin;
1897 
1898  ((TH1F*)listout->FindObject(fillthis))->Fill(d0Prong[iprong],factor[iprong]);
1899 
1900  if(recalcvtx){
1901  fillthis="hd0vmoresB_";
1902  fillthis+=ptbin;
1903  ((TH1F*)listout->FindObject(fillthis))->Fill(d0[iprong],factor[iprong]);
1904  }
1905  }
1906  } //loop on prongs
1907 
1908  if(fSys==0){
1909  fillthis="hd0d0moresB_";
1910  fillthis+=ptbin;
1911  ((TH1F*)listout->FindObject(fillthis))->Fill(part->Prodd0d0(),factor[0]*factor[1]);
1912 
1913  fillthis="hcosthetapointmoresB_";
1914  fillthis+=ptbin;
1915  ((TH1F*)listout->FindObject(fillthis))->Fill(cosPointingAngle,factor[0]*factor[1]);
1916 
1917  if(recalcvtx){
1918  fillthis="hd0d0vmoresB_";
1919  fillthis+=ptbin;
1920  ((TH1F*)listout->FindObject(fillthis))->Fill(d0[0]*d0[1],factor[0]*factor[1]);
1921  }
1922  }
1923  } //readMC
1924 
1925  if(fSys==0){
1926  //normalise pt distr to half afterwards
1927  fillthis="hptB_";
1928  fillthis+=ptbin;
1929  ((TH1F*)listout->FindObject(fillthis))->Fill(ptProng[0]);
1930  ((TH1F*)listout->FindObject(fillthis))->Fill(ptProng[1]);
1931 
1932  fillthis="hcosthetastarB_";
1933  fillthis+=ptbin;
1934  if (!fCutOnDistr || (fCutOnDistr && (fIsSelectedCandidate==1 || fIsSelectedCandidate==3)))((TH1F*)listout->FindObject(fillthis))->Fill(cosThetaStarD0);
1935  if (!fCutOnDistr || (fCutOnDistr && fIsSelectedCandidate>1))((TH1F*)listout->FindObject(fillthis))->Fill(cosThetaStarD0bar);
1936 
1937 
1938  fillthis="hd0p0B_";
1939  fillthis+=ptbin;
1940  ((TH1F*)listout->FindObject(fillthis))->Fill(d0Prong[0]);
1941  fillthis="hd0p1B_";
1942  fillthis+=ptbin;
1943  ((TH1F*)listout->FindObject(fillthis))->Fill(d0Prong[1]);
1944 
1945  fillthis="hcosthpointd0d0B_";
1946  fillthis+=ptbin;
1947  ((TH2F*)listout->FindObject(fillthis))->Fill(cosPointingAngle,part->Prodd0d0());
1948 
1949  fillthis="hcosthpointd0B_";
1950  fillthis+=ptbin;
1951  ((TH1F*)listout->FindObject(fillthis))->Fill(cosPointingAngle,d0Prong[0]);
1952  ((TH1F*)listout->FindObject(fillthis))->Fill(cosPointingAngle,d0Prong[1]);
1953 
1954 
1955  if(recalcvtx){
1956 
1957  fillthis="hd0vp0B_";
1958  fillthis+=ptbin;
1959  ((TH1F*)listout->FindObject(fillthis))->Fill(d0[0]);
1960  fillthis="hd0vp1B_";
1961  fillthis+=ptbin;
1962  ((TH1F*)listout->FindObject(fillthis))->Fill(d0[1]);
1963 
1964  fillthis="hd0vB_";
1965  fillthis+=ptbin;
1966  ((TH1F*)listout->FindObject(fillthis))->Fill(d0[0]);
1967  ((TH1F*)listout->FindObject(fillthis))->Fill(d0[1]);
1968 
1969  }
1970 
1971  }
1972 
1973  fillthis="hdcaB_";
1974  fillthis+=ptbin;
1975  ((TH1F*)listout->FindObject(fillthis))->Fill(part->GetDCA());
1976 
1977  fillthis="hd0d0B_";
1978  fillthis+=ptbin;
1979  ((TH1F*)listout->FindObject(fillthis))->Fill(d0Prong[0]*d0Prong[1]);
1980 
1981  if(recalcvtx){
1982  fillthis="hd0d0vB_";
1983  fillthis+=ptbin;
1984  ((TH1F*)listout->FindObject(fillthis))->Fill(d0[0]*d0[1]);
1985  }
1986 
1987  fillthis="hcosthetapointB_";
1988  fillthis+=ptbin;
1989  ((TH1F*)listout->FindObject(fillthis))->Fill(cosPointingAngle);
1990 
1991  fillthis="hcosthetapointxyB_";
1992  fillthis+=ptbin;
1993  ((TH1F*)listout->FindObject(fillthis))->Fill(part->CosPointingAngleXY());
1994 
1995  fillthis="hdeclB_";
1996  fillthis+=ptbin;
1997  ((TH1F*)listout->FindObject(fillthis))->Fill(decayLength2);
1998 
1999  fillthis="hnormdeclB_";
2000  fillthis+=ptbin;
2001  ((TH1F*)listout->FindObject(fillthis))->Fill(normalizedDecayLength2);
2002 
2003  fillthis="hdeclxyB_";
2004  fillthis+=ptbin;
2005  ((TH1F*)listout->FindObject(fillthis))->Fill(decayLengthxy);
2006 
2007  fillthis="hnormdeclxyB_";
2008  fillthis+=ptbin;
2009  ((TH1F*)listout->FindObject(fillthis))->Fill(normalizedDecayLengthxy);
2010 
2011  fillthis="hdeclxyd0d0B_";
2012  fillthis+=ptbin;
2013  ((TH2F*)listout->FindObject(fillthis))->Fill(decayLengthxy,d0Prong[0]*d0Prong[1]);
2014 
2015  fillthis="hnormdeclxyd0d0B_";
2016  fillthis+=ptbin;
2017  ((TH2F*)listout->FindObject(fillthis))->Fill(normalizedDecayLengthxy,d0Prong[0]*d0Prong[1]);
2018 
2019 
2020  if(recalcvtx) {
2021 
2022  fillthis="hdeclvB_";
2023  fillthis+=ptbin;
2024  ((TH1F*)listout->FindObject(fillthis))->Fill(decl[0]);
2025 
2026  fillthis="hnormdeclvB_";
2027  fillthis+=ptbin;
2028  ((TH1F*)listout->FindObject(fillthis))->Fill(decl[1]);
2029 
2030 
2031  }
2032  }//mass cut
2033  }//else (background)
2034 
2035  return;
2036 }
2037 
2038 //____________________________________________________________________________
2039 void AliAnalysisTaskSED0Mass::FillMassHists(AliAODRecoDecayHF2Prong *part, TClonesArray *arrMC, AliAODMCHeader *mcHeader, AliRDHFCutsD0toKpi* cuts, TList *listout){
2040  //
2042  //
2043 
2044 
2045  Double_t mPDG=TDatabasePDG::Instance()->GetParticle(421)->Mass();
2046 
2047  //cout<<"is selected = "<<fIsSelectedCandidate<<endl;
2048 
2049  // Fill candidate variable Tree (track selection, no candidate selection)
2050  if( fWriteVariableTree && !part->HasBadDaughters()
2051  && fCuts->AreDaughtersSelected(part) && fCuts->IsSelectedPID(part) ){
2052  fCandidateVariables[0] = part->InvMassD0();
2053  fCandidateVariables[1] = part->InvMassD0bar();
2054  fCandidateVariables[2] = part->Pt();
2055  fCandidateVariables[3] = part->GetDCA();
2056  Double_t ctsD0=0. ,ctsD0bar=0.; part->CosThetaStarD0(ctsD0,ctsD0bar);
2057  fCandidateVariables[4] = ctsD0;
2058  fCandidateVariables[5] = ctsD0bar;
2059  fCandidateVariables[6] = part->Pt2Prong(0);
2060  fCandidateVariables[7] = part->Pt2Prong(1);
2061  fCandidateVariables[8] = part->Getd0Prong(0);
2062  fCandidateVariables[9] = part->Getd0Prong(1);
2063  fCandidateVariables[10] = part->Prodd0d0();
2064  fCandidateVariables[11] = part->CosPointingAngle();
2068  fVariablesTree->Fill();
2069  }
2070 
2071  //cout<<"check cuts = "<<endl;
2072  //cuts->PrintAll();
2073  if (!fIsSelectedCandidate){
2074  //cout<<" cut " << cuts << " Rejected because "<<cuts->GetWhy()<<endl;
2075  return;
2076  }
2077 
2078 
2079  if(fDebug>2) cout<<"Candidate selected"<<endl;
2080 
2081  Double_t invmassD0 = part->InvMassD0(), invmassD0bar = part->InvMassD0bar();
2082  //printf("SELECTED\n");
2083  Int_t ptbin=cuts->PtBin(part->Pt());
2084  Double_t pt = part->Pt();
2085  Double_t y = part->YD0();
2086 
2087  Double_t impparXY=part->ImpParXY()*10000.;
2088  Double_t trueImpParXY=0.;
2089  Double_t arrayForSparse[3]={invmassD0,pt,impparXY};
2090  Double_t arrayForSparseTrue[3]={invmassD0,pt,trueImpParXY};
2091 
2092 
2093  // AliAODTrack *prong=(AliAODTrack*)fDaughterTracks.UncheckedAt(0);
2094  // if(!prong) {
2095  // AliDebug(2,"No daughter found");
2096  // return;
2097  // }
2098  // else{
2099  // if(prong->Charge()==1) {
2100  // ((TH1F*)fDistr->FindObject("hpospair"))->Fill(fCuts->GetNPtBins()+ptbin);
2101  // //fTotPosPairs[ptbin]++;
2102  // } else {
2103  // ((TH1F*)fDistr->FindObject("hnegpair"))->Fill(fCuts->GetNPtBins()+ptbin);
2104  // //fTotNegPairs[ptbin]++;
2105  // }
2106  // }
2107 
2108  // for(Int_t it=0;it<2;it++){
2109 
2110  // //request on spd points to be addes
2111  // if(/*nSPD==2 && */part->Pt() > 5. && (TMath::Abs(invmassD0-mPDG)<0.01 || TMath::Abs(invmassD0bar-mPDG)<0.01)){
2112  // FILE *f=fopen("4display.txt","a");
2113  // fprintf(f,"pt: %f \n Rapidity: %f \t Period Number: %x \t Run Number: %d \t BunchCrossNumb: %d \t OrbitNumber: %d\n",part->Pt(),part->Y(421),aod->GetPeriodNumber(),aod->GetRunNumber(),aod->GetBunchCrossNumber(),aod->GetOrbitNumber());
2114  // fclose(f);
2115  // //printf("PrimVtx NContributors: %d \n Prongs Rel Angle: %f \n \n",ncont,relangle);
2116  // }
2117  // }
2118 
2119  TString fillthis="", fillthispt="", fillthismasspt="", fillthismassy="";
2120  Int_t pdgDgD0toKpi[2]={321,211};
2121  Int_t labD0=-1;
2122  Bool_t isPrimary=kTRUE;
2123  if (fReadMC) labD0 = part->MatchToMC(421,arrMC,2,pdgDgD0toKpi); //return MC particle label if the array corresponds to a D0, -1 if not (cf. AliAODRecoDecay.cxx)
2124 
2125  //Define weights for Bayesian (if appropriate)
2126 
2127  Double_t weigD0=1.;
2128  Double_t weigD0bar=1.;
2130  weigD0=fCuts->GetWeightsNegative()[AliPID::kKaon] * fCuts->GetWeightsPositive()[AliPID::kPion];
2131  weigD0bar=fCuts->GetWeightsPositive()[AliPID::kKaon] * fCuts->GetWeightsNegative()[AliPID::kPion];
2132  if (weigD0 > 1.0 || weigD0 < 0.) {weigD0 = 0.;}
2133  if (weigD0bar > 1.0 || weigD0bar < 0.) {weigD0bar = 0.;} //Prevents filling with weight > 1, or < 0
2134  }
2135 
2136  //count candidates selected by cuts
2137  fNentries->Fill(1);
2138  //count true D0 selected by cuts
2139  if (fReadMC && labD0>=0) fNentries->Fill(2);
2140 
2141  if ((fIsSelectedCandidate==1 || fIsSelectedCandidate==3) && fFillOnlyD0D0bar<2) { //D0
2142 
2143  arrayForSparse[0]=invmassD0; arrayForSparse[2]=impparXY;
2144 
2145  if(fReadMC){
2146  if(labD0>=0) {
2147  if(fArray==1) cout<<"LS signal ERROR"<<endl;
2148 
2149  AliAODMCParticle *partD0 = (AliAODMCParticle*)arrMC->At(labD0);
2150  Int_t pdgD0 = partD0->GetPdgCode();
2151  // cout<<"pdg = "<<pdgD0<<endl;
2152 
2153  // old function if(CheckOrigin(arrMC,partD0)==5) isPrimary=kFALSE;
2154  if(AliVertexingHFUtils::CheckOrigin(arrMC,partD0,fUseQuarkTagInKine)==5) isPrimary=kFALSE;
2155  if(!isPrimary)
2156  trueImpParXY=GetTrueImpactParameter(mcHeader,arrMC,partD0)*10000.;
2157  arrayForSparseTrue[0]=invmassD0; arrayForSparseTrue[2]=trueImpParXY;
2158 
2159  if (pdgD0==421){ //D0
2160  // cout<<"Fill S with D0"<<endl;
2161  fillthis="histSgn_";
2162  fillthis+=ptbin;
2163  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0,weigD0);
2164 
2165  if(fFillPtHist){
2166  fillthismasspt="histSgnPt";
2167  ((TH2F*)(fOutputMassPt->FindObject(fillthismasspt)))->Fill(invmassD0,pt,weigD0);
2168  }
2169  if(fFillImpParHist){
2170  if(isPrimary) fHistMassPtImpParTC[1]->Fill(arrayForSparse,weigD0);
2171  else {
2172  fHistMassPtImpParTC[2]->Fill(arrayForSparse,weigD0);
2173  fHistMassPtImpParTC[3]->Fill(arrayForSparseTrue,weigD0);
2174  }
2175  }
2176 
2177  if(fFillYHist){
2178  fillthismassy="histSgnY_";
2179  fillthismassy+=ptbin;
2180  ((TH2F*)(fOutputMassY->FindObject(fillthismassy)))->Fill(invmassD0,y,weigD0);
2181  }
2182 
2183  if(fSys==0){
2184  if(TMath::Abs(invmassD0 - mPDG) < 0.027 && fFillVarHists){
2185  fillthis="histSgn27_";
2186  fillthis+=ptbin;
2187  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0,weigD0);
2188  }
2189  }
2190  } else{ //it was a D0bar
2191  fillthis="histRfl_";
2192  fillthis+=ptbin;
2193  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0,weigD0);
2194 
2195  if(fFillPtHist){
2196  fillthismasspt="histRflPt";
2197  // cout << " Filling "<<fillthismasspt<<" D0bar"<<endl;
2198  ((TH2F*)(fOutputMassPt->FindObject(fillthismasspt)))->Fill(invmassD0,pt,weigD0);
2199  }
2200 
2201  if(fFillYHist){
2202  fillthismassy="histRflY_";
2203  fillthismassy+=ptbin;
2204  // cout << " Filling "<<fillthismassy<<" D0bar"<<endl;
2205  ((TH2F*)(fOutputMassY->FindObject(fillthismassy)))->Fill(invmassD0,y,weigD0);
2206  }
2207 
2208  }
2209  } else {//background
2210  fillthis="histBkg_";
2211  fillthis+=ptbin;
2212  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0,weigD0);
2213 
2214  if(fFillPtHist){
2215  fillthismasspt="histBkgPt";
2216  // cout << " Filling "<<fillthismasspt<<" D0bar"<<endl;
2217  ((TH2F*)(fOutputMassPt->FindObject(fillthismasspt)))->Fill(invmassD0,pt,weigD0);
2218  }
2219  if(fFillImpParHist) fHistMassPtImpParTC[4]->Fill(arrayForSparse,weigD0);
2220 
2221  if(fFillYHist){
2222  fillthismassy="histBkgY_";
2223  fillthismassy+=ptbin;
2224  // cout << " Filling "<<fillthismassy<<" D0bar"<<endl;
2225  ((TH2F*)(fOutputMassY->FindObject(fillthismassy)))->Fill(invmassD0,y,weigD0);
2226  }
2227 
2228  }
2229 
2230  }else{
2231  fillthis="histMass_";
2232  fillthis+=ptbin;
2233  // cout<<"Filling "<<fillthis<<endl;
2234 
2235  // printf("Fill mass with D0");
2236  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0,weigD0);
2237 
2238 
2239  if(fFillPtHist){
2240  fillthismasspt="histMassPt";
2241  // cout<<"Filling "<<fillthismasspt<<endl;
2242  ((TH2F*)(fOutputMassPt->FindObject(fillthismasspt)))->Fill(invmassD0,pt,weigD0);
2243  }
2244  if(fFillImpParHist) {
2245  // cout << "Filling fHistMassPtImpParTC[0]"<<endl;
2246  fHistMassPtImpParTC[0]->Fill(arrayForSparse,weigD0);
2247  }
2248 
2249  if(fFillYHist){
2250  fillthismassy="histMassY_";
2251  fillthismassy+=ptbin;
2252  // cout<<"Filling "<<fillthismassy<<endl;
2253  ((TH2F*)(fOutputMassY->FindObject(fillthismassy)))->Fill(invmassD0,y,weigD0);
2254  }
2255 
2256  }
2257 
2258  }
2259  if (fIsSelectedCandidate>1 && (fFillOnlyD0D0bar==0 || fFillOnlyD0D0bar==2)) { //D0bar
2260 
2261  arrayForSparse[0]=invmassD0bar; arrayForSparse[2]=impparXY;
2262 
2263  if(fReadMC){
2264  if(labD0>=0) {
2265  if(fArray==1) cout<<"LS signal ERROR"<<endl;
2266  AliAODMCParticle *partD0 = (AliAODMCParticle*)arrMC->At(labD0);
2267  Int_t pdgD0 = partD0->GetPdgCode();
2268  // cout<<" pdg = "<<pdgD0<<endl;
2269 
2270  //old function if(CheckOrigin(arrMC,partD0)==5) isPrimary=kFALSE;
2271  if(AliVertexingHFUtils::CheckOrigin(arrMC,partD0,fUseQuarkTagInKine)==5) isPrimary=kFALSE;
2272  if(!isPrimary)
2273  trueImpParXY=GetTrueImpactParameter(mcHeader,arrMC,partD0)*10000.;
2274  arrayForSparseTrue[0]=invmassD0bar; arrayForSparseTrue[2]=trueImpParXY;
2275 
2276  if (pdgD0==-421){ //D0bar
2277  fillthis="histSgn_";
2278  fillthis+=ptbin;
2279  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar,weigD0bar);
2280  // if (TMath::Abs(invmassD0bar - mPDG) < 0.027){
2281  // fillthis="histSgn27_";
2282  // fillthis+=ptbin;
2283  // ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar);
2284  // }
2285 
2286  if(fFillPtHist){
2287  fillthismasspt="histSgnPt";
2288  // cout<<" Filling "<< fillthismasspt << endl;
2289  ((TH2F*)(fOutputMassPt->FindObject(fillthismasspt)))->Fill(invmassD0bar,pt,weigD0bar);
2290  }
2291  if(fFillImpParHist){
2292  // cout << " Filling impact parameter thnsparse"<<endl;
2293  if(isPrimary) fHistMassPtImpParTC[1]->Fill(arrayForSparse,weigD0bar);
2294  else {
2295  fHistMassPtImpParTC[2]->Fill(arrayForSparse,weigD0bar);
2296  fHistMassPtImpParTC[3]->Fill(arrayForSparseTrue,weigD0bar);
2297  }
2298  }
2299 
2300  if(fFillYHist){
2301  fillthismassy="histSgnY_";
2302  fillthismassy+=ptbin;
2303  // cout<<" Filling "<< fillthismassy << endl;
2304  ((TH2F*)(fOutputMassY->FindObject(fillthismassy)))->Fill(invmassD0bar,y,weigD0bar);
2305  }
2306 
2307  } else{
2308  fillthis="histRfl_";
2309  fillthis+=ptbin;
2310  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar,weigD0bar);
2311  if(fFillPtHist){
2312  fillthismasspt="histRflPt";
2313  // cout << " Filling "<<fillthismasspt<<endl;
2314  ((TH2F*)(fOutputMassPt->FindObject(fillthismasspt)))->Fill(invmassD0bar,pt,weigD0bar);
2315  }
2316  if(fFillYHist){
2317  fillthismassy="histRflY_";
2318  fillthismassy+=ptbin;
2319  // cout << " Filling "<<fillthismassy<<endl;
2320  ((TH2F*)(fOutputMassY->FindObject(fillthismassy)))->Fill(invmassD0bar,y,weigD0bar);
2321  }
2322  }
2323  } else {//background or LS
2324  fillthis="histBkg_";
2325  fillthis+=ptbin;
2326  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar,weigD0bar);
2327 
2328  if(fFillPtHist){
2329  fillthismasspt="histBkgPt";
2330  // cout<<" Filling "<< fillthismasspt << endl;
2331  ((TH2F*)(fOutputMassPt->FindObject(fillthismasspt)))->Fill(invmassD0bar,pt,weigD0bar);
2332  }
2333  if(fFillImpParHist) fHistMassPtImpParTC[4]->Fill(arrayForSparse,weigD0bar);
2334  if(fFillYHist){
2335  fillthismassy="histBkgY_";
2336  fillthismassy+=ptbin;
2337  // cout<<" Filling "<< fillthismassy << endl;
2338  ((TH2F*)(fOutputMassY->FindObject(fillthismassy)))->Fill(invmassD0bar,y,weigD0bar);
2339  }
2340  }
2341  }else{
2342  fillthis="histMass_";
2343  fillthis+=ptbin;
2344  // printf("Fill mass with D0bar");
2345 
2346  ((TH1F*)listout->FindObject(fillthis))->Fill(invmassD0bar,weigD0bar);
2347 
2348 
2349  if(fFillPtHist){
2350  fillthismasspt="histMassPt";
2351  // cout<<" Filling "<< fillthismasspt << endl;
2352  ((TH2F*)(fOutputMassPt->FindObject(fillthismasspt)))->Fill(invmassD0bar,pt,weigD0bar);
2353  }
2354  if(fFillImpParHist) fHistMassPtImpParTC[0]->Fill(arrayForSparse,weigD0bar);
2355  if(fFillYHist){
2356  fillthismassy="histMassY_";
2357  fillthismassy+=ptbin;
2358  // cout<<" Filling "<< fillthismassy << endl;
2359  ((TH2F*)(fOutputMassY->FindObject(fillthismassy)))->Fill(invmassD0bar,y,weigD0bar);
2360  }
2361  }
2362  }
2363 
2364  return;
2365 }
2366 
2367 //__________________________________________________________________________
2370 
2371  Int_t skipped[2];
2372  Int_t nTrksToSkip=2;
2373  AliAODTrack *dgTrack = (AliAODTrack*)fDaughterTracks.UncheckedAt(0);
2374  if(!dgTrack){
2375  AliDebug(2,"no daughter found!");
2376  return 0x0;
2377  }
2378  skipped[0]=dgTrack->GetID();
2379  dgTrack = (AliAODTrack*)fDaughterTracks.UncheckedAt(1);
2380  if(!dgTrack){
2381  AliDebug(2,"no daughter found!");
2382  return 0x0;
2383  }
2384  skipped[1]=dgTrack->GetID();
2385 
2386  AliESDVertex *vertexESD=0x0;
2387  AliAODVertex *vertexAOD=0x0;
2388  AliVertexerTracks *vertexer = new AliVertexerTracks(aodev->GetMagneticField());
2389 
2390  //
2391  vertexer->SetSkipTracks(nTrksToSkip,skipped);
2392  vertexer->SetMinClusters(4);
2393  vertexESD = (AliESDVertex*)vertexer->FindPrimaryVertex(aodev);
2394  if(!vertexESD) return vertexAOD;
2395  if(vertexESD->GetNContributors()<=0) {
2396  AliDebug(2,"vertexing failed");
2397  delete vertexESD; vertexESD=NULL;
2398  return vertexAOD;
2399  }
2400 
2401  delete vertexer; vertexer=NULL;
2402 
2403 
2404  // convert to AliAODVertex
2405  Double_t pos[3],cov[6],chi2perNDF;
2406  vertexESD->GetXYZ(pos); // position
2407  vertexESD->GetCovMatrix(cov); //covariance matrix
2408  chi2perNDF = vertexESD->GetChi2toNDF();
2409  delete vertexESD; vertexESD=NULL;
2410 
2411  vertexAOD = new AliAODVertex(pos,cov,chi2perNDF);
2412  return vertexAOD;
2413 
2414 }
2415 
2416 
2417 //________________________________________________________________________
2419 {
2421  //
2422  if(fDebug > 1) printf("AnalysisTaskSED0Mass: Terminate() \n");
2423 
2424 
2425  fOutputMass = dynamic_cast<TList*> (GetOutputData(1));
2426  if (!fOutputMass) {
2427  printf("ERROR: fOutputMass not available\n");
2428  return;
2429  }
2430  fOutputMassPt = dynamic_cast<TList*> (GetOutputData(6));
2431  if ((fFillPtHist || fFillImpParHist) && !fOutputMassPt) {
2432  printf("ERROR: fOutputMass not available\n");
2433  return;
2434  }
2435 
2436  if(fFillVarHists){
2437  fDistr = dynamic_cast<TList*> (GetOutputData(2));
2438  if (!fDistr) {
2439  printf("ERROR: fDistr not available\n");
2440  return;
2441  }
2442  }
2443 
2444  fNentries = dynamic_cast<TH1F*>(GetOutputData(3));
2445 
2446  if(!fNentries){
2447  printf("ERROR: fNEntries not available\n");
2448  return;
2449  }
2450  fCuts = dynamic_cast<AliRDHFCutsD0toKpi*>(GetOutputData(4));
2451  if(!fCuts){
2452  printf("ERROR: fCuts not available\n");
2453  return;
2454  }
2455  fCounter = dynamic_cast<AliNormalizationCounter*>(GetOutputData(5));
2456  if (!fCounter) {
2457  printf("ERROR: fCounter not available\n");
2458  return;
2459  }
2460  if (fDrawDetSignal) {
2461  fDetSignal = dynamic_cast<TList*>(GetOutputData(8));
2462  if (!fDetSignal) {
2463  printf("ERROR: fDetSignal not available\n");
2464  return;
2465  }
2466  }
2467  if(fFillYHist){
2468  fOutputMassY = dynamic_cast<TList*> (GetOutputData(9));
2469  if (fFillYHist && !fOutputMassY) {
2470  printf("ERROR: fOutputMassY not available\n");
2471  return;
2472  }
2473  }
2474 
2476  for(Int_t ipt=0;ipt<nptbins;ipt++){
2477 
2478  if(fArray==1 && fFillVarHists){
2479  fLsNormalization = 2.*TMath::Sqrt(((TH1F*)fOutputMass->FindObject("hpospair"))->Integral(nptbins+ipt+1,nptbins+ipt+2)*((TH1F*)fOutputMass->FindObject("hnegpair"))->Integral(nptbins+ipt+1,nptbins+ipt+2)); //after cuts
2480 
2481 
2482  if(fLsNormalization>1e-6) {
2483 
2484  TString massName="histMass_";
2485  massName+=ipt;
2486  ((TH1F*)fOutputMass->FindObject(massName))->Scale((1/fLsNormalization)*((TH1F*)fOutputMass->FindObject(massName))->GetEntries());
2487 
2488  }
2489 
2490 
2491  fLsNormalization = 2.*TMath::Sqrt(((TH1F*)fOutputMass->FindObject("hpospair"))->Integral(ipt+1,ipt+2)*((TH1F*)fOutputMass->FindObject("hnegpair"))->Integral(ipt+1,ipt+2));
2492  //fLsNormalization = 2.*TMath::Sqrt(fTotPosPairs[4]*fTotNegPairs[4]);
2493 
2494  if(fLsNormalization>1e-6) {
2495 
2496  TString nameDistr="hdcaB_";
2497  nameDistr+=ipt;
2498  ((TH1F*)fDistr->FindObject(nameDistr))->Scale((1/fLsNormalization)*((TH1F*)fDistr->FindObject(nameDistr))->GetEntries());
2499  nameDistr="hd0B_";
2500  nameDistr+=ipt;
2501  ((TH1F*)fDistr->FindObject(nameDistr))->Scale((1/fLsNormalization)*((TH1F*)fDistr->FindObject(nameDistr))->GetEntries());
2502  nameDistr="hd0d0B_";
2503  nameDistr+=ipt;
2504  ((TH1F*)fDistr->FindObject(nameDistr))->Scale((1/fLsNormalization)*((TH1F*)fDistr->FindObject(nameDistr))->GetEntries());
2505  nameDistr="hcosthetapointB_";
2506  nameDistr+=ipt;
2507  ((TH1F*)fDistr->FindObject(nameDistr))->Scale((1/fLsNormalization)*((TH1F*)fDistr->FindObject(nameDistr))->GetEntries());
2508  if(fSys==0){
2509  nameDistr="hptB_";
2510  nameDistr+=ipt;
2511  ((TH1F*)fDistr->FindObject(nameDistr))->Scale((1/fLsNormalization)*((TH1F*)fDistr->FindObject(nameDistr))->GetEntries());
2512  nameDistr="hcosthetastarB_";
2513  nameDistr+=ipt;
2514  ((TH1F*)fDistr->FindObject(nameDistr))->Scale((1/fLsNormalization)*((TH1F*)fDistr->FindObject(nameDistr))->GetEntries());
2515  nameDistr="hcosthpointd0d0B_";
2516  nameDistr+=ipt;
2517  ((TH2F*)fDistr->FindObject(nameDistr))->Scale((1/fLsNormalization)*((TH2F*)fDistr->FindObject(nameDistr))->GetEntries());
2518  }
2519  }
2520  }
2521  }
2522  TString cvname,cstname;
2523 
2524  if (fArray==0){
2525  cvname="D0invmass";
2526  cstname="cstat0";
2527  } else {
2528  cvname="LSinvmass";
2529  cstname="cstat1";
2530  }
2531 
2532  TCanvas *cMass=new TCanvas(cvname,cvname);
2533  cMass->cd();
2534  ((TH1F*)fOutputMass->FindObject("histMass_3"))->Draw();
2535 
2536  TCanvas* cStat=new TCanvas(cstname,Form("Stat%s",fArray ? "LS" : "D0"));
2537  cStat->cd();
2538  cStat->SetGridy();
2539  fNentries->Draw("htext0");
2540 
2541  // TCanvas *ccheck=new TCanvas(Form("cc%d",fArray),Form("cc%d",fArray));
2542  // ccheck->cd();
2543 
2544  return;
2545 }
2546 
2547 
2548 //________________________________________________________________________
2551 
2552  Int_t nmassbins=200;
2553  Double_t fLowmasslimit=1.5648, fUpmasslimit=2.1648;
2554  Int_t fNImpParBins=400;
2555  Double_t fLowerImpPar=-2000., fHigherImpPar=2000.;
2556  Int_t nbins[3]={nmassbins,200,fNImpParBins};
2557  Double_t xmin[3]={fLowmasslimit,0.,fLowerImpPar};
2558  Double_t xmax[3]={fUpmasslimit,20.,fHigherImpPar};
2559 
2560 
2561  fHistMassPtImpParTC[0]=new THnSparseF("hMassPtImpParAll",
2562  "Mass vs. pt vs.imppar - All",
2563  3,nbins,xmin,xmax);
2564  fHistMassPtImpParTC[1]=new THnSparseF("hMassPtImpParPrompt",
2565  "Mass vs. pt vs.imppar - promptD",
2566  3,nbins,xmin,xmax);
2567  fHistMassPtImpParTC[2]=new THnSparseF("hMassPtImpParBfeed",
2568  "Mass vs. pt vs.imppar - DfromB",
2569  3,nbins,xmin,xmax);
2570  fHistMassPtImpParTC[3]=new THnSparseF("hMassPtImpParTrueBfeed",
2571  "Mass vs. pt vs.true imppar -DfromB",
2572  3,nbins,xmin,xmax);
2573  fHistMassPtImpParTC[4]=new THnSparseF("hMassPtImpParBkg",
2574  "Mass vs. pt vs.imppar - backgr.",
2575  3,nbins,xmin,xmax);
2576 
2577  for(Int_t i=0; i<5;i++){
2579  }
2580 }
2581 
2582 //_________________________________________________________________________________________________
2583 Float_t AliAnalysisTaskSED0Mass::GetTrueImpactParameter(AliAODMCHeader *mcHeader, TClonesArray* arrayMC, AliAODMCParticle *partD0) const {
2585 
2586  printf(" AliAnalysisTaskSED0MassV1::GetTrueImpactParameter() \n");
2587 
2588  Double_t vtxTrue[3];
2589  mcHeader->GetVertex(vtxTrue);
2590  Double_t origD[3];
2591  partD0->XvYvZv(origD);
2592  Short_t charge=partD0->Charge();
2593  Double_t pXdauTrue[2],pYdauTrue[2],pZdauTrue[2];
2594  for(Int_t iDau=0; iDau<2; iDau++){
2595  pXdauTrue[iDau]=0.;
2596  pYdauTrue[iDau]=0.;
2597  pZdauTrue[iDau]=0.;
2598  }
2599 
2600  // Int_t nDau=partD0->GetNDaughters();
2601  Int_t labelFirstDau = partD0->GetDaughter(0);
2602 
2603  for(Int_t iDau=0; iDau<2; iDau++){
2604  Int_t ind = labelFirstDau+iDau;
2605  AliAODMCParticle* part = dynamic_cast<AliAODMCParticle*>(arrayMC->At(ind));
2606  if(!part) continue;
2607  Int_t pdgCode = TMath::Abs( part->GetPdgCode() );
2608  if(!part){
2609  AliError("Daughter particle not found in MC array");
2610  return 99999.;
2611  }
2612  if(pdgCode==211 || pdgCode==321){
2613  pXdauTrue[iDau]=part->Px();
2614  pYdauTrue[iDau]=part->Py();
2615  pZdauTrue[iDau]=part->Pz();
2616  }
2617  }
2618 
2619  Double_t d0dummy[2]={0.,0.};
2620  AliAODRecoDecayHF aodDzeroMC(vtxTrue,origD,2,charge,pXdauTrue,pYdauTrue,pZdauTrue,d0dummy);
2621  return aodDzeroMC.ImpParXY();
2622 
2623 }
2624 
2625 //_________________________________________________________________________________________________
2626 Int_t AliAnalysisTaskSED0Mass::CheckOrigin(TClonesArray* arrayMC, AliAODMCParticle *mcPartCandidate) const {
2627  // obsolete method
2629  //
2630  printf(" AliAnalysisTaskSED0Mass V1::CheckOrigin() \n");
2631 
2632  Int_t pdgGranma = 0;
2633  Int_t mother = 0;
2634  mother = mcPartCandidate->GetMother();
2635  Int_t istep = 0;
2636  Int_t abspdgGranma =0;
2637  Bool_t isFromB=kFALSE;
2638  Bool_t isQuarkFound=kFALSE;
2639  while (mother >0 ){
2640  istep++;
2641  AliAODMCParticle* mcGranma = dynamic_cast<AliAODMCParticle*>(arrayMC->At(mother));
2642  if (mcGranma){
2643  pdgGranma = mcGranma->GetPdgCode();
2644  abspdgGranma = TMath::Abs(pdgGranma);
2645  if ((abspdgGranma > 500 && abspdgGranma < 600) || (abspdgGranma > 5000 && abspdgGranma < 6000)){
2646  isFromB=kTRUE;
2647  }
2648  if(abspdgGranma==4 || abspdgGranma==5) isQuarkFound=kTRUE;
2649  mother = mcGranma->GetMother();
2650  }else{
2651  AliError("Failed casting the mother particle!");
2652  break;
2653  }
2654  }
2655 
2656  if(isFromB) return 5;
2657  else return 4;
2658 }
2659 //_______________________________________
2662 
2663  const Int_t nVarPrompt = 2;
2664  const Int_t nVarFD = 3;
2665 
2666  Int_t nbinsPrompt[nVarPrompt]={200,100};
2667  Int_t nbinsFD[nVarFD]={200,100,200};
2668 
2669  Double_t xminPrompt[nVarPrompt] = {0.,-1.};
2670  Double_t xmaxPrompt[nVarPrompt] = {40.,1.};
2671 
2672  Double_t xminFD[nVarFD] = {0.,-1.,0.};
2673  Double_t xmaxFD[nVarFD] = {40.,1.,40.};
2674 
2675  //pt, y
2676  fMCAccPrompt = new THnSparseF("hMCAccPrompt","kStepMCAcceptance pt vs. y - promptD",nVarPrompt,nbinsPrompt,xminPrompt,xmaxPrompt);
2677  fMCAccPrompt->GetAxis(0)->SetTitle("p_{T} (GeV/c)");
2678  fMCAccPrompt->GetAxis(1)->SetTitle("y");
2679 
2680  //pt,y,ptB
2681  fMCAccBFeed = new THnSparseF("hMCAccBFeed","kStepMCAcceptance pt vs. y vs. ptB - DfromB",nVarFD,nbinsFD,xminFD,xmaxFD);
2682  fMCAccBFeed->GetAxis(0)->SetTitle("p_{T} (GeV/c)");
2683  fMCAccBFeed->GetAxis(1)->SetTitle("y");
2684  fMCAccBFeed->GetAxis(2)->SetTitle("p_{T}^{B} (GeV/c)");
2685 
2686  fOutputMass->Add(fMCAccPrompt);
2687  fOutputMass->Add(fMCAccBFeed);
2688 }
2689 //___________________________________________________________________________________________________
2690 void AliAnalysisTaskSED0Mass::FillMCAcceptanceHistos(TClonesArray *arrayMC, AliAODMCHeader *mcHeader){
2692  const Int_t nProng = 2;
2693  Double_t zMCVertex = mcHeader->GetVtxZ(); //vertex MC
2694 
2695  for(Int_t iPart=0; iPart<arrayMC->GetEntriesFast(); iPart++){
2696  AliAODMCParticle* mcPart = dynamic_cast<AliAODMCParticle*>(arrayMC->At(iPart));
2697  if (TMath::Abs(mcPart->GetPdgCode()) == 421){
2698 
2699  Int_t orig=AliVertexingHFUtils::CheckOrigin(arrayMC,mcPart,fUseQuarkTagInKine);//Prompt = 4, FeedDown = 5
2700 
2701  Int_t deca = 0;
2702  Bool_t isGoodDecay=kFALSE;
2703  Int_t labDau[4]={-1,-1,-1,-1};
2704  Bool_t isInAcc = kFALSE;
2705  Bool_t isFidAcc = kFALSE;
2706 
2707  deca=AliVertexingHFUtils::CheckD0Decay(arrayMC,mcPart,labDau);
2708  if(deca > 0) isGoodDecay=kTRUE;
2709 
2710  if(labDau[0]==-1){
2711  continue; //protection against unfilled array of labels
2712  }
2713 
2714  isFidAcc=fCuts->IsInFiducialAcceptance(mcPart->Pt(),mcPart->Y());
2715  isInAcc=CheckAcc(arrayMC,nProng,labDau);
2716 
2717  if(isGoodDecay && TMath::Abs(zMCVertex) < fCuts->GetMaxVtxZ() && isFidAcc && isInAcc) {
2718  //for prompt
2719  if(orig == 4){
2720  //fill histo for prompt
2721  Double_t arrayMCprompt[2] = {mcPart->Pt(),mcPart->Y()};
2722  fMCAccPrompt->Fill(arrayMCprompt);
2723  }
2724  //for FD
2725  else if(orig == 5){
2726  Double_t ptB = AliVertexingHFUtils::GetBeautyMotherPt(arrayMC,mcPart);
2727  //fill histo for FD
2728  Double_t arrayMCFD[3] = {mcPart->Pt(),mcPart->Y(),ptB};
2729  fMCAccBFeed->Fill(arrayMCFD);
2730  }
2731  else
2732  continue;
2733  }
2734  }
2735  }
2736 }
2737 //______________________________________________________________________
2738 Bool_t AliAnalysisTaskSED0Mass::CheckAcc(TClonesArray* arrayMC,Int_t nProng, Int_t *labDau){
2740  for (Int_t iProng = 0; iProng<nProng; iProng++){
2741  AliAODMCParticle* mcPartDaughter=dynamic_cast<AliAODMCParticle*>(arrayMC->At(labDau[iProng]));
2742  if(!mcPartDaughter) return kFALSE;
2743  Double_t eta = mcPartDaughter->Eta();
2744  Double_t pt = mcPartDaughter->Pt();
2745  if (TMath::Abs(eta) > 0.9 || pt < 0.1) return kFALSE;
2746  }
2747  return kTRUE;
2748 }
2749 //________________________________________
2751  Double_t normIP[2];
2752  Double_t pointD0[9];
2753  Double_t pointD0bar[9];
2754  Double_t pointD0MC[9];
2755  Float_t ptB;
2756  Int_t signalType=0;//background
2757 
2758  //fill variables
2759  Double_t diffIP[2], errdiffIP[2];
2760  part->Getd0MeasMinusExpProng(0,aod->GetMagneticField(),diffIP[0],errdiffIP[0]);
2761  part->Getd0MeasMinusExpProng(1,aod->GetMagneticField(),diffIP[1],errdiffIP[1]);
2762  normIP[0]=diffIP[0]/errdiffIP[0];
2763  normIP[1]=diffIP[1]/errdiffIP[1];
2764 
2765  AliAODVertex* secvtx=part->GetSecondaryVtx();
2766  AliAODVertex *primvtx=part->GetPrimaryVtx();
2767  Double_t err2decaylength=secvtx->Error2DistanceXYToVertex(primvtx);
2768  Double_t lxy=part->AliAODRecoDecay::DecayLengthXY(primvtx);
2769  Bool_t ispid=fCuts->GetIsUsePID();
2770 
2771  if(fReadMC){
2772  // MC THnSparse
2773  Int_t pdgDgD0toKpi[2]={321,211};
2774  Int_t lab=-9999;
2775  lab=part->MatchToMC(421,arrMC,2,pdgDgD0toKpi); //return MC particle label if the array corresponds to a D0, -1 if not (cf. AliAODRecoDecay.cxx)
2776  if(lab>=0){
2777  signalType=1;//signal
2778  AliAODMCParticle *partD0 = (AliAODMCParticle*)arrMC->At(lab);
2780  if(orig==4)signalType=1;//prompt
2781  else if(orig==5)signalType=2;//feed down
2782  //pt, ptB, normImpParTrk1, normImpParTrk2, decLXY, normDecLXY, iscut, ispid
2783  if(ispid)fCuts->SetUsePID(kFALSE);// if PID on, switch it off
2784  Int_t isCuts=fCuts->IsSelected(part,AliRDHFCuts::kAll,aod);
2785  if(ispid)fCuts->SetUsePID(kTRUE);//if PID was on, switch it on
2786  if(!isCuts) return;
2787  Int_t isPid= fCuts->IsSelectedPID(part);
2788  ptB=AliVertexingHFUtils::GetBeautyMotherPt(arrMC,partD0);
2789  Double_t arrayMC[8]={part->Pt(), ptB,normIP[0], normIP[1], lxy, lxy/TMath::Sqrt(err2decaylength),(Double_t)isCuts, (Double_t)isPid};
2790 
2791  //fill pointD0MC
2792  for(Int_t i=0; i<8; i++){
2793  pointD0MC[i]=arrayMC[i];
2794  }
2795  if(signalType==1)fhStudyImpParSingleTrackSign->Fill(pointD0MC);
2796  if(signalType==2){
2797  fhStudyImpParSingleTrackFd->Fill(pointD0MC);
2798  }
2799  }
2800  }else{
2801  // data
2802  if(ispid)fCuts->SetUsePID(kFALSE);// if PID on, switch it off
2803  Int_t IsSelectedPIDoff=fCuts->IsSelected(part,AliRDHFCuts::kAll,aod);
2804  if(ispid)fCuts->SetUsePID(kTRUE);//if PID was on, switch it on
2805  if(!IsSelectedPIDoff) return;
2806  Int_t pid=fCuts->IsSelectedPID(part);
2807  if((IsSelectedPIDoff==1 || IsSelectedPIDoff==3) && fFillOnlyD0D0bar<2){
2808  Double_t array[9]={part->Pt(),normIP[0],normIP[1],lxy,lxy/TMath::Sqrt(err2decaylength),part->InvMassD0(),(Double_t)IsSelectedPIDoff,(Double_t)pid,0.};
2809  for(Int_t i=0;i<9;i++){
2810  pointD0[i]=array[i];
2811  }
2812  fhStudyImpParSingleTrackCand->Fill(pointD0);
2813  }
2814  if (IsSelectedPIDoff>1 && (fFillOnlyD0D0bar==0 || fFillOnlyD0D0bar==2)){
2815  Double_t arrayD0bar[9]={part->Pt(),normIP[0],normIP[1],lxy,lxy/TMath::Sqrt(err2decaylength),part->InvMassD0bar(),(Double_t)IsSelectedPIDoff,(Double_t)pid,1.};
2816  for(Int_t i=0;i<9;i++){
2817  pointD0bar[i]=arrayD0bar[i];
2818  }
2819  fhStudyImpParSingleTrackCand->Fill(pointD0bar);
2820  }
2821  }
2822 }
Int_t charge
THnSparseF * fMCAccBFeed
!histo for StepMCAcc for D0 FD (pt,y,ptB)
Double_t NormalizedDecayLengthXY() const
TList * fOutputMassPt
! list send on output slot 6
Bool_t fDrawDetSignal
flag to decide whether to fill "PID = x" bins in fNentrie
double Double_t
Definition: External.C:58
void Draw(const char *filename, const char *title="", const char *others="ALL", const char *options="DEFAULT", const char *outFlg="ALL", UShort_t rebin=5, Float_t eff=0, const char *base="")
Definition: DrawdNdeta.C:3603
Bool_t fFillYHist
flag to fill Pt and Impact Parameter Histograms
void StoreCandidates(AliVEvent *, Int_t nCand=0, Bool_t flagFilter=kTRUE)
Definition: External.C:236
void Getd0MeasMinusExpProng(Int_t ip, Double_t magf, Double_t &d0diff, Double_t &errd0diff) const
void FillVarHists(AliAODEvent *aodev, AliAODRecoDecayHF2Prong *part, TClonesArray *arrMC, AliRDHFCutsD0toKpi *cuts, TList *listout)
Bool_t GetCombPID() const
Double_t * fCandidateVariables
! variables to be written to the tree
Bool_t HasSelectionBit(Int_t i) const
Double_t NormalizedDecayLength2() const
Double_t ImpParXY() const
Int_t fAODProtection
flag to check or not the selection bit
static Int_t CheckMatchingAODdeltaAODevents()
TList * fDistr
! list send on output slot 2
void NormIPvar(AliAODEvent *aod, AliAODRecoDecayHF2Prong *part, TClonesArray *arrMC)
Int_t GetWhyRejection() const
Definition: AliRDHFCuts.h:297
TList * fDetSignal
!Detector signal histograms (on output slot 8)
Double_t CosPointingAngleXY() const
Bool_t fUseSelectionBit
flag to fill Pt and Impact Parameter Histograms
Bool_t FillRecoCand(AliVEvent *event, AliAODRecoDecayHF3Prong *rd3)
void FillMCAcceptanceHistos(TClonesArray *arrayMC, AliAODMCHeader *mcHeader)
Int_t fSys
flag to enable filling variable histos
ULong_t GetSelectionMap() const
TList * fOutputMass
! list send on output slot 1
virtual Bool_t IsInFiducialAcceptance(Double_t pt, Double_t y) const
Bool_t CheckAcc(TClonesArray *arrayMC, Int_t nProng, Int_t *labDau)
Bool_t fReadMC
can be D0 or Like Sign candidates
static Int_t CheckOrigin(TClonesArray *arrayMC, AliAODMCParticle *mcPart, Bool_t searchUpToQuark=kTRUE)
Double_t GetMaxVtxZ() const
Definition: AliRDHFCuts.h:241
Int_t GetBayesianStrategy() const
Bool_t HasBadDaughters() const
THnSparseF * fHistMassPtImpParTC[5]
! histograms for impact paramter studies
TH1F * fNentries
! histogram with number of events on output slot 3
virtual void UserCreateOutputObjects()
Implementation of interface methods.
Bool_t fFillPtHist
flag to reject events with SDD clusters
Bool_t fFillImpParHist
flag to fill Y Histograms
int Int_t
Definition: External.C:63
Definition: External.C:204
unsigned int UInt_t
Definition: External.C:33
Int_t GetIsFilled() const
Double_t CosThetaStarD0bar() const
angle of K
float Float_t
Definition: External.C:68
Bool_t fUseQuarkTagInKine
flag to decide whether to draw the TPC dE/dx and TOF signal before/after PID
Bool_t fIsRejectSDDClusters
fSys=0 -> p-p; fSys=1 ->PbPb (in this case fFillVarHists=kFALSE by default: set it to kTRUE after if ...
static Int_t CheckD0Decay(AliStack *stack, Int_t label, Int_t *arrayDauLab)
Int_t CheckOrigin(TClonesArray *arrayMC, AliAODMCParticle *mcPartCandidate) const
virtual void Terminate(Option_t *option)
Float_t GetTrueImpactParameter(AliAODMCHeader *mcHeader, TClonesArray *arrayMC, AliAODMCParticle *partD0) const
void DrawDetSignal(AliAODRecoDecayHF2Prong *part, TList *ListDetSignal)
Double_t fLsNormalization
number of pt bins
Bool_t fUsePid4Distr
flag to decide if apply cut also on distributions: 0 no cuts, 1 looser cuts, 2 tighter/ cuts ...
Double_t DecayLength2() const
kinematics & topology
THnSparseF * fhStudyImpParSingleTrackFd
! sparse with imp par residual cuts for MC
Bool_t Data(TH1F *h, Double_t *rangefit, Bool_t writefit, Double_t &sgn, Double_t &errsgn, Double_t &bkg, Double_t &errbkg, Double_t &sgnf, Double_t &errsgnf, Double_t &sigmafit, Int_t &status)
Bool_t fCutOnDistr
flag for MC array: kTRUE = read it, kFALSE = do not read it
AliAODVertex * GetPrimaryVtxSkipped(AliAODEvent *aodev)
short Short_t
Definition: External.C:23
Int_t IsSelectedSpecialCuts(AliAODRecoDecayHF *d) const
THnSparseF * fhStudyImpParSingleTrackSign
! sparse with imp par residual cuts for MC
void SetOwnPrimaryVtx(const AliAODVertex *vtx)
Double_t DecayLengthXY() const
static Double_t GetBeautyMotherPt(TClonesArray *arrayMC, AliAODMCParticle *mcPart)
void FillMassHists(AliAODRecoDecayHF2Prong *part, TClonesArray *arrMC, AliAODMCHeader *mcHeader, AliRDHFCutsD0toKpi *cuts, TList *listout)
Int_t fIsSelectedCandidate
keeps the daughter tracks
Bool_t GetIsPrimaryWithoutDaughters() const
Definition: AliRDHFCuts.h:262
Double_t * GetWeightsPositive() const
Bool_t AreDaughtersSelected(AliAODRecoDecayHF *rd, const AliAODEvent *aod=0x0) const
Bool_t IsEventSelected(AliVEvent *event)
void StoreEvent(AliVEvent *, AliRDHFCuts *, Bool_t mc=kFALSE, Int_t multiplicity=-9999, Double_t spherocity=-99.)
void SetUsePID(Bool_t flag=kTRUE)
Definition: AliRDHFCuts.h:204
THnSparseF * fhStudyImpParSingleTrackCand
! sparse with imp par residual cuts for Data
const Int_t nVar
virtual Int_t IsSelectedPID(AliAODRecoDecayHF *rd)
TList * fOutputMassY
! list send on output slot 9
TTree * fVariablesTree
flag to decide whether to write the candidate variables on a tree variables
TH2 * Scale(TH2 *h, TH1 *g)
virtual Int_t IsSelected(TObject *obj, Int_t selectionLevel)
Bool_t GetIsUsePID() const
Definition: AliRDHFCuts.h:254
const char Option_t
Definition: External.C:48
Int_t GetNPtBins() const
Definition: AliRDHFCuts.h:233
AliAODVertex * GetPrimaryVtx() const
virtual void UserExec(Option_t *option)
TObjArray fDaughterTracks
flag to fill mass histogram with D0/D0bar only (0 = fill with both, 1 = fill with D0 only...
const Int_t nbins
bool Bool_t
Definition: External.C:53
Double_t CosPointingAngle() const
Double_t DecayLengthXYError() const
Int_t fFillOnlyD0D0bar
normalization
THnSparseF * fMCAccPrompt
!histo for StepMCAcc for D0 prompt (pt,y,ptB)
Int_t PtBin(Double_t pt) const
Int_t nptbins
Bool_t fFillVarHists
selection outcome
AliNormalizationCounter * fCounter
! AliNormalizationCounter on output slot 5
Double_t * GetWeightsNegative() const
ClassImp(AliAnalysisTaskCascadeTester) AliAnalysisTaskCascadeTester