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