]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGHF/correlationHF/AliAnalysisTaskSED0Correlations.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Usage of D0 efficiency introduced
[u/mrichter/AliRoot.git] / PWGHF / correlationHF / AliAnalysisTaskSED0Correlations.cxx
CommitLineData
a2ad7da1 1/**************************************************************************
2 * Copyright(c) 1998-2012, 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
18/////////////////////////////////////////////////////////////
19//
20// AliAnalysisTaskSE for D0 candidates (2Prongs)
21// and hadrons correlations
22//
23// Authors: Chiara Bianchin, chiara.bianchin@pd.infn.it (invariant mass)
24// Fabio Colamaria, fabio.colamaria@ba.infn.it (correlations)
25/////////////////////////////////////////////////////////////
26
27#include <Riostream.h>
28#include <TClonesArray.h>
29#include <TCanvas.h>
30#include <TNtuple.h>
31#include <TTree.h>
32#include <TList.h>
33#include <TH1F.h>
34#include <TH2F.h>
35#include <THnSparse.h>
36#include <TDatabasePDG.h>
37
38#include <AliAnalysisDataSlot.h>
39#include <AliAnalysisDataContainer.h>
40#include "AliAnalysisManager.h"
41#include "AliESDtrack.h"
42#include "AliVertexerTracks.h"
43#include "AliAODHandler.h"
44#include "AliAODEvent.h"
45#include "AliAODVertex.h"
46#include "AliAODTrack.h"
47#include "AliAODMCHeader.h"
48#include "AliAODMCParticle.h"
49#include "AliAODRecoDecayHF2Prong.h"
50#include "AliAODRecoCascadeHF.h"
51#include "AliAnalysisVertexingHF.h"
52#include "AliAnalysisTaskSE.h"
53#include "AliAnalysisTaskSED0Correlations.h"
54#include "AliNormalizationCounter.h"
55
56using std::cout;
57using std::endl;
58
59ClassImp(AliAnalysisTaskSED0Correlations)
60
61
62//________________________________________________________________________
63AliAnalysisTaskSED0Correlations::AliAnalysisTaskSED0Correlations():
64AliAnalysisTaskSE(),
65 fNPtBinsCorr(0),
66 fBinLimsCorr(),
67 fPtThreshLow(),
68 fPtThreshUp(),
7f221b36 69 fD0Eff(),
a2ad7da1 70 fEvents(0),
71 fAlreadyFilled(kFALSE),
72 fOutputMass(0),
73 fOutputCorr(0),
74 fOutputStudy(0),
75 fNentries(0),
76 fCutsD0(0),
77 fCutsTracks(0),
bce70c96 78 fCorrelatorTr(0),
79 fCorrelatorKc(0),
80 fCorrelatorK0(0),
a2ad7da1 81 fReadMC(0),
4d978a4a 82 fRecoTr(kTRUE),
83 fRecoD0(kTRUE),
84 fSelEvType(kFALSE),
bce70c96 85 fMixing(kFALSE),
a2ad7da1 86 fCounter(0),
87 fNPtBins(1),
88 fFillOnlyD0D0bar(0),
89 fIsSelectedCandidate(0),
90 fSys(0),
8c2d7467 91 fEtaForCorrel(0),
a2ad7da1 92 fIsRejectSDDClusters(0),
93 fFillGlobal(kTRUE)
94{
95 // Default constructor
96
97}
98
99//________________________________________________________________________
100AliAnalysisTaskSED0Correlations::AliAnalysisTaskSED0Correlations(const char *name,AliRDHFCutsD0toKpi* cutsD0, AliHFAssociatedTrackCuts* cutsTrk):
101 AliAnalysisTaskSE(name),
102 fNPtBinsCorr(0),
103 fBinLimsCorr(),
104 fPtThreshLow(),
105 fPtThreshUp(),
7f221b36 106 fD0Eff(),
a2ad7da1 107 fEvents(0),
108 fAlreadyFilled(kFALSE),
109 fOutputMass(0),
110 fOutputCorr(0),
111 fOutputStudy(0),
112 fNentries(0),
113 fCutsD0(0),
114 fCutsTracks(cutsTrk),
bce70c96 115 fCorrelatorTr(0),
116 fCorrelatorKc(0),
117 fCorrelatorK0(0),
a2ad7da1 118 fReadMC(0),
4d978a4a 119 fRecoTr(kTRUE),
120 fRecoD0(kTRUE),
121 fSelEvType(kFALSE),
bce70c96 122 fMixing(kFALSE),
a2ad7da1 123 fCounter(0),
124 fNPtBins(1),
125 fFillOnlyD0D0bar(0),
126 fIsSelectedCandidate(0),
127 fSys(0),
8c2d7467 128 fEtaForCorrel(0),
a2ad7da1 129 fIsRejectSDDClusters(0),
130 fFillGlobal(kTRUE)
131{
132 // Default constructor
133
134 fNPtBins=cutsD0->GetNPtBins();
135
136 fCutsD0=cutsD0;
137
138 // Output slot #1 writes into a TList container (mass with cuts)
139 DefineOutput(1,TList::Class()); //My private output
140 // Output slot #2 writes into a TH1F container (number of events)
141 DefineOutput(2,TH1F::Class()); //My private output
142 // Output slot #3 writes into a AliRDHFD0toKpi container (cuts)
143 DefineOutput(3,AliRDHFCutsD0toKpi::Class()); //My private output
144 // Output slot #4 writes Normalization Counter
145 DefineOutput(4,AliNormalizationCounter::Class());
146 // Output slot #5 writes into a TList container (correl output)
147 DefineOutput(5,TList::Class()); //My private output
148 // Output slot #6 writes into a TList container (correl advanced)
149 DefineOutput(6,TList::Class()); //My private output
150 // Output slot #7 writes into a AliHFAssociatedTrackCuts container (cuts)
151 DefineOutput(7,AliHFAssociatedTrackCuts::Class()); //My private output
152}
153
154//________________________________________________________________________
155AliAnalysisTaskSED0Correlations::AliAnalysisTaskSED0Correlations(const AliAnalysisTaskSED0Correlations &source):
156 AliAnalysisTaskSE(source),
157 fNPtBinsCorr(source.fNPtBinsCorr),
158 fBinLimsCorr(source.fBinLimsCorr),
159 fPtThreshLow(source.fPtThreshLow),
160 fPtThreshUp(source.fPtThreshUp),
7f221b36 161 fD0Eff(source.fD0Eff),
a2ad7da1 162 fEvents(source.fEvents),
163 fAlreadyFilled(source.fAlreadyFilled),
164 fOutputMass(source.fOutputMass),
165 fOutputCorr(source.fOutputCorr),
166 fOutputStudy(source.fOutputStudy),
167 fNentries(source.fNentries),
168 fCutsD0(source.fCutsD0),
169 fCutsTracks(source.fCutsTracks),
bce70c96 170 fCorrelatorTr(source.fCorrelatorTr),
171 fCorrelatorKc(source.fCorrelatorKc),
172 fCorrelatorK0(source.fCorrelatorK0),
a2ad7da1 173 fReadMC(source.fReadMC),
4d978a4a 174 fRecoTr(source.fRecoTr),
175 fRecoD0(source.fRecoD0),
176 fSelEvType(source.fSelEvType),
bce70c96 177 fMixing(source.fMixing),
a2ad7da1 178 fCounter(source.fCounter),
179 fNPtBins(source.fNPtBins),
180 fFillOnlyD0D0bar(source.fFillOnlyD0D0bar),
181 fIsSelectedCandidate(source.fIsSelectedCandidate),
182 fSys(source.fSys),
8c2d7467 183 fEtaForCorrel(source.fEtaForCorrel),
a2ad7da1 184 fIsRejectSDDClusters(source.fIsRejectSDDClusters),
185 fFillGlobal(source.fFillGlobal)
186{
187 // Copy constructor
188}
189
190//________________________________________________________________________
191AliAnalysisTaskSED0Correlations::~AliAnalysisTaskSED0Correlations()
192{
193 if (fOutputMass) {
194 delete fOutputMass;
195 fOutputMass = 0;
196 }
197 if (fOutputCorr) {
198 delete fOutputCorr;
199 fOutputCorr = 0;
200 }
201 if (fOutputStudy) {
202 delete fOutputStudy;
203 fOutputStudy = 0;
204 }
205 if (fCutsD0) {
206 delete fCutsD0;
207 fCutsD0 = 0;
208 }
209 if (fNentries){
210 delete fNentries;
211 fNentries = 0;
212 }
bce70c96 213 if (fCorrelatorTr) {
214 delete fCorrelatorTr;
215 fCorrelatorTr = 0;
216 }
217 if (fCorrelatorKc) {
218 delete fCorrelatorKc;
219 fCorrelatorKc = 0;
220 }
221 if (fCorrelatorK0) {
222 delete fCorrelatorK0;
223 fCorrelatorK0 = 0;
224 }
225 if (fCounter){
a2ad7da1 226 delete fCounter;
227 fCounter=0;
228 }
229}
230
231//______________________________________________________________________________
232AliAnalysisTaskSED0Correlations& AliAnalysisTaskSED0Correlations::operator=(const AliAnalysisTaskSED0Correlations& orig)
233{
234// Assignment
235 if (&orig == this) return *this; //if address is the same (same object), returns itself
236
237 AliAnalysisTaskSE::operator=(orig); //Uses the AliAnalysisTaskSE operator to assign the inherited part of the class
238 fNPtBinsCorr = orig.fNPtBinsCorr;
239 fBinLimsCorr = orig.fBinLimsCorr;
240 fPtThreshLow = orig.fPtThreshLow;
241 fPtThreshUp = orig.fPtThreshUp;
7f221b36 242 fD0Eff = orig.fD0Eff;
a2ad7da1 243 fEvents = orig.fEvents;
244 fAlreadyFilled = orig.fAlreadyFilled;
245 fOutputMass = orig.fOutputMass;
246 fOutputCorr = orig.fOutputCorr;
247 fOutputStudy = orig.fOutputStudy;
248 fNentries = orig.fNentries;
249 fCutsD0 = orig.fCutsD0;
250 fCutsTracks = orig.fCutsTracks;
bce70c96 251 fCorrelatorTr = orig.fCorrelatorTr;
252 fCorrelatorKc = orig.fCorrelatorKc;
253 fCorrelatorK0 = orig.fCorrelatorK0;
a2ad7da1 254 fReadMC = orig.fReadMC;
4d978a4a 255 fRecoTr = orig.fRecoTr;
256 fRecoD0 = orig.fRecoD0;
257 fSelEvType = orig.fSelEvType;
bce70c96 258 fMixing = orig.fMixing;
a2ad7da1 259 fCounter = orig.fCounter;
260 fNPtBins = orig.fNPtBins;
261 fFillOnlyD0D0bar = orig.fFillOnlyD0D0bar;
262 fIsSelectedCandidate = orig.fIsSelectedCandidate;
263 fSys = orig.fSys;
8c2d7467 264 fEtaForCorrel = orig.fEtaForCorrel;
a2ad7da1 265 fIsRejectSDDClusters = orig.fIsRejectSDDClusters;
266 fFillGlobal = orig.fFillGlobal;
267
268 return *this; //returns pointer of the class
269}
270
271//________________________________________________________________________
272void AliAnalysisTaskSED0Correlations::Init()
273{
274 // Initialization
275
276 if(fDebug > 1) printf("AnalysisTaskSED0Correlations::Init() \n");
277
278 //Copy of cuts objects
279 AliRDHFCutsD0toKpi* copyfCutsD0 = new AliRDHFCutsD0toKpi(*fCutsD0);
280 const char* nameoutput=GetOutputSlot(3)->GetContainer()->GetName();
281 copyfCutsD0->SetName(nameoutput);
282
bce70c96 283 //needer to clear completely the objects inside with Clear() method
a2ad7da1 284 // Post the data
285 PostData(3,copyfCutsD0);
286 PostData(7,fCutsTracks);
287
288 return;
289}
290
291//________________________________________________________________________
292void AliAnalysisTaskSED0Correlations::UserCreateOutputObjects()
293{
294
295 // Create the output container
296 //
297 if(fDebug > 1) printf("AnalysisTaskSED0Correlations::UserCreateOutputObjects() \n");
298
bce70c96 299 //HFCorrelator creation and definition
300 fCorrelatorTr = new AliHFCorrelator("CorrelatorTr",fCutsTracks,fSys);
301 fCorrelatorKc = new AliHFCorrelator("CorrelatorKc",fCutsTracks,fSys);
302 fCorrelatorK0 = new AliHFCorrelator("CorrelatorK0",fCutsTracks,fSys);
7f221b36 303 fCorrelatorTr->SetDeltaPhiInterval(-TMath::Pi()/2+TMath::Pi()/32.,3*TMath::Pi()/2+TMath::Pi()/32.);// set the Delta Phi Interval you want
304 fCorrelatorKc->SetDeltaPhiInterval(-TMath::Pi()/2+TMath::Pi()/32.,3*TMath::Pi()/2+TMath::Pi()/32.);
305 fCorrelatorK0->SetDeltaPhiInterval(-TMath::Pi()/2+TMath::Pi()/32.,3*TMath::Pi()/2+TMath::Pi()/32.);
bce70c96 306 fCorrelatorTr->SetEventMixing(fMixing);// sets the analysis on a single event (kFALSE) or mixed events (kTRUE)
307 fCorrelatorKc->SetEventMixing(fMixing);
308 fCorrelatorK0->SetEventMixing(fMixing);
309 fCorrelatorTr->SetAssociatedParticleType(1);// set 1 for correlations with hadrons, 2 with kaons, 3 with KZeros
310 fCorrelatorKc->SetAssociatedParticleType(2);// set 1 for correlations with hadrons, 2 with kaons, 3 with KZeros
311 fCorrelatorK0->SetAssociatedParticleType(3);// set 1 for correlations with hadrons, 2 with kaons, 3 with KZeros
312 fCorrelatorTr->SetApplyDisplacementCut(2); //0: don't calculate d0; 1: return d0; 2: return d0/d0err
313 fCorrelatorKc->SetApplyDisplacementCut(2);
314 fCorrelatorK0->SetApplyDisplacementCut(0);
315 fCorrelatorTr->SetUseMC(fReadMC);// sets Montecarlo flag
316 fCorrelatorKc->SetUseMC(fReadMC);
317 fCorrelatorK0->SetUseMC(fReadMC);
4d978a4a 318 fCorrelatorTr->SetUseReco(fRecoTr);// sets (if MC analysis) wheter to analyze Reco or Kinem tracks
319 fCorrelatorKc->SetUseReco(fRecoTr);
320 fCorrelatorK0->SetUseReco(fRecoTr);
bce70c96 321 fCorrelatorKc->SetPIDmode(2); //switch for K+/- PID option
322 Bool_t pooldefTr = fCorrelatorTr->DefineEventPool();// method that defines the properties ot the event mixing (zVtx and Multipl. bins)
323 Bool_t pooldefKc = fCorrelatorKc->DefineEventPool();// method that defines the properties ot the event mixing (zVtx and Multipl. bins)
324 Bool_t pooldefK0 = fCorrelatorK0->DefineEventPool();// method that defines the properties ot the event mixing (zVtx and Multipl. bins)
325 if(!pooldefTr) AliInfo("Warning:: Event pool not defined properly");
326 if(!pooldefKc) AliInfo("Warning:: Event pool not defined properly");
327 if(!pooldefK0) AliInfo("Warning:: Event pool not defined properly");
328
a2ad7da1 329 // Several histograms are more conveniently managed in a TList
330 fOutputMass = new TList();
331 fOutputMass->SetOwner();
332 fOutputMass->SetName("listMass");
333
334 fOutputCorr = new TList();
335 fOutputCorr->SetOwner();
336 fOutputCorr->SetName("correlationslist");
337
338 fOutputStudy = new TList();
339 fOutputStudy->SetOwner();
340 fOutputStudy->SetName("MCstudyplots");
341
7f221b36 342 TString nameMass=" ",nameSgn=" ", nameBkg=" ", nameRfl=" ",nameMassWg=" ",nameSgnWg=" ", nameBkgWg=" ", nameRflWg=" ";
a2ad7da1 343
344 for(Int_t i=0;i<fCutsD0->GetNPtBins();i++){
345
346 nameMass="histMass_";
347 nameMass+=i;
7f221b36 348 nameMassWg="histMass_WeigD0Eff_";
349 nameMassWg+=i;
a2ad7da1 350 nameSgn="histSgn_";
351 nameSgn+=i;
7f221b36 352 nameSgnWg="histSgn_WeigD0Eff_";
353 nameSgnWg+=i;
a2ad7da1 354 nameBkg="histBkg_";
355 nameBkg+=i;
7f221b36 356 nameBkgWg="histBkg_WeigD0Eff_";
357 nameBkgWg+=i;
a2ad7da1 358 nameRfl="histRfl_";
359 nameRfl+=i;
7f221b36 360 nameRflWg="histRfl_WeigD0Eff_";
361 nameRflWg+=i;
a2ad7da1 362
363 //histograms of invariant mass distributions
364
365 //MC signal
366 if(fReadMC){
367 TH1F* tmpSt = new TH1F(nameSgn.Data(), "D^{0} invariant mass - MC; M [GeV]; Entries",120,1.5648,2.1648);
7f221b36 368 TH1F* tmpStWg = new TH1F(nameSgnWg.Data(), "D^{0} invariant mass - MC; M [GeV] - weight 1/D0eff; Entries",120,1.5648,2.1648);
a2ad7da1 369 tmpSt->Sumw2();
7f221b36 370 tmpStWg->Sumw2();
a2ad7da1 371
372 //Reflection: histo filled with D0Mass which pass the cut (also) as D0bar and with D0bar which pass (also) the cut as D0
373 TH1F* tmpRt = new TH1F(nameRfl.Data(), "Reflected signal invariant mass - MC; M [GeV]; Entries",120,1.5648,2.1648);
7f221b36 374 TH1F* tmpRtWg = new TH1F(nameRflWg.Data(), "Reflected signal invariant mass - MC - weight 1/D0eff; M [GeV]; Entries",120,1.5648,2.1648);
a2ad7da1 375 TH1F* tmpBt = new TH1F(nameBkg.Data(), "Background invariant mass - MC; M [GeV]; Entries",120,1.5648,2.1648);
7f221b36 376 TH1F* tmpBtWg = new TH1F(nameBkgWg.Data(), "Background invariant mass - MC - weight 1/D0eff; M [GeV]; Entries",120,1.5648,2.1648);
a2ad7da1 377 tmpBt->Sumw2();
7f221b36 378 tmpBtWg->Sumw2();
a2ad7da1 379 tmpRt->Sumw2();
7f221b36 380 tmpRtWg->Sumw2();
a2ad7da1 381 fOutputMass->Add(tmpSt);
7f221b36 382 fOutputMass->Add(tmpStWg);
a2ad7da1 383 fOutputMass->Add(tmpRt);
7f221b36 384 fOutputMass->Add(tmpRtWg);
a2ad7da1 385 fOutputMass->Add(tmpBt);
7f221b36 386 fOutputMass->Add(tmpBtWg);
a2ad7da1 387 }
388
389 //mass
390 TH1F* tmpMt = new TH1F(nameMass.Data(),"D^{0} invariant mass; M [GeV]; Entries",120,1.5648,2.1648);
391 tmpMt->Sumw2();
392 fOutputMass->Add(tmpMt);
7f221b36 393 //mass weighted by 1/D0eff
394 TH1F* tmpMtwg = new TH1F(nameMassWg.Data(),"D^{0} invariant mass - weight 1/D0eff; M [GeV]; Entries",120,1.5648,2.1648);
395 tmpMtwg->Sumw2();
396 fOutputMass->Add(tmpMtwg);
a2ad7da1 397 }
398
399 const char* nameoutput=GetOutputSlot(2)->GetContainer()->GetName();
400
4d978a4a 401 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", 20,-0.5,19.5);
a2ad7da1 402
403 fNentries->GetXaxis()->SetBinLabel(1,"nEventsAnal");
404 fNentries->GetXaxis()->SetBinLabel(2,"nCandSel(Cuts)");
405 fReadMC ? fNentries->GetXaxis()->SetBinLabel(3,"nD0Selected") : fNentries->GetXaxis()->SetBinLabel(3,"Dstar<-D0");
406 fNentries->GetXaxis()->SetBinLabel(4,"nEventsGoodVtxS");
407 fNentries->GetXaxis()->SetBinLabel(5,"ptbin = -1");
408 fNentries->GetXaxis()->SetBinLabel(6,"no daughter");
409 if(fSys==0) fNentries->GetXaxis()->SetBinLabel(7,"nCandSel(Tr)");
410 if(fReadMC && fSys==0){
411 fNentries->GetXaxis()->SetBinLabel(12,"K");
412 fNentries->GetXaxis()->SetBinLabel(13,"Lambda");
413 }
414 fNentries->GetXaxis()->SetBinLabel(14,"Pile-up Rej");
415 fNentries->GetXaxis()->SetBinLabel(15,"N. of 0SMH");
416 if(fSys==1) fNentries->GetXaxis()->SetBinLabel(16,"Nev in centr");
417 if(fIsRejectSDDClusters) fNentries->GetXaxis()->SetBinLabel(17,"SDD-Cls Rej");
418 fNentries->GetXaxis()->SetBinLabel(18,"Phys.Sel.Rej");
4d978a4a 419 fNentries->GetXaxis()->SetBinLabel(19,"nEventsSelected");
420 if(fReadMC) fNentries->GetXaxis()->SetBinLabel(20,"nEvsWithProdMech");
a2ad7da1 421 fNentries->GetXaxis()->SetNdivisions(1,kFALSE);
422
423 fCounter = new AliNormalizationCounter(Form("%s",GetOutputSlot(4)->GetContainer()->GetName()));
424 fCounter->Init();
425
426 CreateCorrelationsObjs(); //creates histos for correlations analysis
427
428 // Post the data
429 PostData(1,fOutputMass);
430 PostData(2,fNentries);
431 PostData(4,fCounter);
432 PostData(5,fOutputCorr);
433 PostData(6,fOutputStudy);
434
435 return;
436}
437
438//________________________________________________________________________
439void AliAnalysisTaskSED0Correlations::UserExec(Option_t */*option*/)
440{
441 // Execute analysis for current event:
442 // heavy flavor candidates association to MC truth
443 //cout<<"I'm in UserExec"<<endl;
444
445
446 //cuts order
447 // printf(" |M-MD0| [GeV] < %f\n",fD0toKpiCuts[0]);
448 // printf(" dca [cm] < %f\n",fD0toKpiCuts[1]);
449 // printf(" cosThetaStar < %f\n",fD0toKpiCuts[2]);
450 // printf(" pTK [GeV/c] > %f\n",fD0toKpiCuts[3]);
451 // printf(" pTpi [GeV/c] > %f\n",fD0toKpiCuts[4]);
452 // printf(" |d0K| [cm] < %f\n",fD0toKpiCuts[5]);
453 // printf(" |d0pi| [cm] < %f\n",fD0toKpiCuts[6]);
454 // printf(" d0d0 [cm^2] < %f\n",fD0toKpiCuts[7]);
455 // printf(" cosThetaPoint > %f\n",fD0toKpiCuts[8]);
456
a2ad7da1 457 AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
458 fEvents++;
459
460 TString bname="D0toKpi";
461
462 TClonesArray *inputArray=0;
463
464 if(!aod && AODEvent() && IsStandardAOD()) {
465 // In case there is an AOD handler writing a standard AOD, use the AOD
466 // event in memory rather than the input (ESD) event.
467 aod = dynamic_cast<AliAODEvent*> (AODEvent());
468 // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
469 // have to taken from the AOD event hold by the AliAODExtension
470 AliAODHandler* aodHandler = (AliAODHandler*)
471 ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
472
473 if(aodHandler->GetExtensions()) {
474 AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
475 AliAODEvent* aodFromExt = ext->GetAOD();
476 inputArray=(TClonesArray*)aodFromExt->GetList()->FindObject(bname.Data());
477 }
478 } else if(aod) {
479 inputArray=(TClonesArray*)aod->GetList()->FindObject(bname.Data());
480 }
481
482 if(!inputArray || !aod) {
483 printf("AliAnalysisTaskSED0Correlations::UserExec: input branch not found!\n");
484 return;
485 }
486
487 // fix for temporary bug in ESDfilter
488 // the AODs with null vertex pointer didn't pass the PhysSel
489 if(!aod->GetPrimaryVertex() || TMath::Abs(aod->GetMagneticField())<0.001) return;
490
491 TClonesArray *mcArray = 0;
492 AliAODMCHeader *mcHeader = 0;
493
494 if(fReadMC) {
495 // load MC particles
496 mcArray = (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName());
497 if(!mcArray) {
498 printf("AliAnalysisTaskSED0Correlations::UserExec: MC particles branch not found!\n");
499 return;
500 }
501
502 // load MC header
503 mcHeader = (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
504 if(!mcHeader) {
505 printf("AliAnalysisTaskSED0Correlations::UserExec: MC header branch not found!\n");
506 return;
507 }
508 }
bce70c96 509
a2ad7da1 510 //histogram filled with 1 for every AOD
511 fNentries->Fill(0);
512 fCounter->StoreEvent(aod,fCutsD0,fReadMC);
513
514 // trigger class for PbPb C0SMH-B-NOPF-ALLNOTRD, C0SMH-B-NOPF-ALL
515 TString trigclass=aod->GetFiredTriggerClasses();
516 if(trigclass.Contains("C0SMH-B-NOPF-ALLNOTRD") || trigclass.Contains("C0SMH-B-NOPF-ALL")) fNentries->Fill(14);
517
518 if(!fCutsD0->IsEventSelected(aod)) {
519 if(fCutsD0->GetWhyRejection()==1) // rejected for pileup
520 fNentries->Fill(13);
521 if(fSys==1 && (fCutsD0->GetWhyRejection()==2 || fCutsD0->GetWhyRejection()==3)) fNentries->Fill(15);
522 if(fCutsD0->GetWhyRejection()==7) fNentries->Fill(17);
523 return;
524 }
525
4d978a4a 526 fNentries->Fill(18); //event selected after selection
527
8c2d7467 528 //Setting PIDResponse for associated tracks
529 fCorrelatorTr->SetPidAssociated();
530 fCorrelatorKc->SetPidAssociated();
531 fCorrelatorK0->SetPidAssociated();
532
4d978a4a 533 //Selection on production type (MC)
534 if(fReadMC && fSelEvType){
535
536 Bool_t isMCeventgood = kFALSE;
537
538 Int_t eventType = mcHeader->GetEventType();
539 Int_t NMCevents = fCutsTracks->GetNofMCEventType();
540
541 for(Int_t k=0; k<NMCevents; k++){
542 Int_t * MCEventType = fCutsTracks->GetMCEventType();
543
544 if(eventType == MCEventType[k]) isMCeventgood= kTRUE;
545 ((TH1D*)fOutputStudy->FindObject("EventTypeMC"))->Fill(eventType);
546 }
547
548 if(NMCevents && !isMCeventgood){
549 std::cout << "The MC event " << eventType << " not interesting for this analysis: skipping" << std::endl;
550 return;
551 }
552 fNentries->Fill(19); //event with particular production type
553
554 } //end of selection
555
556
a2ad7da1 557 // Check the Nb of SDD clusters
558 if (fIsRejectSDDClusters) {
559 Bool_t skipEvent = kFALSE;
560 Int_t ntracks = 0;
561 if (aod) ntracks = aod->GetNTracks();
562 for(Int_t itrack=0; itrack<ntracks; itrack++) { // loop on tacks
563 // ... get the track
564 AliAODTrack * track = aod->GetTrack(itrack);
565 if(TESTBIT(track->GetITSClusterMap(),2) || TESTBIT(track->GetITSClusterMap(),3) ){
566 skipEvent=kTRUE;
567 fNentries->Fill(16);
568 break;
569 }
570 }
571 if (skipEvent) return;
572 }
573
bce70c96 574 //HFCorrelators initialization (for this event)
575 fCorrelatorTr->SetAODEvent(aod); // set the AOD event from which you are processing
576 fCorrelatorKc->SetAODEvent(aod);
577 fCorrelatorK0->SetAODEvent(aod);
578 Bool_t correlatorONTr = fCorrelatorTr->Initialize(); // initialize the pool for event mixing
579 Bool_t correlatorONKc = fCorrelatorKc->Initialize();
580 Bool_t correlatorONK0 = fCorrelatorK0->Initialize();
581 if(!correlatorONTr) {AliInfo("AliHFCorrelator (tracks) didn't initialize the pool correctly or processed a bad event"); return;}
582 if(!correlatorONKc) {AliInfo("AliHFCorrelator (charged K) didn't initialize the pool correctly or processed a bad event"); return;}
583 if(!correlatorONK0) {AliInfo("AliHFCorrelator (K0) didn't initialize the pool correctly or processed a bad event"); return;}
584
585 if(fReadMC) {
586 fCorrelatorTr->SetMCArray(mcArray); // set the TClonesArray *fmcArray for analysis on monte carlo
587 fCorrelatorKc->SetMCArray(mcArray);
588 fCorrelatorK0->SetMCArray(mcArray);
589 }
590
a2ad7da1 591 // AOD primary vertex
592 AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex();
a2ad7da1 593
594 //vtx1->Print();
595 TString primTitle = vtx1->GetTitle();
596 if(primTitle.Contains("VertexerTracks") && vtx1->GetNContributors()>0) {
a2ad7da1 597 fNentries->Fill(3);
598 }
599
600 //Reset flag for tracks distributions fill
601 fAlreadyFilled=kFALSE;
602
bce70c96 603 //***** Loop over D0 candidates *****
a2ad7da1 604 Int_t nInD0toKpi = inputArray->GetEntriesFast();
605 if(fDebug>2) printf("Number of D0->Kpi: %d\n",nInD0toKpi);
606
607 if(fFillGlobal) { //loop on V0 and tracks for each event, to fill Pt distr. and InvMass distr.
608
609 TClonesArray *v0array = (TClonesArray*)aod->GetList()->FindObject("v0s");
610 Int_t pdgCodes[2] = {211,211};
611 Int_t idArrayV0[v0array->GetEntriesFast()][2];
612 for(int iV0=0; iV0<v0array->GetEntriesFast(); iV0++) {
613 for(int j=0; j<2; j++) {idArrayV0[iV0][j]=-2;}
614 AliAODv0 *v0 = (AliAODv0*)v0array->UncheckedAt(iV0);
615 if(SelectV0(v0,vtx1,2,idArrayV0)) { //option 2 = for mass inv plots only
4d978a4a 616 if(fReadMC && fRecoTr && (v0->MatchToMC(310,mcArray,2,pdgCodes)<0)) continue; //310 = K0s, 311 = K0 generico!!
a2ad7da1 617 ((TH2F*)fOutputStudy->FindObject("hK0MassInv"))->Fill(v0->MassK0Short(),v0->Pt()); //invariant mass plot
f80e7bba 618 ((TH1F*)fOutputStudy->FindObject("hist_Pt_K0_AllEv"))->Fill(v0->Pt()); //pT distribution (in all events), K0 case
a2ad7da1 619 }
620 }
621
622 for(Int_t itrack=0; itrack<aod->GetNTracks(); itrack++) { // loop on tacks
623 AliAODTrack * track = aod->GetTrack(itrack);
624 //rejection of tracks
625 if(track->GetID() < 0) continue; //discard negative ID tracks
626 if(track->Pt() < fPtThreshLow.at(0) || track->Pt() > fPtThreshUp.at(0)) continue; //discard tracks outside pt range for hadrons/K
bce70c96 627 if(!fCutsTracks->IsHadronSelected(track) || !fCutsTracks->CheckHadronKinematic(track->Pt(),0.1)) continue; //0.1 = dummy (d0 value, no cut on it for me)
a2ad7da1 628 //pT distribution (in all events), charg and hadr cases
629 ((TH1F*)fOutputStudy->FindObject("hist_Pt_Charg_AllEv"))->Fill(track->Pt());
630 if(fCutsTracks->CheckKaonCompatibility(track,kFALSE,0,2)) ((TH1F*)fOutputStudy->FindObject("hist_Pt_Kcharg_AllEv"))->Fill(track->Pt());
631 }
632
633 } //end of loops for global plot fill
634
635 Int_t nSelectedloose=0,nSelectedtight=0;
4d978a4a 636
637 //RecoD0 case ************************************************
638 if(fRecoD0) {
639
640 for (Int_t iD0toKpi = 0; iD0toKpi < nInD0toKpi; iD0toKpi++) {
641 AliAODRecoDecayHF2Prong *d = (AliAODRecoDecayHF2Prong*)inputArray->UncheckedAt(iD0toKpi);
7f221b36 642
643 if(d->Pt()<2.) continue; //to save time and merging memory...
644
4d978a4a 645 if(d->GetSelectionMap()) if(!d->HasSelectionBit(AliRDHFCuts::kD0toKpiCuts)){
646 fNentries->Fill(2);
647 continue; //skip the D0 from Dstar
a2ad7da1 648 }
649
4d978a4a 650 if (fCutsD0->IsInFiducialAcceptance(d->Pt(),d->Y(421)) ) {
651 nSelectedloose++;
652 nSelectedtight++;
653 if(fSys==0){
654 if(fCutsD0->IsSelected(d,AliRDHFCuts::kTracks,aod))fNentries->Fill(6);
655 }
656 Int_t ptbin=fCutsD0->PtBin(d->Pt());
657 if(ptbin==-1) {fNentries->Fill(4); continue;} //out of bounds
658
659 fIsSelectedCandidate=fCutsD0->IsSelected(d,AliRDHFCuts::kAll,aod); //D0 selected
660 if(!fIsSelectedCandidate) continue;
661
662 //D0 infos
663 Double_t phiD0 = fCorrelatorTr->SetCorrectPhiRange(d->Phi());
664 phiD0 = fCorrelatorKc->SetCorrectPhiRange(d->Phi()); //bad usage, but returns a Double_t...
665 phiD0 = fCorrelatorK0->SetCorrectPhiRange(d->Phi());
666 fCorrelatorTr->SetTriggerParticleProperties(d->Pt(),phiD0,d->Eta()); // sets the parameters of the trigger particles that are needed
667 fCorrelatorKc->SetTriggerParticleProperties(d->Pt(),phiD0,d->Eta());
668 fCorrelatorK0->SetTriggerParticleProperties(d->Pt(),phiD0,d->Eta());
669 fCorrelatorTr->SetD0Properties(d,fIsSelectedCandidate); //sets special properties for D0
670 fCorrelatorKc->SetD0Properties(d,fIsSelectedCandidate);
671 fCorrelatorK0->SetD0Properties(d,fIsSelectedCandidate);
672
673 if(!fReadMC) {
674 if (TMath::Abs(d->Eta())<fEtaForCorrel) CalculateCorrelations(d); //correlations on real data
675 } else { //correlations on MC -> association of selected D0 to MCinfo with MCtruth
676 if (TMath::Abs(d->Eta())<fEtaForCorrel) {
677 Int_t pdgDgD0toKpi[2]={321,211};
678 Int_t labD0 = d->MatchToMC(421,mcArray,2,pdgDgD0toKpi); //return MC particle label if the array corresponds to a D0, -1 if not
679 if (labD0>-1) CalculateCorrelations(d,labD0,mcArray);
680 }
681 }
682
683 FillMassHists(d,mcArray,fCutsD0,fOutputMass);
a2ad7da1 684 }
4d978a4a 685 }
686 }
687 //End RecoD0 case ************************************************
688
689 //MCKineD0 case ************************************************
690 if(fReadMC && !fRecoD0) {
691
692 for (Int_t iPart=0; iPart<mcArray->GetEntriesFast(); iPart++) { //Loop over all the tracks of MCArray
693 AliAODMCParticle* mcPart = dynamic_cast<AliAODMCParticle*>(mcArray->At(iPart));
694 if (!mcPart) {
695 AliWarning("Particle not found in tree, skipping");
696 continue;
697 }
698
699 if(TMath::Abs(mcPart->GetPdgCode()) == 421){ // THIS IS A D0
700 if (fCutsD0->IsInFiducialAcceptance(mcPart->Pt(),mcPart->Y()) ) {
701 nSelectedloose++;
702 nSelectedtight++;
703 if(fSys==0) fNentries->Fill(6);
704 Int_t ptbin=fCutsD0->PtBin(mcPart->Pt());
705 if(ptbin==-1) {fNentries->Fill(4); continue;} //out of bounds
706
707 //D0 infos
708 Double_t phiD0 = fCorrelatorTr->SetCorrectPhiRange(mcPart->Phi());
709 phiD0 = fCorrelatorKc->SetCorrectPhiRange(mcPart->Phi()); //bad usage, but returns a Double_t...
710 phiD0 = fCorrelatorK0->SetCorrectPhiRange(mcPart->Phi());
711 fCorrelatorTr->SetTriggerParticleProperties(mcPart->Pt(),phiD0,mcPart->Eta()); // sets the parameters of the trigger particles that are needed
712 fCorrelatorKc->SetTriggerParticleProperties(mcPart->Pt(),phiD0,mcPart->Eta());
713 fCorrelatorK0->SetTriggerParticleProperties(mcPart->Pt(),phiD0,mcPart->Eta());
714 //fCorrelatorTr->SetD0Properties(mcPart,fIsSelectedCandidate); //needed for D* soft pions rejection, useless in MCKine
715 //fCorrelatorKc->SetD0Properties(mcPart,fIsSelectedCandidate);
716 //fCorrelatorK0->SetD0Properties(mcPart,fIsSelectedCandidate);
717
718 if (TMath::Abs(mcPart->Eta())<fEtaForCorrel) {
719
720 //Removal of D0 from D* feeddown! This solves also the problem of soft pions, now excluded
721 /* Int_t mother = mcPart->GetMother();
722 AliAODMCParticle* mcMoth = dynamic_cast<AliAODMCParticle*>(mcArray->At(mother));
723 if(!mcMoth) continue;
724 if(TMath::Abs(mcMoth->GetPdgCode())==413) continue;
725 */
726 if (mcPart->GetPdgCode()==421) fIsSelectedCandidate = 1;
727 else fIsSelectedCandidate = 2;
728
729 TString fillthis="histSgn_"; fillthis+=ptbin;
730 ((TH1F*)(fOutputMass->FindObject(fillthis)))->Fill(1.864);
731
732 CalculateCorrelationsMCKine(mcPart,mcArray);
733 }
bce70c96 734 }
a2ad7da1 735 }
a2ad7da1 736 }
4d978a4a 737 }
738 //End MCKineD0 case ************************************************
bce70c96 739
740 if(fMixing /* && fAlreadyFilled*/) { // update the pool for Event Mixing, if: enabled, event is ok, at least a SelD0 found! (fAlreadyFilled's role!)
741 Bool_t updatedTr = fCorrelatorTr->PoolUpdate();
742 Bool_t updatedKc = fCorrelatorKc->PoolUpdate();
743 Bool_t updatedK0 = fCorrelatorK0->PoolUpdate();
744 if(!updatedTr || !updatedKc || !updatedK0) AliInfo("Pool was not updated");
745 }
746
a2ad7da1 747 fCounter->StoreCandidates(aod,nSelectedloose,kTRUE);
748 fCounter->StoreCandidates(aod,nSelectedtight,kFALSE);
749
750 // Post the data
751 PostData(1,fOutputMass);
752 PostData(2,fNentries);
753 PostData(4,fCounter);
754 PostData(5,fOutputCorr);
755 PostData(6,fOutputStudy);
756
757 return;
758}
759
760//____________________________________________________________________________
761void AliAnalysisTaskSED0Correlations::FillMassHists(AliAODRecoDecayHF2Prong *part, TClonesArray *arrMC, AliRDHFCutsD0toKpi* cuts, TList *listout){
762 //
763 // function used in UserExec to fill mass histograms:
764 //
765 if (!fIsSelectedCandidate) return;
766
767 if(fDebug>2) cout<<"Candidate selected"<<endl;
768
769 Double_t invmassD0 = part->InvMassD0(), invmassD0bar = part->InvMassD0bar();
770 Int_t ptbin = cuts->PtBin(part->Pt());
771
772 TString fillthis="";
773 Int_t pdgDgD0toKpi[2]={321,211};
774 Int_t labD0=-1;
775 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)
776
777 //count candidates selected by cuts
778 fNentries->Fill(1);
779 //count true D0 selected by cuts
780 if (fReadMC && labD0>=0) fNentries->Fill(2);
781
782 if ((fIsSelectedCandidate==1 || fIsSelectedCandidate==3) && fFillOnlyD0D0bar<2) { //D0
783
784 if(fReadMC){
785 if(labD0>=0) {
786 AliAODMCParticle *partD0 = (AliAODMCParticle*)arrMC->At(labD0);
787 Int_t pdgD0 = partD0->GetPdgCode();
788 if (pdgD0==421){ //D0
789 fillthis="histSgn_";
790 fillthis+=ptbin;
791 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0);
7f221b36 792 fillthis="histSgn_WeigD0Eff_";
793 fillthis+=ptbin;
794 Double_t effD0 = fD0Eff.at(ptbin);
795 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0,1./effD0);
a2ad7da1 796 } else{ //it was a D0bar
797 fillthis="histRfl_";
798 fillthis+=ptbin;
7f221b36 799 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0);
800 fillthis="histRfl_WeigD0Eff_";
801 fillthis+=ptbin;
802 Double_t effD0 = fD0Eff.at(ptbin);
803 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0,1./effD0);
a2ad7da1 804 }
805 } else {//background
806 fillthis="histBkg_";
807 fillthis+=ptbin;
808 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0);
7f221b36 809 fillthis="histBkg_WeigD0Eff_";
810 fillthis+=ptbin;
811 Double_t effD0 = fD0Eff.at(ptbin);
812 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0,1./effD0);
a2ad7da1 813 }
814 }else{
815 fillthis="histMass_";
816 fillthis+=ptbin;
817 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0);
7f221b36 818 fillthis="histMass_WeigD0Eff_";
819 fillthis+=ptbin;
820 Double_t effD0 = fD0Eff.at(ptbin);
821 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0,1./effD0);
a2ad7da1 822 }
823
824 }
825 if (fIsSelectedCandidate>1 && (fFillOnlyD0D0bar==0 || fFillOnlyD0D0bar==2)) { //D0bar
826
827 if(fReadMC){
828 if(labD0>=0) {
829 AliAODMCParticle *partD0 = (AliAODMCParticle*)arrMC->At(labD0);
830 Int_t pdgD0 = partD0->GetPdgCode();
831
832 if (pdgD0==-421){ //D0bar
833 fillthis="histSgn_";
834 fillthis+=ptbin;
835 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar);
7f221b36 836 fillthis="histSgn_WeigD0Eff_";
837 fillthis+=ptbin;
838 Double_t effD0 = fD0Eff.at(ptbin);
839 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar,1./effD0);
a2ad7da1 840 } else{
841 fillthis="histRfl_";
842 fillthis+=ptbin;
843 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar);
7f221b36 844 fillthis="histRfl_WeigD0Eff_";
845 fillthis+=ptbin;
846 Double_t effD0 = fD0Eff.at(ptbin);
847 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar,1./effD0);
a2ad7da1 848 }
849 } else {//background or LS
850 fillthis="histBkg_";
851 fillthis+=ptbin;
852 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar);
7f221b36 853 fillthis="histBkg_WeigD0Eff_";
854 fillthis+=ptbin;
855 Double_t effD0 = fD0Eff.at(ptbin);
856 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar,1./effD0);
a2ad7da1 857 }
858 }else{
859 fillthis="histMass_";
860 fillthis+=ptbin;
861 ((TH1F*)listout->FindObject(fillthis))->Fill(invmassD0bar);
7f221b36 862 fillthis="histMass_WeigD0Eff_";
863 fillthis+=ptbin;
864 Double_t effD0 = fD0Eff.at(ptbin);
865 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar,1./effD0);
a2ad7da1 866 }
867 }
868
869 return;
870}
871
872//________________________________________________________________________
873void AliAnalysisTaskSED0Correlations::Terminate(Option_t */*option*/)
874{
875 // Terminate analysis
876 //
877 if(fDebug > 1) printf("AnalysisTaskSED0Correlations: Terminate() \n");
878
879 fOutputMass = dynamic_cast<TList*> (GetOutputData(1));
880 if (!fOutputMass) {
881 printf("ERROR: fOutputMass not available\n");
882 return;
883 }
884
885 fNentries = dynamic_cast<TH1F*>(GetOutputData(2));
886
887 if(!fNentries){
888 printf("ERROR: fNEntries not available\n");
889 return;
890 }
891
892 fCutsD0 = dynamic_cast<AliRDHFCutsD0toKpi*>(GetOutputData(3));
893 if(!fCutsD0){
894 printf("ERROR: fCuts not available\n");
895 return;
896 }
897
898 fCounter = dynamic_cast<AliNormalizationCounter*>(GetOutputData(4));
899 if (!fCounter) {
900 printf("ERROR: fCounter not available\n");
901 return;
902 }
903 fOutputCorr = dynamic_cast<TList*> (GetOutputData(5));
904 if (!fOutputCorr) {
905 printf("ERROR: fOutputCorr not available\n");
906 return;
907 }
908 fOutputStudy = dynamic_cast<TList*> (GetOutputData(6));
909 if (!fOutputStudy) {
910 printf("ERROR: fOutputStudy not available\n");
911 return;
912 }
913 fCutsTracks = dynamic_cast<AliHFAssociatedTrackCuts*>(GetOutputData(7));
914 if(!fCutsTracks){
915 printf("ERROR: fCutsTracks not available\n");
916 return;
917 }
918
919 return;
920}
921
922//_________________________________________________________________________________________________
923Int_t AliAnalysisTaskSED0Correlations::CheckD0Origin(TClonesArray* arrayMC, AliAODMCParticle *mcPartCandidate) const {
924 //
925 // checking whether the mother of the particles come from a charm or a bottom quark
926 //
bce70c96 927 printf("AliAnalysisTaskSED0Correlations::CheckD0Origin() \n");
a2ad7da1 928
929 Int_t pdgGranma = 0;
930 Int_t mother = 0;
931 mother = mcPartCandidate->GetMother();
932 Int_t abspdgGranma =0;
933 Bool_t isFromB=kFALSE;
934 Bool_t isQuarkFound=kFALSE;
935
936 while (mother > 0){
4d978a4a 937 AliAODMCParticle* mcMoth = dynamic_cast<AliAODMCParticle*>(arrayMC->At(mother));
938 if (mcMoth){
939 pdgGranma = mcMoth->GetPdgCode();
a2ad7da1 940 abspdgGranma = TMath::Abs(pdgGranma);
941 if ((abspdgGranma > 500 && abspdgGranma < 600) || (abspdgGranma > 5000 && abspdgGranma < 6000)){
942 isFromB=kTRUE;
943 }
944 if(abspdgGranma==4 || abspdgGranma==5) isQuarkFound=kTRUE;
4d978a4a 945 mother = mcMoth->GetMother();
a2ad7da1 946 }else{
947 AliError("Failed casting the mother particle!");
948 break;
949 }
950 }
951
952 if(isQuarkFound) {
953 if(isFromB) return 5;
954 else return 4;
955 }
956 else return 1;
957}
958
959//________________________________________________________________________
960void AliAnalysisTaskSED0Correlations::CreateCorrelationsObjs() {
961//
962
963 TString namePlot = "";
964
965 //These for limits in THnSparse (one per bin, same limits).
bce70c96 966 //Vars: DeltaPhi, InvMass, PtTrack, Displacement, DeltaEta
bec72d8c 967 Int_t nBinsPhi[5] = {32,150,6,3,16};
7f221b36 968 Double_t binMinPhi[5] = {-TMath::Pi()/2.+TMath::Pi()/32.,1.6,0.,0.,-1.6}; //is the minimum for all the bins
969 Double_t binMaxPhi[5] = {3*TMath::Pi()/2.+TMath::Pi()/32.,2.2,3.0,3.,1.6}; //is the maximum for all the bins
a2ad7da1 970
bce70c96 971 //Vars: DeltaPhi, InvMass, DeltaEta
972 Int_t nBinsMix[3] = {32,150,16};
7f221b36 973 Double_t binMinMix[3] = {-TMath::Pi()/2.+TMath::Pi()/32.,1.6,-1.6}; //is the minimum for all the bins
974 Double_t binMaxMix[3] = {3*TMath::Pi()/2.+TMath::Pi()/32.,2.2,1.6}; //is the maximum for all the bins
a2ad7da1 975
bce70c96 976 for(Int_t i=0;i<fNPtBinsCorr;i++){
a2ad7da1 977
bce70c96 978 if(!fMixing) {
979 //THnSparse plots: correlations for various invariant mass (MC and data)
980 namePlot="hPhi_K0_Bin";
a2ad7da1 981 namePlot+=i;
982
4d978a4a 983 THnSparseF *hPhiK = new THnSparseF(namePlot.Data(), "Azimuthal correlation; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
bce70c96 984 hPhiK->Sumw2();
985 fOutputCorr->Add(hPhiK);
a2ad7da1 986
bce70c96 987 namePlot="hPhi_Kcharg_Bin";
a2ad7da1 988 namePlot+=i;
989
4d978a4a 990 THnSparseF *hPhiH = new THnSparseF(namePlot.Data(), "Azimuthal correlation; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
bce70c96 991 hPhiH->Sumw2();
992 fOutputCorr->Add(hPhiH);
a2ad7da1 993
bce70c96 994 namePlot="hPhi_Charg_Bin";
a2ad7da1 995 namePlot+=i;
996
4d978a4a 997 THnSparseF *hPhiC = new THnSparseF(namePlot.Data(), "Azimuthal correlation; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
bce70c96 998 hPhiC->Sumw2();
999 fOutputCorr->Add(hPhiC);
a2ad7da1 1000
bce70c96 1001 //histos for c/b origin for D0 (MC only)
1002 if (fReadMC) {
a2ad7da1 1003
bce70c96 1004 //generic origin for tracks
1005 namePlot="hPhi_K0_From_c_Bin";
1006 namePlot+=i;
a2ad7da1 1007
4d978a4a 1008 THnSparseF *hPhiK_c = new THnSparseF(namePlot.Data(), "Azimuthal correlation - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
bce70c96 1009 hPhiK_c->Sumw2();
1010 fOutputCorr->Add(hPhiK_c);
a2ad7da1 1011
bce70c96 1012 namePlot="hPhi_Kcharg_From_c_Bin";
1013 namePlot+=i;
a2ad7da1 1014
4d978a4a 1015 THnSparseF *hPhiH_c = new THnSparseF(namePlot.Data(), "Azimuthal correlation - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
bce70c96 1016 hPhiH_c->Sumw2();
1017 fOutputCorr->Add(hPhiH_c);
a2ad7da1 1018
bce70c96 1019 namePlot="hPhi_Charg_From_c_Bin";
1020 namePlot+=i;
a2ad7da1 1021
4d978a4a 1022 THnSparseF *hPhiC_c = new THnSparseF(namePlot.Data(), "Azimuthal correlation - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
bce70c96 1023 hPhiC_c->Sumw2();
1024 fOutputCorr->Add(hPhiC_c);
1025
1026 namePlot="hPhi_K0_From_b_Bin";
1027 namePlot+=i;
a2ad7da1 1028
4d978a4a 1029 THnSparseF *hPhiK_b = new THnSparseF(namePlot.Data(), "Azimuthal correlation - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
bce70c96 1030 hPhiK_b->Sumw2();
1031 fOutputCorr->Add(hPhiK_b);
a2ad7da1 1032
bce70c96 1033 namePlot="hPhi_Kcharg_From_b_Bin";
1034 namePlot+=i;
a2ad7da1 1035
4d978a4a 1036 THnSparseF *hPhiH_b = new THnSparseF(namePlot.Data(), "Azimuthal correlation - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
bce70c96 1037 hPhiH_b->Sumw2();
1038 fOutputCorr->Add(hPhiH_b);
a2ad7da1 1039
bce70c96 1040 namePlot="hPhi_Charg_From_b_Bin";
1041 namePlot+=i;
a2ad7da1 1042
4d978a4a 1043 THnSparseF *hPhiC_b = new THnSparseF(namePlot.Data(), "Azimuthal correlation - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
bce70c96 1044 hPhiC_b->Sumw2();
1045 fOutputCorr->Add(hPhiC_b);
a2ad7da1 1046
bce70c96 1047 //HF-only tracks (c for c->D0, b for b->D0)
1048 namePlot="hPhi_K0_HF_From_c_Bin";
1049 namePlot+=i;
a2ad7da1 1050
4d978a4a 1051 THnSparseF *hPhiK_HF_c = new THnSparseF(namePlot.Data(), "Azimuthal correlation HF - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
bce70c96 1052 hPhiK_HF_c->Sumw2();
1053 fOutputCorr->Add(hPhiK_HF_c);
a2ad7da1 1054
bce70c96 1055 namePlot="hPhi_Kcharg_HF_From_c_Bin";
1056 namePlot+=i;
a2ad7da1 1057
4d978a4a 1058 THnSparseF *hPhiH_HF_c = new THnSparseF(namePlot.Data(), "Azimuthal correlation HF - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
bce70c96 1059 hPhiH_HF_c->Sumw2();
1060 fOutputCorr->Add(hPhiH_HF_c);
a2ad7da1 1061
bce70c96 1062 namePlot="hPhi_Charg_HF_From_c_Bin";
1063 namePlot+=i;
a2ad7da1 1064
4d978a4a 1065 THnSparseF *hPhiC_HF_c = new THnSparseF(namePlot.Data(), "Azimuthal correlation HF - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
bce70c96 1066 hPhiC_HF_c->Sumw2();
1067 fOutputCorr->Add(hPhiC_HF_c);
a2ad7da1 1068
bce70c96 1069 namePlot="hPhi_K0_HF_From_b_Bin";
1070 namePlot+=i;
a2ad7da1 1071
4d978a4a 1072 THnSparseF *hPhiK_HF_b = new THnSparseF(namePlot.Data(), "Azimuthal correlation HF - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
bce70c96 1073 hPhiK_HF_b->Sumw2();
1074 fOutputCorr->Add(hPhiK_HF_b);
a2ad7da1 1075
bce70c96 1076 namePlot="hPhi_Kcharg_HF_From_b_Bin";
1077 namePlot+=i;
a2ad7da1 1078
4d978a4a 1079 THnSparseF *hPhiH_HF_b = new THnSparseF(namePlot.Data(), "Azimuthal correlation HF - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
bce70c96 1080 hPhiH_HF_b->Sumw2();
1081 fOutputCorr->Add(hPhiH_HF_b);
a2ad7da1 1082
bce70c96 1083 namePlot="hPhi_Charg_HF_From_b_Bin";
1084 namePlot+=i;
a2ad7da1 1085
4d978a4a 1086 THnSparseF *hPhiC_HF_b = new THnSparseF(namePlot.Data(), "Azimuthal correlation HF - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
bce70c96 1087 hPhiC_HF_b->Sumw2();
1088 fOutputCorr->Add(hPhiC_HF_b);
4d978a4a 1089
1090 namePlot="hPhi_K0_NonHF_Bin";
1091 namePlot+=i;
1092
1093 THnSparseF *hPhiK_Non = new THnSparseF(namePlot.Data(), "Azimuthal correlation - Non HF; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
1094 hPhiK_Non->Sumw2();
1095 fOutputCorr->Add(hPhiK_Non);
1096
1097 namePlot="hPhi_Kcharg_NonHF_Bin";
1098 namePlot+=i;
1099
1100 THnSparseF *hPhiH_Non = new THnSparseF(namePlot.Data(), "Azimuthal correlation - Non HF; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
1101 hPhiH_Non->Sumw2();
1102 fOutputCorr->Add(hPhiH_Non);
1103
1104 namePlot="hPhi_Charg_NonHF_Bin";
1105 namePlot+=i;
1106
1107 THnSparseF *hPhiC_Non = new THnSparseF(namePlot.Data(), "Azimuthal correlation - Non HF; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
1108 hPhiC_Non->Sumw2();
1109 fOutputCorr->Add(hPhiC_Non);
bce70c96 1110 }
a2ad7da1 1111
bce70c96 1112 //leading hadron correlations
1113 namePlot="hPhi_Lead_Bin";
a2ad7da1 1114 namePlot+=i;
1115
4d978a4a 1116 THnSparseF *hCorrLead = new THnSparseF(namePlot.Data(), "Leading particle correlations; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix);
bce70c96 1117 hCorrLead->Sumw2();
1118 fOutputCorr->Add(hCorrLead);
a2ad7da1 1119
bce70c96 1120 if (fReadMC) {
1121 namePlot="hPhi_Lead_From_c_Bin";
1122 namePlot+=i;
a2ad7da1 1123
4d978a4a 1124 THnSparseF *hCorrLead_c = new THnSparseF(namePlot.Data(), "Leading particle correlations - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix);
bce70c96 1125 hCorrLead_c->Sumw2();
1126 fOutputCorr->Add(hCorrLead_c);
1127
1128 namePlot="hPhi_Lead_From_b_Bin";
1129 namePlot+=i;
1130
4d978a4a 1131 THnSparseF *hCorrLead_b = new THnSparseF(namePlot.Data(), "Leading particle correlations - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix);
bce70c96 1132 hCorrLead_b->Sumw2();
1133 fOutputCorr->Add(hCorrLead_b);
1134
1135 namePlot="hPhi_Lead_HF_From_c_Bin";
1136 namePlot+=i;
1137
4d978a4a 1138 THnSparseF *hCorrLead_HF_c = new THnSparseF(namePlot.Data(), "Leading particle correlations HF - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix);
bce70c96 1139 hCorrLead_HF_c->Sumw2();
1140 fOutputCorr->Add(hCorrLead_HF_c);
1141
1142 namePlot="hPhi_Lead_HF_From_b_Bin";
1143 namePlot+=i;
1144
4d978a4a 1145 THnSparseF *hCorrLead_HF_b = new THnSparseF(namePlot.Data(), "Leading particle correlations HF - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix);
bce70c96 1146 hCorrLead_HF_b->Sumw2();
1147 fOutputCorr->Add(hCorrLead_HF_b);
4d978a4a 1148
1149 namePlot="hPhi_Lead_NonHF_Bin";
1150 namePlot+=i;
1151
1152 THnSparseF *hCorrLead_Non = new THnSparseF(namePlot.Data(), "Leading particle correlations - Non HF; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix);
1153 hCorrLead_Non->Sumw2();
1154 fOutputCorr->Add(hCorrLead_Non);
bce70c96 1155 }
1156
1157 //pT weighted correlations
1158 namePlot="hPhi_Weig_Bin";
a2ad7da1 1159 namePlot+=i;
bce70c96 1160
4d978a4a 1161 THnSparseF *hCorrWeig = new THnSparseF(namePlot.Data(), "Charged particle correlations (pT weighted); #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix);
bce70c96 1162 fOutputCorr->Add(hCorrWeig);
1163
1164 if (fReadMC) {
1165 namePlot="hPhi_Weig_From_c_Bin";
1166 namePlot+=i;
1167
4d978a4a 1168 THnSparseF *hCorrWeig_c = new THnSparseF(namePlot.Data(), "Charged particle correlations (pT weighted) - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix);
bce70c96 1169 fOutputCorr->Add(hCorrWeig_c);
1170
1171 namePlot="hPhi_Weig_From_b_Bin";
1172 namePlot+=i;
1173
4d978a4a 1174 THnSparseF *hCorrWeig_b = new THnSparseF(namePlot.Data(), "Charged particle correlations (pT weighted) - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix);
bce70c96 1175 fOutputCorr->Add(hCorrWeig_b);
1176
1177 namePlot="hPhi_Weig_HF_From_c_Bin";
1178 namePlot+=i;
1179
4d978a4a 1180 THnSparseF *hCorrWeig_HF_c = new THnSparseF(namePlot.Data(), "Charged particle correlations (pT weighted) HF - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix);
bce70c96 1181 fOutputCorr->Add(hCorrWeig_HF_c);
1182
1183 namePlot="hPhi_Weig_HF_From_b_Bin";
1184 namePlot+=i;
1185
4d978a4a 1186 THnSparseF *hCorrWeig_HF_b = new THnSparseF(namePlot.Data(), "Charged particle correlations (pT weighted) HF - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix);
bce70c96 1187 fOutputCorr->Add(hCorrWeig_HF_b);
4d978a4a 1188
1189 namePlot="hPhi_Weig_NonHF_Bin";
1190 namePlot+=i;
1191
1192 THnSparseF *hCorrWeig_Non = new THnSparseF(namePlot.Data(), "Charged particle correlations (pT weighted) - Non HF; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix);
1193 fOutputCorr->Add(hCorrWeig_Non);
bce70c96 1194 }
a2ad7da1 1195
1196 //pT distribution histos
1197 namePlot = "hist_Pt_Charg_Bin"; namePlot+=i;
1198 TH1F *hPtC = new TH1F(namePlot.Data(), "Charged track pT (in D0 evs); p_{T} (GeV/c)",240,0.,12.);
1199 hPtC->SetMinimum(0);
1200 fOutputStudy->Add(hPtC);
1201
1202 namePlot = "hist_Pt_Kcharg_Bin"; namePlot+=i;
1203 TH1F *hPtH = new TH1F(namePlot.Data(), "Hadrons pT (in D0 evs); p_{T} (GeV/c)",240,0.,12.);
1204 hPtH->SetMinimum(0);
1205 fOutputStudy->Add(hPtH);
1206
f80e7bba 1207 namePlot = "hist_Pt_K0_Bin"; namePlot+=i;
a2ad7da1 1208 TH1F *hPtK = new TH1F(namePlot.Data(), "Kaons pT (in D0 evs); p_{T} (GeV/c)",240,0.,12.);
1209 hPtK->SetMinimum(0);
1210 fOutputStudy->Add(hPtK);
1211
1212 //D* feeddown pions rejection histos
1213 namePlot = "hDstarPions_Bin"; namePlot+=i;
4d978a4a 1214 TH2F *hDstarPions = new TH2F(namePlot.Data(), "Tracks rejected for D* inv.mass cut; # Tracks",2,0.,2.,150,1.6,2.2);
a2ad7da1 1215 hDstarPions->GetXaxis()->SetBinLabel(1,"Not rejected");
1216 hDstarPions->GetXaxis()->SetBinLabel(2,"Rejected");
1217 hDstarPions->SetMinimum(0);
bce70c96 1218 fOutputStudy->Add(hDstarPions);
1219
1220 }
1221
1222 if(fMixing) {
1223 //THnSparse plots for event mixing!
1224 namePlot="hPhi_K0_Bin";
1225 namePlot+=i;namePlot+="_EvMix";
1226
4d978a4a 1227 THnSparseF *hPhiK_EvMix = new THnSparseF(namePlot.Data(), "Az. corr. EvMix; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix);
bce70c96 1228 hPhiK_EvMix->Sumw2();
1229 fOutputCorr->Add(hPhiK_EvMix);
1230
1231 namePlot="hPhi_Kcharg_Bin";
1232 namePlot+=i;namePlot+="_EvMix";
1233
4d978a4a 1234 THnSparseF *hPhiH_EvMix = new THnSparseF(namePlot.Data(), "Az. corr. EvMix; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix);
bce70c96 1235 hPhiH_EvMix->Sumw2();
1236 fOutputCorr->Add(hPhiH_EvMix);
1237
1238 namePlot="hPhi_Charg_Bin";
1239 namePlot+=i;namePlot+="_EvMix";
1240
4d978a4a 1241 THnSparseF *hPhiC_EvMix = new THnSparseF(namePlot.Data(), "Az. corr. EvMix; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix);
bce70c96 1242 hPhiC_EvMix->Sumw2();
1243 fOutputCorr->Add(hPhiC_EvMix);
1244 }
a2ad7da1 1245 }
1246
1247 //out of bin loop
da975030 1248 if(!fMixing) {
1249 TH1F *hCountC = new TH1F("hist_Count_Charg", "Charged track counter; # Tracks",100,0.,100.);
1250 hCountC->SetMinimum(0);
1251 fOutputStudy->Add(hCountC);
bce70c96 1252
da975030 1253 TH1F *hCountH = new TH1F("hist_Count_Kcharg", "Hadrons counter; # Tracks",20,0.,20.);
1254 hCountH->SetMinimum(0);
1255 fOutputStudy->Add(hCountH);
bce70c96 1256
da975030 1257 TH1F *hCountK = new TH1F("hist_Count_K0", "Kaons counter; # Tracks",20,0.,20.);
1258 hCountK->SetMinimum(0);
1259 fOutputStudy->Add(hCountK);
1260 }
bce70c96 1261
4d978a4a 1262 if (fReadMC) {
1263 TH1D *hEventTypeMC = new TH1D("EventTypeMC","EventTypeMC",100,-0.5,99.5);
1264 fOutputStudy->Add(hEventTypeMC);
1265 }
1266
a2ad7da1 1267 if (fFillGlobal) { //all-events plots
1268 //pt distributions
1269 TH1F *hPtCAll = new TH1F("hist_Pt_Charg_AllEv", "Charged track pT (All); p_{T} (GeV/c)",240,0.,12.);
1270 hPtCAll->SetMinimum(0);
1271 fOutputStudy->Add(hPtCAll);
1272
1273 TH1F *hPtHAll = new TH1F("hist_Pt_Kcharg_AllEv", "Hadrons pT (All); p_{T} (GeV/c)",240,0.,12.);
1274 hPtHAll->SetMinimum(0);
1275 fOutputStudy->Add(hPtHAll);
1276
bce70c96 1277 TH1F *hPtKAll = new TH1F("hist_Pt_K0_AllEv", "Kaons pT (All); p_{T} (GeV/c)",240,0.,12.);
a2ad7da1 1278 hPtKAll->SetMinimum(0);
1279 fOutputStudy->Add(hPtKAll);
1280
da975030 1281 if(!fMixing) {
1282 //phi distributions
bec72d8c 1283 TH1F *hPhiDistCAll = new TH1F("hist_PhiDistr_Charg", "Charged track phi distr. (All); p_{T} (GeV/c)",64,0,6.283);
da975030 1284 hPhiDistCAll->SetMinimum(0);
1285 fOutputStudy->Add(hPhiDistCAll);
1286
bec72d8c 1287 TH1F *hPhiDistHAll = new TH1F("hist_PhiDistr_Kcharg", "Hadrons phi distr. (All); p_{T} (GeV/c)",64,0,6.283);
da975030 1288 hPhiDistHAll->SetMinimum(0);
1289 fOutputStudy->Add(hPhiDistHAll);
1290
bec72d8c 1291 TH1F *hPhiDistKAll = new TH1F("hist_PhiDistr_K0", "Kaons phi distr. (All); p_{T} (GeV/c)",64,0,6.283);
da975030 1292 hPhiDistKAll->SetMinimum(0);
1293 fOutputStudy->Add(hPhiDistKAll);
1294
bec72d8c 1295 TH1F *hPhiDistDAll = new TH1F("hist_PhiDistr_D0", "D^{0} phi distr. (All); p_{T} (GeV/c)",64,0,6.283);
da975030 1296 hPhiDistDAll->SetMinimum(0);
1297 fOutputStudy->Add(hPhiDistDAll);
1298 }
bce70c96 1299
a2ad7da1 1300 //K0 Invariant Mass plots
1301 TH2F *hK0MassInv = new TH2F("hK0MassInv", "K0 invariant mass; Invariant mass (MeV/c^{2}); pT (GeV/c)",200,0.4,0.6,100,0.,10.);
1302 hK0MassInv->SetMinimum(0);
1303 fOutputStudy->Add(hK0MassInv);
1304 }
1305
1306 //for MC analysis only
da975030 1307 for(Int_t i=0;i<fNPtBinsCorr;i++) {
a2ad7da1 1308
da975030 1309 if (fReadMC && !fMixing) {
a2ad7da1 1310
1311 //displacement histos
f80e7bba 1312 namePlot="histDispl_K0_Bin"; namePlot+=i;
a2ad7da1 1313 TH1F *hDisplK = new TH1F(namePlot.Data(), "Kaons Displacement; DCA",150,0.,0.15);
1314 hDisplK->SetMinimum(0);
1315 fOutputStudy->Add(hDisplK);
1316
f80e7bba 1317 namePlot="histDispl_K0_HF_Bin"; namePlot+=i;
a2ad7da1 1318 TH1F *hDisplK_HF = new TH1F(namePlot.Data(), "Kaons Displacement (from HF decay only); DCA",150,0.,0.15);
1319 hDisplK_HF->SetMinimum(0);
1320 fOutputStudy->Add(hDisplK_HF);
1321
1322 namePlot="histDispl_Kcharg_Bin"; namePlot+=i;
1323 TH1F *hDisplHadr = new TH1F(namePlot.Data(), "Hadrons Displacement; DCA",150,0.,0.15);
1324 hDisplHadr->SetMinimum(0);
1325 fOutputStudy->Add(hDisplHadr);
1326
1327 namePlot="histDispl_Kcharg_HF_Bin"; namePlot+=i;
1328 TH1F *hDisplHadr_HF = new TH1F(namePlot.Data(), "Hadrons Displacement (from HF decay only); DCA",150,0.,0.15);
1329 hDisplHadr_HF->SetMinimum(0);
1330 fOutputStudy->Add(hDisplHadr_HF);
1331
1332 namePlot="histDispl_Charg_Bin"; namePlot+=i;
1333 TH1F *hDisplCharg = new TH1F(namePlot.Data(), "Charged tracks Displacement; DCA",150,0.,0.15);
1334 hDisplCharg->SetMinimum(0);
1335 fOutputStudy->Add(hDisplCharg);
1336
1337 namePlot="histDispl_Charg_HF_Bin"; namePlot+=i;
1338 TH1F *hDisplCharg_HF = new TH1F(namePlot.Data(), "Charged tracks Displacement (from HF decay only); DCA",150,0.,0.15);
1339 hDisplCharg_HF->SetMinimum(0);
1340 fOutputStudy->Add(hDisplCharg_HF);
1341
f80e7bba 1342 namePlot="histDispl_K0_From_c_Bin"; namePlot+=i;
a2ad7da1 1343 TH1F *hDisplK_c = new TH1F(namePlot.Data(), "Kaons Displacement - c origin; DCA",150,0.,0.15);
1344 hDisplK_c->SetMinimum(0);
1345 fOutputStudy->Add(hDisplK_c);
1346
f80e7bba 1347 namePlot="histDispl_K0_HF_From_c_Bin"; namePlot+=i;
a2ad7da1 1348 TH1F *hDisplK_HF_c = new TH1F(namePlot.Data(), "Kaons Displacement (from HF decay only) - c origin; DCA",150,0.,0.15);
1349 hDisplK_HF_c->SetMinimum(0);
1350 fOutputStudy->Add(hDisplK_HF_c);
1351
1352 namePlot="histDispl_Kcharg_From_c_Bin"; namePlot+=i;
1353 TH1F *hDisplHadr_c = new TH1F(namePlot.Data(), "Hadrons Displacement - c origin; DCA",150,0.,0.15);
1354 hDisplHadr_c->SetMinimum(0);
1355 fOutputStudy->Add(hDisplHadr_c);
1356
1357 namePlot="histDispl_Kcharg_HF_From_c_Bin"; namePlot+=i;
1358 TH1F *hDisplHadr_HF_c = new TH1F(namePlot.Data(), "Hadrons Displacement (from HF decay only) - c origin; DCA",150,0.,0.15);
1359 hDisplHadr_HF_c->SetMinimum(0);
1360 fOutputStudy->Add(hDisplHadr_HF_c);
1361
1362 namePlot="histDispl_Charg_From_c_Bin"; namePlot+=i;
1363 TH1F *hDisplCharg_c = new TH1F(namePlot.Data(), "Charged tracks Displacement - c origin; DCA",150,0.,0.15);
1364 hDisplCharg_c->Sumw2();
1365 hDisplCharg_c->SetMinimum(0);
1366 fOutputStudy->Add(hDisplCharg_c);
1367
1368 namePlot="histDispl_Charg_HF_From_c_Bin"; namePlot+=i;
1369 TH1F *hDisplCharg_HF_c = new TH1F(namePlot.Data(), "Charged tracks Displacement (from HF decay only) - c origin; DCA",150,0.,0.15);
1370 hDisplCharg_HF_c->SetMinimum(0);
1371 fOutputStudy->Add(hDisplCharg_HF_c);
1372
f80e7bba 1373 namePlot="histDispl_K0_From_b_Bin"; namePlot+=i;
a2ad7da1 1374 TH1F *hDisplK_b = new TH1F(namePlot.Data(), "Kaons Displacement - b origin; DCA",150,0.,0.15);
1375 hDisplK_b->SetMinimum(0);
1376 fOutputStudy->Add(hDisplK_b);
1377
f80e7bba 1378 namePlot="histDispl_K0_HF_From_b_Bin"; namePlot+=i;
a2ad7da1 1379 TH1F *hDisplK_HF_b = new TH1F(namePlot.Data(), "Kaons Displacement (from HF decay only) - b origin; DCA",150,0.,0.15);
1380 hDisplK_HF_b->SetMinimum(0);
1381 fOutputStudy->Add(hDisplK_HF_b);
1382
1383 namePlot="histDispl_Kcharg_From_b_Bin"; namePlot+=i;
1384 TH1F *hDisplHadr_b = new TH1F(namePlot.Data(), "Hadrons Displacement - b origin; DCA",150,0.,0.15);
1385 hDisplHadr_b->SetMinimum(0);
1386 fOutputStudy->Add(hDisplHadr_b);
1387
1388 namePlot="histDispl_Kcharg_HF_From_b_Bin"; namePlot+=i;
1389 TH1F *hDisplHadr_HF_b = new TH1F(namePlot.Data(), "Hadrons Displacement (from HF decay only) - b origin; DCA",150,0.,0.15);
1390 hDisplHadr_HF_b->SetMinimum(0);
1391 fOutputStudy->Add(hDisplHadr_HF_b);
1392
1393 namePlot="histDispl_Charg_From_b_Bin"; namePlot+=i;
1394 TH1F *hDisplCharg_b = new TH1F(namePlot.Data(), "Charged tracks Displacement - b origin; DCA",150,0.,0.15);
1395 hDisplCharg_b->SetMinimum(0);
1396 fOutputStudy->Add(hDisplCharg_b);
1397
1398 namePlot="histDispl_Charg_HF_From_b_Bin"; namePlot+=i;
1399 TH1F *hDisplCharg_HF_b = new TH1F(namePlot.Data(), "Charged tracks Displacement (from HF decay only) - b origin; DCA",150,0.,0.15);
1400 hDisplCharg_HF_b->SetMinimum(0);
1401 fOutputStudy->Add(hDisplCharg_HF_b);
1402
1403 //origin of tracks histos
1404 namePlot="histOrig_Charg_Bin"; namePlot+=i;
1405 TH1F *hOrigin_Charm = new TH1F(namePlot.Data(), "Origin of charged tracks",9,0.,9.);
1406 hOrigin_Charm->SetMinimum(0);
1407 hOrigin_Charm->GetXaxis()->SetBinLabel(1,"Not HF");
1408 hOrigin_Charm->GetXaxis()->SetBinLabel(2,"D->#");
1409 hOrigin_Charm->GetXaxis()->SetBinLabel(3,"D->X->#");
4d978a4a 1410 hOrigin_Charm->GetXaxis()->SetBinLabel(4,"c hadr.");
1411 hOrigin_Charm->GetXaxis()->SetBinLabel(5,"B->#");
1412 hOrigin_Charm->GetXaxis()->SetBinLabel(6,"B->X-># (X!=D)");
1413 hOrigin_Charm->GetXaxis()->SetBinLabel(7,"B->D->#");
1414 hOrigin_Charm->GetXaxis()->SetBinLabel(8,"B->D->X->#");
a2ad7da1 1415 hOrigin_Charm->GetXaxis()->SetBinLabel(9,"b hadr.");
1416 fOutputStudy->Add(hOrigin_Charm);
1417
1418 namePlot="histOrig_Kcharg_Bin"; namePlot+=i;
1419 TH1F *hOrigin_Kcharg = new TH1F(namePlot.Data(), "Origin of hadrons",9,0.,9.);
1420 hOrigin_Kcharg->SetMinimum(0);
1421 hOrigin_Kcharg->GetXaxis()->SetBinLabel(1,"Not HF");
1422 hOrigin_Kcharg->GetXaxis()->SetBinLabel(2,"D->#");
1423 hOrigin_Kcharg->GetXaxis()->SetBinLabel(3,"D->X->#");
4d978a4a 1424 hOrigin_Kcharg->GetXaxis()->SetBinLabel(4,"c hadr.");
1425 hOrigin_Kcharg->GetXaxis()->SetBinLabel(5,"B->#");
1426 hOrigin_Kcharg->GetXaxis()->SetBinLabel(6,"B->X-># (X!=D)");
1427 hOrigin_Kcharg->GetXaxis()->SetBinLabel(7,"B->D->#");
1428 hOrigin_Kcharg->GetXaxis()->SetBinLabel(8,"B->D->X->#");
a2ad7da1 1429 hOrigin_Kcharg->GetXaxis()->SetBinLabel(9,"b hadr.");
1430 fOutputStudy->Add(hOrigin_Kcharg);
1431
f80e7bba 1432 namePlot="histOrig_K0_Bin"; namePlot+=i;
a2ad7da1 1433 TH1F *hOrigin_K = new TH1F(namePlot.Data(), "Origin of kaons",9,0.,9.);
1434 hOrigin_K->SetMinimum(0);
1435 hOrigin_K->GetXaxis()->SetBinLabel(1,"Not HF");
1436 hOrigin_K->GetXaxis()->SetBinLabel(2,"D->#");
1437 hOrigin_K->GetXaxis()->SetBinLabel(3,"D->X->#");
4d978a4a 1438 hOrigin_K->GetXaxis()->SetBinLabel(4,"c hadr.");
1439 hOrigin_K->GetXaxis()->SetBinLabel(5,"B->#");
1440 hOrigin_K->GetXaxis()->SetBinLabel(6,"B->X-># (X!=D)");
1441 hOrigin_K->GetXaxis()->SetBinLabel(7,"B->D->#");
1442 hOrigin_K->GetXaxis()->SetBinLabel(8,"B->D->X->#");
a2ad7da1 1443 hOrigin_K->GetXaxis()->SetBinLabel(9,"b hadr.");
1444 fOutputStudy->Add(hOrigin_K);
da975030 1445 }
a2ad7da1 1446
da975030 1447 if (fReadMC) {
a2ad7da1 1448 //origin of D0 histos
1449 namePlot="histOrig_D0_Bin"; namePlot+=i;
1450 TH1F *hOrigin_D0 = new TH1F(namePlot.Data(), "Origin of D0",2,0.,2.);
1451 hOrigin_D0->SetMinimum(0);
1452 hOrigin_D0->GetXaxis()->SetBinLabel(1,"From c");
1453 hOrigin_D0->GetXaxis()->SetBinLabel(2,"From b");
1454 fOutputStudy->Add(hOrigin_D0);
1455 }
1456 }
a2ad7da1 1457}
1458
1459//________________________________________________________________________
bce70c96 1460void AliAnalysisTaskSED0Correlations::CalculateCorrelations(AliAODRecoDecayHF2Prong* d, Int_t labD0, TClonesArray* mcArray) {
a2ad7da1 1461//
1462// Method for correlations D0-hadrons study
1463//
bce70c96 1464 Int_t N_Charg = 0, N_KCharg = 0, N_Kaons = 0;
1465 Double_t mD0, mD0bar;
4d978a4a 1466 Int_t origD0 = 0, PDGD0 = 0, ptbin = 0;
a2ad7da1 1467 d->InvMassD0(mD0,mD0bar);
4d978a4a 1468 Double_t mInv[2] = {mD0, mD0bar};
1469 ptbin = PtBinCorr(d->Pt());
1470
a2ad7da1 1471 if(ptbin < 0) return;
a2ad7da1 1472
bec72d8c 1473 //Fill of D0 phi distribution
4d978a4a 1474 if (!fMixing) ((TH1F*)fOutputStudy->FindObject("hist_PhiDistr_D0"))->Fill(d->Phi());
bec72d8c 1475
bce70c96 1476 //Origin of D0
1477 TString orig="";
1478 if(fReadMC) {
1479 origD0=CheckD0Origin(mcArray,(AliAODMCParticle*)mcArray->At(labD0));
1480 PDGD0 = ((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode();
1481 switch (CheckD0Origin(mcArray,(AliAODMCParticle*)mcArray->At(labD0))) {
1482 case 4:
1483 orig = "_From_c";
1484 ((TH1F*)fOutputStudy->FindObject(Form("histOrig_D0_Bin%d",ptbin)))->Fill(0.);
1485 break;
1486 case 5:
1487 orig = "_From_b";
1488 ((TH1F*)fOutputStudy->FindObject(Form("histOrig_D0_Bin%d",ptbin)))->Fill(1.);
1489 break;
1490 default:
1491 return;
1492 }
1493 }
1494
4d978a4a 1495 Double_t highPt = 0; Double_t lead[4] = {0,0,0,1}; //infos for leading particle (pt,deltaphi)
a2ad7da1 1496
bce70c96 1497 //loop over the tracks in the pool
ad8d2dfd 1498 Bool_t execPoolTr = fCorrelatorTr->ProcessEventPool(); //pool is ready? (only in ME, in SE returns kFALSE)
1499 Bool_t execPoolKc = fCorrelatorKc->ProcessEventPool(); //pool is ready? (only in ME, in SE returns kFALSE)
1500 Bool_t execPoolK0 = fCorrelatorK0->ProcessEventPool(); //pool is ready? (only in ME, in SE returns kFALSE)
bce70c96 1501
1502 Int_t NofEventsinPool = 1;
ad8d2dfd 1503 if(fMixing) {
1504 NofEventsinPool = fCorrelatorTr->GetNofEventsInPool();
1505 if(!execPoolTr) {
1506 AliInfo("Mixed event analysis: track pool is not ready");
1507 NofEventsinPool = 0;
1508 }
1509 }
1510
1511 //Charged tracks
bce70c96 1512 for (Int_t jMix =0; jMix < NofEventsinPool; jMix++) {// loop on events in the pool; if it is SE analysis, stops at one (index not needed there)
1513 Bool_t analyzetracksTr = fCorrelatorTr->ProcessAssociatedTracks(jMix);// process all the tracks in the aodEvent, by applying the selection cuts
ad8d2dfd 1514 if(!analyzetracksTr) {
bce70c96 1515 AliInfo("AliHFCorrelator::Cannot process the track array");
1516 continue;
1517 }
ad8d2dfd 1518
bce70c96 1519 for(Int_t iTrack = 0; iTrack<fCorrelatorTr->GetNofTracks(); iTrack++){ // looping on track candidates
1520
1521 Bool_t runcorrelation = fCorrelatorTr->Correlate(iTrack);
1522 if(!runcorrelation) continue;
1523
1524 AliReducedParticle* track = fCorrelatorTr->GetAssociatedParticle();
1525
1526 if(!fMixing) {
7f221b36 1527 if(!track->CheckSoftPi()) { //removal of soft pions
bce70c96 1528 if (fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(1.,mD0);
1529 if (fIsSelectedCandidate >= 2) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(1.,mD0bar);
1530 continue;
1531 } else { //not a soft pion
1532 if (fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(0.,mD0);
1533 if (fIsSelectedCandidate >= 2) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(0.,mD0bar);
7f221b36 1534 }
4d978a4a 1535 Int_t idDaughs[2] = {((AliVTrack*)d->GetDaughter(0))->GetID(),((AliVTrack*)d->GetDaughter(1))->GetID()}; //IDs of daughters to be skipped
1536 if(track->GetID() == idDaughs[0] || track->GetID() == idDaughs[1]) continue; //discards daughters of candidate
ad8d2dfd 1537 }
1538 if(track->Pt() < fPtThreshLow.at(ptbin) || track->Pt() > fPtThreshUp.at(ptbin)) continue; //discard tracks outside pt range for hadrons/K
4d978a4a 1539
7f221b36 1540 Double_t effTr = track->GetWeight(); //extract track efficiency
1541 Double_t effD0 = fD0Eff.at(ptbin);
1542 Double_t eff = effTr*effD0;
1543
4d978a4a 1544 FillSparsePlots(mcArray,mInv,origD0,PDGD0,track,ptbin,kTrack,1./eff); //fills for charged tracks, weight = 1./eff
bce70c96 1545
1546 if(!fMixing) N_Charg++;
1547
1548 //retrieving leading info...
1549 if(track->Pt() > highPt) {
da975030 1550 if(fReadMC && track->GetLabel()<1) continue;
8c2d7467 1551 if(fReadMC && !(AliAODMCParticle*)mcArray->At(track->GetLabel())) continue;
bce70c96 1552 lead[0] = fCorrelatorTr->GetDeltaPhi();
1553 lead[1] = fCorrelatorTr->GetDeltaEta();
1554 lead[2] = fReadMC ? CheckTrackOrigin(mcArray,(AliAODMCParticle*)mcArray->At(track->GetLabel())) : 0;
4d978a4a 1555 lead[3] = 1./track->GetWeight(); //weight is 1./efficiency
bce70c96 1556 highPt = track->Pt();
a2ad7da1 1557 }
bce70c96 1558
1559 } // end of tracks loop
ad8d2dfd 1560 } //end of event loop for fCorrelatorTr
1561
1562 if(fMixing) {
1563 NofEventsinPool = fCorrelatorKc->GetNofEventsInPool();
1564 if(!execPoolKc) {
1565 AliInfo("Mixed event analysis: K+/- pool is not ready");
1566 NofEventsinPool = 0;
1567 }
1568 }
1569
1570 //Charged Kaons loop
4d978a4a 1571 for (Int_t jMix = 0; jMix < NofEventsinPool; jMix++) {// loop on events in the pool; if it is SE analysis, stops at one (index not needed there)
ad8d2dfd 1572 Bool_t analyzetracksKc = fCorrelatorKc->ProcessAssociatedTracks(jMix);
1573 if(!analyzetracksKc) {
1574 AliInfo("AliHFCorrelator::Cannot process the K+/- array");
1575 continue;
1576 }
bce70c96 1577
bce70c96 1578 for(Int_t iTrack = 0; iTrack<fCorrelatorKc->GetNofTracks(); iTrack++){ // looping on charged kaons candidates
1579
1580 Bool_t runcorrelation = fCorrelatorKc->Correlate(iTrack);
1581 if(!runcorrelation) continue;
1582
1583 AliReducedParticle* kCharg = fCorrelatorKc->GetAssociatedParticle();
1584
1585 if(!fMixing) {
ad8d2dfd 1586 if(!kCharg->CheckSoftPi()) { //removal of soft pions
bce70c96 1587 if (fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(1.,mD0);
1588 if (fIsSelectedCandidate >= 2) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(1.,mD0bar);
1589 continue;
1590 } else {
1591 if (fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(0.,mD0);
1592 if (fIsSelectedCandidate >= 2) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(0.,mD0bar);
a2ad7da1 1593 }
4d978a4a 1594 Int_t idDaughs[2] = {((AliVTrack*)d->GetDaughter(0))->GetID(),((AliVTrack*)d->GetDaughter(1))->GetID()}; //IDs of daughters to be skipped
1595 if(kCharg->GetID() == idDaughs[0] || kCharg->GetID() == idDaughs[1]) continue; //discards daughters of candidate
ad8d2dfd 1596 }
1597 if(kCharg->Pt() < fPtThreshLow.at(ptbin) || kCharg->Pt() > fPtThreshUp.at(ptbin)) continue; //discard tracks outside pt range for hadrons/K
4d978a4a 1598
1599 FillSparsePlots(mcArray,mInv,origD0,PDGD0,kCharg,ptbin,kKCharg); //fills for charged tracks
a2ad7da1 1600
bce70c96 1601 if(!fMixing) N_KCharg++;
a2ad7da1 1602
bce70c96 1603 } // end of charged kaons loop
ad8d2dfd 1604 } //end of event loop for fCorrelatorKc
1605
1606 if(fMixing) {
1607 NofEventsinPool = fCorrelatorK0->GetNofEventsInPool();
1608 if(!execPoolK0) {
1609 AliInfo("Mixed event analysis: K0 pool is not ready");
1610 NofEventsinPool = 0;
1611 }
1612 }
1613
1614 //K0 loop
1615 for (Int_t jMix =0; jMix < NofEventsinPool; jMix++) {// loop on events in the pool; if it is SE analysis, stops at one (index not needed there)
1616 Bool_t analyzetracksK0 = fCorrelatorK0->ProcessAssociatedTracks(jMix);
1617 if(!analyzetracksK0) {
1618 AliInfo("AliHFCorrelator::Cannot process the K0 array");
1619 continue;
1620 }
a2ad7da1 1621
bce70c96 1622 for(Int_t iTrack = 0; iTrack<fCorrelatorK0->GetNofTracks(); iTrack++){ // looping on k0 candidates
1623
1624 Bool_t runcorrelation = fCorrelatorK0->Correlate(iTrack);
1625 if(!runcorrelation) continue;
1626
1627 AliReducedParticle* k0 = fCorrelatorK0->GetAssociatedParticle();
a2ad7da1 1628
bce70c96 1629 if(k0->Pt() < fPtThreshLow.at(ptbin) || k0->Pt() > fPtThreshUp.at(ptbin)) continue; //discard tracks outside pt range for hadrons/K
4d978a4a 1630
1631 FillSparsePlots(mcArray,mInv,origD0,PDGD0,k0,ptbin,kK0); //fills for charged tracks
bce70c96 1632
1633 if(!fMixing) N_Kaons++;
1634
1635 } // end of charged kaons loop
ad8d2dfd 1636 } //end of event loop for fCorrelatorK0
bce70c96 1637
4d978a4a 1638 Double_t fillSpLeadD0[3] = {lead[0],mD0,lead[1]};
1639 Double_t fillSpLeadD0bar[3] = {lead[0],mD0bar,lead[1]};
1640
bce70c96 1641 //leading track correlations fill
ad8d2dfd 1642 if(!fMixing) {
1643 if(fReadMC) {
1644 if(((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode()==421) { //D0
4d978a4a 1645 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(fillSpLeadD0,lead[3]); //c and b D0
1646 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadD0,lead[3]); //c or b D0
1647 if(origD0==4&&(int)lead[2]>=1&&(int)lead[2]<=3) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadD0,lead[3]);
1648 if(origD0==5&&(int)lead[2]>=4&&(int)lead[2]<=8) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadD0,lead[3]);
1649 if((int)lead[2]==0) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_NonHF_Bin%d",ptbin)))->Fill(fillSpLeadD0,lead[3]); //non HF
ad8d2dfd 1650 }
1651 if(((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode()==-421) { //D0bar
4d978a4a 1652 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(fillSpLeadD0bar,lead[3]);
1653 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadD0bar,lead[3]); //c or b D0
1654 if(origD0==4&&(int)lead[2]>=1&&(int)lead[2]<=3) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadD0bar,lead[3]);
1655 if(origD0==5&&(int)lead[2]>=4&&(int)lead[2]<=8) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadD0bar,lead[3]);
1656 if((int)lead[2]==0) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_NonHF_Bin%d",ptbin)))->Fill(fillSpLeadD0bar,lead[3]); //non HF
ad8d2dfd 1657 }
1658 } else {
4d978a4a 1659 if(fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(fillSpLeadD0,lead[3]);
1660 if(fIsSelectedCandidate == 2 || fIsSelectedCandidate == 3) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(fillSpLeadD0bar,lead[3]);
a2ad7da1 1661 }
1662
ad8d2dfd 1663 //Fill of count histograms
1664 if (!fAlreadyFilled) {
1665 ((TH1F*)fOutputStudy->FindObject("hist_Count_Charg"))->Fill(N_Charg);
1666 ((TH1F*)fOutputStudy->FindObject("hist_Count_Kcharg"))->Fill(N_KCharg);
1667 ((TH1F*)fOutputStudy->FindObject("hist_Count_K0"))->Fill(N_Kaons);
1668 }
bce70c96 1669 }
a2ad7da1 1670
bce70c96 1671 fAlreadyFilled=kTRUE; //at least a D0 analyzed in the event; distribution plots already filled
a2ad7da1 1672
bce70c96 1673}
1674
1675//________________________________________________________________________
4d978a4a 1676void AliAnalysisTaskSED0Correlations::CalculateCorrelationsMCKine(AliAODMCParticle* d, TClonesArray* mcArray) {
1677//
1678// Method for correlations D0-hadrons study
1679//
1680 Int_t N_Charg = 0, N_KCharg = 0, N_Kaons = 0;
1681 Double_t mD0 = 1.864, mD0bar = 1.864;
1682 Double_t mInv[2] = {mD0, mD0bar};
1683 Int_t origD0 = 0, PDGD0 = 0;
1684 Int_t ptbin = PtBinCorr(d->Pt());
1685
1686 if(ptbin < 0) return;
1687
1688 //Fill of D0 phi distribution
1689 if (!fMixing) ((TH1F*)fOutputStudy->FindObject("hist_PhiDistr_D0"))->Fill(d->Phi());
1690
1691 //Origin of D0
1692 TString orig="";
1693 origD0=CheckD0Origin(mcArray,d);
1694 PDGD0 = d->GetPdgCode();
1695 switch (CheckD0Origin(mcArray,d)) {
1696 case 4:
1697 orig = "_From_c";
1698 ((TH1F*)fOutputStudy->FindObject(Form("histOrig_D0_Bin%d",ptbin)))->Fill(0.);
1699 break;
1700 case 5:
1701 orig = "_From_b";
1702 ((TH1F*)fOutputStudy->FindObject(Form("histOrig_D0_Bin%d",ptbin)))->Fill(1.);
1703 break;
1704 default:
1705 return;
1706 }
1707
1708 Double_t highPt = 0; Double_t lead[3] = {0,0,0}; //infos for leading particle (pt,deltaphi)
1709
1710 //loop over the tracks in the pool
1711 Bool_t execPoolTr = fCorrelatorTr->ProcessEventPool(); //pool is ready? (only in ME, in SE returns kFALSE)
1712 Bool_t execPoolKc = fCorrelatorKc->ProcessEventPool(); //pool is ready? (only in ME, in SE returns kFALSE)
1713 Bool_t execPoolK0 = fCorrelatorK0->ProcessEventPool(); //pool is ready? (only in ME, in SE returns kFALSE)
1714
1715 Int_t NofEventsinPool = 1;
1716 if(fMixing) {
1717 NofEventsinPool = fCorrelatorTr->GetNofEventsInPool();
1718 if(!execPoolTr) {
1719 AliInfo("Mixed event analysis: track pool is not ready");
1720 NofEventsinPool = 0;
1721 }
1722 }
1723
1724 //Charged tracks
1725 for (Int_t jMix =0; jMix < NofEventsinPool; jMix++) {// loop on events in the pool; if it is SE analysis, stops at one (index not needed there)
1726
1727 Bool_t analyzetracksTr = fCorrelatorTr->ProcessAssociatedTracks(jMix);// process all the tracks in the aodEvent, by applying the selection cuts
1728 if(!analyzetracksTr) {
1729 AliInfo("AliHFCorrelator::Cannot process the track array");
1730 continue;
1731 }
1732
1733 for(Int_t iTrack = 0; iTrack<fCorrelatorTr->GetNofTracks(); iTrack++){ // looping on track candidates
1734
1735 Bool_t runcorrelation = fCorrelatorTr->Correlate(iTrack);
1736 if(!runcorrelation) continue;
1737
1738 AliReducedParticle* track = fCorrelatorTr->GetAssociatedParticle();
1739 if(track->GetLabel()<0) continue;
1740 if(track->Pt() < fPtThreshLow.at(ptbin) || track->Pt() > fPtThreshUp.at(ptbin)) continue; //discard tracks outside pt range for hadrons/K
1741 if(track->Pt() < 0.3 || TMath::Abs(track->Eta())>0.8) continue; //discard tracks outside barrel (since it's kinematic MC and produces tracks all over rapidity region
1742 if(!fMixing) N_Charg++;
1743
1744 AliAODMCParticle *trkMC = (AliAODMCParticle*)mcArray->At(track->GetLabel());
1745 if(!trkMC) continue;
1746
1747 if (!trkMC->IsPhysicalPrimary()) continue; //reject material budget, or other fake tracks
1748 if (IsDDaughter(d,trkMC,mcArray)) continue;
1749 FillSparsePlots(mcArray,mInv,origD0,PDGD0,track,ptbin,kTrack); //fills for charged tracks
1750
1751 //retrieving leading info...
1752 if(track->Pt() > highPt) {
1753 lead[0] = fCorrelatorTr->GetDeltaPhi();
1754 lead[1] = fCorrelatorTr->GetDeltaEta();
1755 lead[2] = fReadMC ? CheckTrackOrigin(mcArray,trkMC) : 0;
1756 highPt = track->Pt();
1757 }
1758
1759 } // end of tracks loop
1760 } //end of event loop for fCorrelatorTr
1761
1762 if(fMixing) {
1763 NofEventsinPool = fCorrelatorKc->GetNofEventsInPool();
1764 if(!execPoolKc) {
1765 AliInfo("Mixed event analysis: K+/- pool is not ready");
1766 NofEventsinPool = 0;
1767 }
1768 }
1769
1770 //Charged Kaons loop
1771 for (Int_t jMix =0; jMix < NofEventsinPool; jMix++) {// loop on events in the pool; if it is SE analysis, stops at one (index not needed there)
1772 Bool_t analyzetracksKc = fCorrelatorKc->ProcessAssociatedTracks(jMix);
1773 if(!analyzetracksKc) {
1774 AliInfo("AliHFCorrelator::Cannot process the K+/- array");
1775 continue;
1776 }
1777
1778 for(Int_t iTrack = 0; iTrack<fCorrelatorKc->GetNofTracks(); iTrack++){ // looping on charged kaons candidates
1779
1780 Bool_t runcorrelation = fCorrelatorKc->Correlate(iTrack);
1781 if(!runcorrelation) continue;
1782
1783 AliReducedParticle* kCharg = fCorrelatorKc->GetAssociatedParticle();
1784 if(kCharg->GetLabel()<1) continue;
1785 if(kCharg->Pt() < fPtThreshLow.at(ptbin) || kCharg->Pt() > fPtThreshUp.at(ptbin)) continue; //discard tracks outside pt range for hadrons/K
1786 if(TMath::Abs(kCharg->Eta())>0.8) continue; //discard tracks outside barrel (since it's kinematic MC and produces tracks all over rapidity region
1787 if(!fMixing) N_KCharg++;
1788
1789 AliAODMCParticle *kChargMC = (AliAODMCParticle*)mcArray->At(kCharg->GetLabel());
1790 if(!kChargMC) continue;
1791
1792 if (IsDDaughter(d,kChargMC,mcArray)) continue;
1793 FillSparsePlots(mcArray,mInv,origD0,PDGD0,kCharg,ptbin,kKCharg); //fills for charged tracks
1794
1795 } // end of charged kaons loop
1796 } //end of event loop for fCorrelatorKc
1797
1798 if(fMixing) {
1799 NofEventsinPool = fCorrelatorK0->GetNofEventsInPool();
1800 if(!execPoolK0) {
1801 AliInfo("Mixed event analysis: K0 pool is not ready");
1802 NofEventsinPool = 0;
1803 }
1804 }
1805
1806 //K0 loop
1807 for (Int_t jMix =0; jMix < NofEventsinPool; jMix++) {// loop on events in the pool; if it is SE analysis, stops at one (index not needed there)
1808 Bool_t analyzetracksK0 = fCorrelatorK0->ProcessAssociatedTracks(jMix);
1809 if(!analyzetracksK0) {
1810 AliInfo("AliHFCorrelator::Cannot process the K0 array");
1811 continue;
1812 }
1813
1814 for(Int_t iTrack = 0; iTrack<fCorrelatorK0->GetNofTracks(); iTrack++){ // looping on k0 candidates
1815
1816 Bool_t runcorrelation = fCorrelatorK0->Correlate(iTrack);
1817 if(!runcorrelation) continue;
1818
1819 AliReducedParticle* k0 = fCorrelatorK0->GetAssociatedParticle();
1820 if(k0->GetLabel()<1) continue;
1821 if(k0->Pt() < fPtThreshLow.at(ptbin) || k0->Pt() > fPtThreshUp.at(ptbin)) continue; //discard tracks outside pt range for hadrons/K
1822 if(TMath::Abs(k0->Eta())>0.8) continue; //discard tracks outside barrel (since it's kinematic MC and produces tracks all over rapidity region
1823
1824 AliAODMCParticle *k0MC = (AliAODMCParticle*)mcArray->At(k0->GetLabel());
1825 if(!k0MC) continue;
1826
1827 if (IsDDaughter(d,k0MC,mcArray)) continue;
1828 FillSparsePlots(mcArray,mInv,origD0,PDGD0,k0,ptbin,kK0); //fills for charged tracks
1829
1830 if(!fMixing) N_Kaons++;
1831
1832 } // end of charged kaons loop
1833 } //end of event loop for fCorrelatorK0
1834
1835 Double_t fillSpLeadMC[3] = {lead[0],mD0,lead[1]}; //mD0 = mD0bar = 1.864
1836
1837 //leading track correlations fill
1838 if(!fMixing) {
1839 if(d->GetPdgCode()==421) { //D0
1840 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(fillSpLeadMC); //c and b D0
1841 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadMC); //c or b D0
1842 if(origD0==4&&(int)lead[2]>=1&&(int)lead[2]<=3) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadMC);
1843 if(origD0==5&&(int)lead[2]>=4&&(int)lead[2]<=8) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadMC);
1844 if((int)lead[2]==0) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_NonHF_Bin%d",ptbin)))->Fill(fillSpLeadMC); //non HF
1845 }
1846 if(d->GetPdgCode()==-421) { //D0bar
1847 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(fillSpLeadMC);
1848 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadMC); //c or b D0
1849 if(origD0==4&&(int)lead[2]>=1&&(int)lead[2]<=3) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadMC);
1850 if(origD0==5&&(int)lead[2]>=4&&(int)lead[2]<=8) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadMC);
1851 if((int)lead[2]==0) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_NonHF_Bin%d",ptbin)))->Fill(fillSpLeadMC); //non HF
1852 }
1853
1854 //Fill of count histograms
1855 if (!fAlreadyFilled) {
1856 ((TH1F*)fOutputStudy->FindObject("hist_Count_Charg"))->Fill(N_Charg);
1857 ((TH1F*)fOutputStudy->FindObject("hist_Count_Kcharg"))->Fill(N_KCharg);
1858 ((TH1F*)fOutputStudy->FindObject("hist_Count_K0"))->Fill(N_Kaons);
1859 }
1860
1861 fAlreadyFilled=kTRUE; //at least a D0 analyzed in the event; distribution plots already filled
1862 }
1863
1864}
1865
1866//________________________________________________________________________
1867void AliAnalysisTaskSED0Correlations::FillSparsePlots(TClonesArray* mcArray, Double_t mInv[], Int_t origD0, Int_t PdgD0, AliReducedParticle* track, Int_t ptbin, Int_t type, Double_t wg) {
bce70c96 1868 //
1869 //fills the THnSparse for correlations, calculating the variables
1870 //
1871
1872 //Initialization of variables
bce70c96 1873 Double_t mD0, mD0bar, deltaphi = 0., deltaeta = 0.;
4d978a4a 1874 mD0 = mInv[0];
1875 mD0bar = mInv[1];
bce70c96 1876
da975030 1877 if (fReadMC && track->GetLabel()<1) return;
1878 if (fReadMC && !(AliAODMCParticle*)mcArray->At(track->GetLabel())) return;
bce70c96 1879 Double_t ptTrack = track->Pt();
1880 Double_t d0Track = type!=kK0 ? track->GetImpPar() : 0.;
1881 Double_t phiTr = track->Phi();
1882 Double_t origTr = fReadMC ? CheckTrackOrigin(mcArray,(AliAODMCParticle*)mcArray->At(track->GetLabel())) : 0;
1883
1884 TString part = "", orig = "";
1885
1886 switch (type) {
1887 case(kTrack): {
1888 part = "Charg";
1889 deltaphi = fCorrelatorTr->GetDeltaPhi();
1890 deltaeta = fCorrelatorTr->GetDeltaEta();
1891 break;
1892 }
1893 case(kKCharg): {
1894 part = "Kcharg";
1895 deltaphi = fCorrelatorKc->GetDeltaPhi();
1896 deltaeta = fCorrelatorKc->GetDeltaEta();
1897 break;
1898 }
1899 case(kK0): {
1900 part = "K0";
1901 deltaphi = fCorrelatorK0->GetDeltaPhi();
1902 deltaeta = fCorrelatorK0->GetDeltaEta();
1903 break;
a2ad7da1 1904 }
bce70c96 1905 }
bce70c96 1906
1907 if(fMixing == kSE) {
1908
ad8d2dfd 1909 //Fixes limits; needed to include overflow into THnSparse projections!
1910 Double_t pTorig = track->Pt();
1911 Double_t d0orig = track->GetImpPar();
4d978a4a 1912 Double_t ptLim_Sparse = ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->GetAxis(2)->GetXmax(); //all plots have same axes...
1913 Double_t displLim_Sparse = ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->GetAxis(3)->GetXmax();
1914 Double_t EtaLim_Sparse = ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->GetAxis(4)->GetXmax();
ad8d2dfd 1915 if(ptTrack > ptLim_Sparse) ptTrack = ptLim_Sparse-0.01;
1916 if(d0Track > displLim_Sparse) d0Track = (displLim_Sparse-0.001);
1917 if(deltaeta > EtaLim_Sparse) deltaeta = EtaLim_Sparse-0.01;
1918 if(deltaeta < -EtaLim_Sparse) deltaeta = -EtaLim_Sparse+0.01;
1919
1920 //variables for filling histos
1921 Double_t fillSpPhiD0[5] = {deltaphi,mD0,ptTrack,d0Track,deltaeta};
1922 Double_t fillSpPhiD0bar[5] = {deltaphi,mD0bar,ptTrack,d0Track,deltaeta};
4d978a4a 1923 Double_t fillSpWeigD0[3] = {deltaphi,mD0,deltaeta};
1924 Double_t fillSpWeigD0bar[3] = {deltaphi,mD0bar,deltaeta};
a2ad7da1 1925
bce70c96 1926 if(fReadMC == 0) {
1927 //sparse fill for data (tracks, K+-, K0) + weighted
1928 if(fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) { //D0
4d978a4a 1929 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0,wg);
1930 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig_Bin%d",ptbin)))->Fill(fillSpWeigD0,pTorig*wg);
a2ad7da1 1931 }
bce70c96 1932 if(fIsSelectedCandidate == 2 || fIsSelectedCandidate == 3) { //D0bar
4d978a4a 1933 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0bar,wg);
1934 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig_Bin%d",ptbin)))->Fill(fillSpWeigD0bar,pTorig*wg);
a2ad7da1 1935 }
1936 if(!fAlreadyFilled) {
bce70c96 1937 ((TH1F*)fOutputStudy->FindObject(Form("hist_Pt_%s_Bin%d",part.Data(),ptbin)))->Fill(pTorig);
1938 ((TH1F*)fOutputStudy->FindObject(Form("hist_PhiDistr_%s",part.Data())))->Fill(phiTr);
a2ad7da1 1939 }
bce70c96 1940 }
a2ad7da1 1941
bce70c96 1942 if(fReadMC) {
a2ad7da1 1943
bce70c96 1944 if(origD0==4) {orig = "_From_c";} else {orig = "_From_b";}
a2ad7da1 1945
bce70c96 1946 //sparse fill for data (tracks, K+-, K0) + weighted
1947 if(PdgD0==421) { //D0 (from MCTruth)
4d978a4a 1948 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0,wg);
1949 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0,wg);
1950 if(origD0==4&&origTr>=1&&origTr<=3) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_HF%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0,wg);
1951 if(origD0==5&&origTr>=4&&origTr<=8) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_HF%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0,wg);
1952 if(origTr==0) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_NonHF_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0,wg);
1953 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig_Bin%d",ptbin)))->Fill(fillSpWeigD0,pTorig*wg);
1954 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpWeigD0,pTorig*wg);
1955 if(origD0==4&&origTr>=1&&origTr<=3) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpWeigD0,pTorig*wg);
1956 if(origD0==5&&origTr>=4&&origTr<=8) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpWeigD0,pTorig*wg);
1957 if(origTr==0) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig_NonHF_Bin%d",ptbin)))->Fill(fillSpWeigD0,pTorig*wg);
a2ad7da1 1958 }
bce70c96 1959 if(PdgD0==-421) { //D0bar
4d978a4a 1960 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0bar,wg);
1961 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0bar,wg);
1962 if(origD0==4&&origTr>=1&&origTr<=3) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_HF%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0bar,wg);
1963 if(origD0==5&&origTr>=4&&origTr<=8) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_HF%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0bar,wg);
1964 if(origTr==0) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_NonHF_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0bar,wg);
1965 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig_Bin%d",ptbin)))->Fill(fillSpWeigD0bar,pTorig*wg);
1966 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpWeigD0bar,pTorig*wg);
1967 if(origD0==4&&origTr>=1&&origTr<=3) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpWeigD0bar,pTorig*wg);
1968 if(origD0==5&&origTr>=4&&origTr<=8) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpWeigD0bar,pTorig*wg);
1969 if(origTr==0) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig_NonHF_Bin%d",ptbin)))->Fill(fillSpWeigD0bar,pTorig*wg);
bce70c96 1970 }
a2ad7da1 1971 if(!fAlreadyFilled) {
bce70c96 1972 ((TH1F*)fOutputStudy->FindObject(Form("histDispl_%s_Bin%d",part.Data(),ptbin)))->Fill(d0orig); //Fills displacement histos
4d978a4a 1973 if (origTr>=1&&origTr<=8) ((TH1F*)fOutputStudy->FindObject(Form("histDispl_%s_HF_Bin%d",part.Data(),ptbin)))->Fill(d0orig);
1974 if (origTr>=1&&origTr<=8) ((TH1F*)fOutputStudy->FindObject(Form("histDispl_%s_HF%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(d0orig);
bce70c96 1975 ((TH1F*)fOutputStudy->FindObject(Form("histDispl_%s%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(d0orig); //Fills displacement histos
1976 ((TH1F*)fOutputStudy->FindObject(Form("hist_Pt_%s_Bin%d",part.Data(),ptbin)))->Fill(pTorig);
1977 ((TH1F*)fOutputStudy->FindObject(Form("histOrig_%s_Bin%d",part.Data(),ptbin)))->Fill(origTr);
1978 ((TH1F*)fOutputStudy->FindObject(Form("hist_PhiDistr_%s",part.Data())))->Fill(phiTr);
a2ad7da1 1979 }
bce70c96 1980 }//end MC case
a2ad7da1 1981
bce70c96 1982 } //end of SE fill
a2ad7da1 1983
bce70c96 1984 if(fMixing == kME) {
a2ad7da1 1985
ad8d2dfd 1986 //Fixes limits; needed to include overflow into THnSparse projections!
4d978a4a 1987 Double_t EtaLim_Sparse = ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d_EvMix",ptbin)))->GetAxis(2)->GetXmax();
ad8d2dfd 1988 if(deltaeta > EtaLim_Sparse) deltaeta = EtaLim_Sparse-0.01;
1989 if(deltaeta < -EtaLim_Sparse) deltaeta = -EtaLim_Sparse+0.01;
1990
1991 //variables for filling histos
1992 Double_t fillSpPhiD0[3] = {deltaphi,mD0,deltaeta};
1993 Double_t fillSpPhiD0bar[3] = {deltaphi,mD0bar,deltaeta};
a2ad7da1 1994
bce70c96 1995 if(fReadMC == 0) {
1996 //sparse fill for data (tracks, K+-, K0)
4d978a4a 1997 if(fIsSelectedCandidate == 1||fIsSelectedCandidate == 3) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d_EvMix",part.Data(),ptbin)))->Fill(fillSpPhiD0,wg);
1998 if(fIsSelectedCandidate == 2||fIsSelectedCandidate == 3) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d_EvMix",part.Data(),ptbin)))->Fill(fillSpPhiD0bar,wg);
bce70c96 1999 }
2000 if(fReadMC == 1) {
2001 //sparse fill for data (tracks, K+-, K0)
4d978a4a 2002 if(PdgD0==421) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d_EvMix",part.Data(),ptbin)))->Fill(fillSpPhiD0,wg);
2003 if(PdgD0==-421) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d_EvMix",part.Data(),ptbin)))->Fill(fillSpPhiD0bar,wg);
bce70c96 2004 }//end MC case
a2ad7da1 2005
bce70c96 2006 } //end of ME fill
2007
2008 return;
a2ad7da1 2009}
2010
2011//_________________________________________________________________________________________________
2012Int_t AliAnalysisTaskSED0Correlations::CheckTrackOrigin(TClonesArray* arrayMC, AliAODMCParticle *mcPartCandidate) const {
2013 //
2014 // checks on particle (#) origin:
2015 // 0) Not HF
2016 // 1) D->#
2017 // 2) D->X->#
4d978a4a 2018 // 3) c hadronization
2019 // 4) B->#
2020 // 5) B->X-># (X!=D)
2021 // 6) B->D->#
2022 // 7) B->D->X->#
a2ad7da1 2023 // 8) b hadronization
2024 //
bce70c96 2025 if(fDebug>2) printf("AliAnalysisTaskSED0Correlations::CheckTrkOrigin() \n");
a2ad7da1 2026
2027 Int_t pdgGranma = 0;
2028 Int_t mother = 0;
2029 mother = mcPartCandidate->GetMother();
2030 Int_t istep = 0;
2031 Int_t abspdgGranma =0;
2032 Bool_t isFromB=kFALSE;
2033 Bool_t isDdaugh=kFALSE;
2034 Bool_t isDchaindaugh=kFALSE;
2035 Bool_t isBdaugh=kFALSE;
2036 Bool_t isBchaindaugh=kFALSE;
2037 Bool_t isQuarkFound=kFALSE;
2038
4d978a4a 2039 if (mother<0) return -1;
2040 while (mother >= 0){
a2ad7da1 2041 istep++;
4d978a4a 2042 AliAODMCParticle* mcMoth = dynamic_cast<AliAODMCParticle*>(arrayMC->At(mother));
2043 if (mcMoth){
2044 pdgGranma = mcMoth->GetPdgCode();
a2ad7da1 2045 abspdgGranma = TMath::Abs(pdgGranma);
2046 if ((abspdgGranma > 500 && abspdgGranma < 600) || (abspdgGranma > 5000 && abspdgGranma < 6000)){
2047 isBchaindaugh=kTRUE;
2048 if(istep==1) isBdaugh=kTRUE;
2049 }
2050 if ((abspdgGranma > 400 && abspdgGranma < 500) || (abspdgGranma > 4000 && abspdgGranma < 5000)){
2051 isDchaindaugh=kTRUE;
2052 if(istep==1) isDdaugh=kTRUE;
2053 }
2054 if(abspdgGranma==4 || abspdgGranma==5) {isQuarkFound=kTRUE; if(abspdgGranma==5) isFromB = kTRUE;}
4d978a4a 2055 mother = mcMoth->GetMother();
a2ad7da1 2056 }else{
2057 AliError("Failed casting the mother particle!");
4d978a4a 2058 return -1;
a2ad7da1 2059 }
2060 }
2061
2062 //decides what to return based on the flag status
2063 if(isQuarkFound) {
2064 if(!isFromB) { //charm
2065 if(isDdaugh) return 1; //charm immediate
2066 else if(isDchaindaugh) return 2; //charm chain
4d978a4a 2067 else return 3; //charm hadronization
a2ad7da1 2068 }
2069 else { //beauty
4d978a4a 2070 if(isBdaugh) return 4; //b immediate
a2ad7da1 2071 else if(isBchaindaugh) { //b chain
2072 if(isDchaindaugh) {
4d978a4a 2073 if(isDdaugh) return 6; //d immediate
2074 return 7; //d chain
a2ad7da1 2075 }
4d978a4a 2076 else return 5; //b, not d
a2ad7da1 2077 }
2078 else return 8; //b hadronization
2079 }
2080 }
4d978a4a 2081 else if(!isDdaugh && !isDchaindaugh && !isBdaugh && !isBchaindaugh) return 0; //no HF decay
2082 //in this case, it's !isQuarkFound, but not in 100% cases it's a non HF particle!
2083 //rarely you can find a D/B meson which comes from a -1! It isn't a Non-HF, in that case! And I'll return -1...
2084
2085 return -1; //some problem spotted
2086}
2087//________________________________________________________________________
2088Bool_t AliAnalysisTaskSED0Correlations::IsDDaughter(AliAODMCParticle* d, AliAODMCParticle* track, TClonesArray* mcArray) const {
2089
2090 //Daughter removal in MCKine case
2091 Bool_t isDaughter = kFALSE;
2092 Int_t labelD0 = d->GetLabel();
2093
2094 Int_t mother = track->GetMother();
2095
2096 //Loop on the mothers to find the D0 label (it must be the trigger D0, not a generic D0!)
2097 while (mother > 0){
2098 AliAODMCParticle* mcMoth = dynamic_cast<AliAODMCParticle*>(mcArray->At(mother)); //it's the mother of the track!
2099 if (mcMoth){
2100 if (mcMoth->GetLabel() == labelD0) isDaughter = kTRUE;
2101 mother = mcMoth->GetMother(); //goes back by one
2102 } else{
2103 AliError("Failed casting the mother particle!");
2104 break;
2105 }
2106 }
2107
2108 return isDaughter;
a2ad7da1 2109}
2110
a2ad7da1 2111//________________________________________________________________________
2112Int_t AliAnalysisTaskSED0Correlations::PtBinCorr(Double_t pt) const {
2113 //
2114 //give the pt bin where the pt lies.
2115 //
2116 Int_t ptbin=-1;
2117 if(pt<fBinLimsCorr.at(0)) return ptbin; //out of bounds
2118
2119 Int_t i = 0;
2120 while(pt>fBinLimsCorr.at(i)) {ptbin=i; i++;}
2121
2122 return ptbin;
2123}
2124
a2ad7da1 2125//---------------------------------------------------------------------------
2126Bool_t AliAnalysisTaskSED0Correlations::SelectV0(AliAODv0* v0, AliAODVertex *vtx, Int_t opt, Int_t idArrayV0[][2]) const
2127{
2128 //
2129 // Selection for K0 hypotheses
2130 // options: 1 = selects mass invariant about 3 sigma inside the peak + threshold of 0.3 GeV
2131 // 2 = no previous selections
2132
2133 if(!fCutsTracks->IsKZeroSelected(v0,vtx)) return kFALSE;
2134
2135 AliAODTrack *v0Daug1 = (AliAODTrack*)v0->GetDaughter(0);
2136 AliAODTrack *v0Daug2 = (AliAODTrack*)v0->GetDaughter(1);
2137
2138 if(opt==1) { //additional cuts for correlations (V0 has to be closer than 3 sigma from K0 mass)
bce70c96 2139 if(TMath::Abs(v0->MassK0Short()-0.4976) > 3*0.004) return kFALSE;
a2ad7da1 2140 }
2141
2142 //This part removes double counting for swapped tracks!
2143 Int_t i = 0; //while loop (until the last-written entry pair of ID!
2144 while(idArrayV0[i][0]!=-2 && idArrayV0[i][1]!=-2) {
2145 if((v0Daug1->GetID()==idArrayV0[i][0] && v0Daug2->GetID()==idArrayV0[i][1])||
2146 (v0Daug1->GetID()==idArrayV0[i][1] && v0Daug2->GetID()==idArrayV0[i][0])) return kFALSE;
2147 i++;
2148 }
2149 idArrayV0[i][0]=v0Daug1->GetID();
2150 idArrayV0[i][1]=v0Daug2->GetID();
2151
2152 return kTRUE;
2153}
2154
2155//________________________________________________________________________
2156void AliAnalysisTaskSED0Correlations::PrintBinsAndLimits() {
2157
2158 cout << "--------------------------\n";
2159 cout << "PtBins = " << fNPtBinsCorr << "\n";
2160 cout << "PtBin limits--------------\n";
2161 for (int i=0; i<fNPtBinsCorr; i++) {
2162 cout << "Bin "<<i+1<<" = "<<fBinLimsCorr.at(i)<<" to "<<fBinLimsCorr.at(i+1)<<"\n";
2163 }
2164 cout << "\n--------------------------\n";
2165 cout << "PtBin tresh. tracks low---\n";
2166 for (int i=0; i<fNPtBinsCorr; i++) {
2167 cout << fPtThreshLow.at(i) << "; ";
2168 }
2169 cout << "PtBin tresh. tracks up----\n";
2170 for (int i=0; i<fNPtBinsCorr; i++) {
2171 cout << fPtThreshUp.at(i) << "; ";
2172 }
7f221b36 2173 cout << "D0 efficiencies:----\n";
2174 for (int i=0; i<fNPtBinsCorr; i++) {
2175 cout << fD0Eff.at(i) << "; ";
2176 }
a2ad7da1 2177 cout << "\n--------------------------\n";
8c2d7467 2178 cout << "D0 Eta cut for Correl = "<<fEtaForCorrel<<"\n";
2179 cout << "--------------------------\n";
4d978a4a 2180 cout << "MC Truth = "<<fReadMC<<" - MC Reco: Trk = "<<fRecoTr<<", D0 = "<<fRecoD0<<"\n";
2181 cout << "--------------------------\n";
2182 cout << "Sel of Event tpye (PP,GS,FE,...)= "<<fSelEvType<<"\n";
a2ad7da1 2183 cout << "--------------------------\n";
bce70c96 2184 cout << "Ev Mixing = "<<fMixing<<"\n";
2185 cout << "--------------------------\n";
a2ad7da1 2186}
2187