]> 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
Event Mixing for Balance Function Psi Task (PhiCorrelation style)
[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(),
69 fEvents(0),
70 fAlreadyFilled(kFALSE),
71 fOutputMass(0),
72 fOutputCorr(0),
73 fOutputStudy(0),
74 fNentries(0),
75 fCutsD0(0),
76 fCutsTracks(0),
bce70c96 77 fCorrelatorTr(0),
78 fCorrelatorKc(0),
79 fCorrelatorK0(0),
a2ad7da1 80 fReadMC(0),
bce70c96 81 fMixing(kFALSE),
a2ad7da1 82 fCounter(0),
83 fNPtBins(1),
84 fFillOnlyD0D0bar(0),
85 fIsSelectedCandidate(0),
86 fSys(0),
87 fIsRejectSDDClusters(0),
88 fFillGlobal(kTRUE)
89{
90 // Default constructor
91
92}
93
94//________________________________________________________________________
95AliAnalysisTaskSED0Correlations::AliAnalysisTaskSED0Correlations(const char *name,AliRDHFCutsD0toKpi* cutsD0, AliHFAssociatedTrackCuts* cutsTrk):
96 AliAnalysisTaskSE(name),
97 fNPtBinsCorr(0),
98 fBinLimsCorr(),
99 fPtThreshLow(),
100 fPtThreshUp(),
101 fEvents(0),
102 fAlreadyFilled(kFALSE),
103 fOutputMass(0),
104 fOutputCorr(0),
105 fOutputStudy(0),
106 fNentries(0),
107 fCutsD0(0),
108 fCutsTracks(cutsTrk),
bce70c96 109 fCorrelatorTr(0),
110 fCorrelatorKc(0),
111 fCorrelatorK0(0),
a2ad7da1 112 fReadMC(0),
bce70c96 113 fMixing(kFALSE),
a2ad7da1 114 fCounter(0),
115 fNPtBins(1),
116 fFillOnlyD0D0bar(0),
117 fIsSelectedCandidate(0),
118 fSys(0),
119 fIsRejectSDDClusters(0),
120 fFillGlobal(kTRUE)
121{
122 // Default constructor
123
124 fNPtBins=cutsD0->GetNPtBins();
125
126 fCutsD0=cutsD0;
127
128 // Output slot #1 writes into a TList container (mass with cuts)
129 DefineOutput(1,TList::Class()); //My private output
130 // Output slot #2 writes into a TH1F container (number of events)
131 DefineOutput(2,TH1F::Class()); //My private output
132 // Output slot #3 writes into a AliRDHFD0toKpi container (cuts)
133 DefineOutput(3,AliRDHFCutsD0toKpi::Class()); //My private output
134 // Output slot #4 writes Normalization Counter
135 DefineOutput(4,AliNormalizationCounter::Class());
136 // Output slot #5 writes into a TList container (correl output)
137 DefineOutput(5,TList::Class()); //My private output
138 // Output slot #6 writes into a TList container (correl advanced)
139 DefineOutput(6,TList::Class()); //My private output
140 // Output slot #7 writes into a AliHFAssociatedTrackCuts container (cuts)
141 DefineOutput(7,AliHFAssociatedTrackCuts::Class()); //My private output
142}
143
144//________________________________________________________________________
145AliAnalysisTaskSED0Correlations::AliAnalysisTaskSED0Correlations(const AliAnalysisTaskSED0Correlations &source):
146 AliAnalysisTaskSE(source),
147 fNPtBinsCorr(source.fNPtBinsCorr),
148 fBinLimsCorr(source.fBinLimsCorr),
149 fPtThreshLow(source.fPtThreshLow),
150 fPtThreshUp(source.fPtThreshUp),
151 fEvents(source.fEvents),
152 fAlreadyFilled(source.fAlreadyFilled),
153 fOutputMass(source.fOutputMass),
154 fOutputCorr(source.fOutputCorr),
155 fOutputStudy(source.fOutputStudy),
156 fNentries(source.fNentries),
157 fCutsD0(source.fCutsD0),
158 fCutsTracks(source.fCutsTracks),
bce70c96 159 fCorrelatorTr(source.fCorrelatorTr),
160 fCorrelatorKc(source.fCorrelatorKc),
161 fCorrelatorK0(source.fCorrelatorK0),
a2ad7da1 162 fReadMC(source.fReadMC),
bce70c96 163 fMixing(source.fMixing),
a2ad7da1 164 fCounter(source.fCounter),
165 fNPtBins(source.fNPtBins),
166 fFillOnlyD0D0bar(source.fFillOnlyD0D0bar),
167 fIsSelectedCandidate(source.fIsSelectedCandidate),
168 fSys(source.fSys),
169 fIsRejectSDDClusters(source.fIsRejectSDDClusters),
170 fFillGlobal(source.fFillGlobal)
171{
172 // Copy constructor
173}
174
175//________________________________________________________________________
176AliAnalysisTaskSED0Correlations::~AliAnalysisTaskSED0Correlations()
177{
178 if (fOutputMass) {
179 delete fOutputMass;
180 fOutputMass = 0;
181 }
182 if (fOutputCorr) {
183 delete fOutputCorr;
184 fOutputCorr = 0;
185 }
186 if (fOutputStudy) {
187 delete fOutputStudy;
188 fOutputStudy = 0;
189 }
190 if (fCutsD0) {
191 delete fCutsD0;
192 fCutsD0 = 0;
193 }
194 if (fNentries){
195 delete fNentries;
196 fNentries = 0;
197 }
bce70c96 198 if (fCorrelatorTr) {
199 delete fCorrelatorTr;
200 fCorrelatorTr = 0;
201 }
202 if (fCorrelatorKc) {
203 delete fCorrelatorKc;
204 fCorrelatorKc = 0;
205 }
206 if (fCorrelatorK0) {
207 delete fCorrelatorK0;
208 fCorrelatorK0 = 0;
209 }
210 if (fCounter){
a2ad7da1 211 delete fCounter;
212 fCounter=0;
213 }
214}
215
216//______________________________________________________________________________
217AliAnalysisTaskSED0Correlations& AliAnalysisTaskSED0Correlations::operator=(const AliAnalysisTaskSED0Correlations& orig)
218{
219// Assignment
220 if (&orig == this) return *this; //if address is the same (same object), returns itself
221
222 AliAnalysisTaskSE::operator=(orig); //Uses the AliAnalysisTaskSE operator to assign the inherited part of the class
223 fNPtBinsCorr = orig.fNPtBinsCorr;
224 fBinLimsCorr = orig.fBinLimsCorr;
225 fPtThreshLow = orig.fPtThreshLow;
226 fPtThreshUp = orig.fPtThreshUp;
227 fEvents = orig.fEvents;
228 fAlreadyFilled = orig.fAlreadyFilled;
229 fOutputMass = orig.fOutputMass;
230 fOutputCorr = orig.fOutputCorr;
231 fOutputStudy = orig.fOutputStudy;
232 fNentries = orig.fNentries;
233 fCutsD0 = orig.fCutsD0;
234 fCutsTracks = orig.fCutsTracks;
bce70c96 235 fCorrelatorTr = orig.fCorrelatorTr;
236 fCorrelatorKc = orig.fCorrelatorKc;
237 fCorrelatorK0 = orig.fCorrelatorK0;
a2ad7da1 238 fReadMC = orig.fReadMC;
bce70c96 239 fMixing = orig.fMixing;
a2ad7da1 240 fCounter = orig.fCounter;
241 fNPtBins = orig.fNPtBins;
242 fFillOnlyD0D0bar = orig.fFillOnlyD0D0bar;
243 fIsSelectedCandidate = orig.fIsSelectedCandidate;
244 fSys = orig.fSys;
245 fIsRejectSDDClusters = orig.fIsRejectSDDClusters;
246 fFillGlobal = orig.fFillGlobal;
247
248 return *this; //returns pointer of the class
249}
250
251//________________________________________________________________________
252void AliAnalysisTaskSED0Correlations::Init()
253{
254 // Initialization
255
256 if(fDebug > 1) printf("AnalysisTaskSED0Correlations::Init() \n");
257
258 //Copy of cuts objects
259 AliRDHFCutsD0toKpi* copyfCutsD0 = new AliRDHFCutsD0toKpi(*fCutsD0);
260 const char* nameoutput=GetOutputSlot(3)->GetContainer()->GetName();
261 copyfCutsD0->SetName(nameoutput);
262
bce70c96 263 //needer to clear completely the objects inside with Clear() method
a2ad7da1 264 // Post the data
265 PostData(3,copyfCutsD0);
266 PostData(7,fCutsTracks);
267
268 return;
269}
270
271//________________________________________________________________________
272void AliAnalysisTaskSED0Correlations::UserCreateOutputObjects()
273{
274
275 // Create the output container
276 //
277 if(fDebug > 1) printf("AnalysisTaskSED0Correlations::UserCreateOutputObjects() \n");
278
bce70c96 279 //HFCorrelator creation and definition
280 fCorrelatorTr = new AliHFCorrelator("CorrelatorTr",fCutsTracks,fSys);
281 fCorrelatorKc = new AliHFCorrelator("CorrelatorKc",fCutsTracks,fSys);
282 fCorrelatorK0 = new AliHFCorrelator("CorrelatorK0",fCutsTracks,fSys);
283 fCorrelatorTr->SetDeltaPhiInterval(-1.57,4.71);// set the Delta Phi Interval you want (in this case -0.5Pi to 1.5 Pi)
284 fCorrelatorKc->SetDeltaPhiInterval(-1.57,4.71);
285 fCorrelatorK0->SetDeltaPhiInterval(-1.57,4.71);
286 fCorrelatorTr->SetEventMixing(fMixing);// sets the analysis on a single event (kFALSE) or mixed events (kTRUE)
287 fCorrelatorKc->SetEventMixing(fMixing);
288 fCorrelatorK0->SetEventMixing(fMixing);
289 fCorrelatorTr->SetAssociatedParticleType(1);// set 1 for correlations with hadrons, 2 with kaons, 3 with KZeros
290 fCorrelatorKc->SetAssociatedParticleType(2);// set 1 for correlations with hadrons, 2 with kaons, 3 with KZeros
291 fCorrelatorK0->SetAssociatedParticleType(3);// set 1 for correlations with hadrons, 2 with kaons, 3 with KZeros
292 fCorrelatorTr->SetApplyDisplacementCut(2); //0: don't calculate d0; 1: return d0; 2: return d0/d0err
293 fCorrelatorKc->SetApplyDisplacementCut(2);
294 fCorrelatorK0->SetApplyDisplacementCut(0);
295 fCorrelatorTr->SetUseMC(fReadMC);// sets Montecarlo flag
296 fCorrelatorKc->SetUseMC(fReadMC);
297 fCorrelatorK0->SetUseMC(fReadMC);
298 fCorrelatorKc->SetPIDmode(2); //switch for K+/- PID option
299 Bool_t pooldefTr = fCorrelatorTr->DefineEventPool();// method that defines the properties ot the event mixing (zVtx and Multipl. bins)
300 Bool_t pooldefKc = fCorrelatorKc->DefineEventPool();// method that defines the properties ot the event mixing (zVtx and Multipl. bins)
301 Bool_t pooldefK0 = fCorrelatorK0->DefineEventPool();// method that defines the properties ot the event mixing (zVtx and Multipl. bins)
302 if(!pooldefTr) AliInfo("Warning:: Event pool not defined properly");
303 if(!pooldefKc) AliInfo("Warning:: Event pool not defined properly");
304 if(!pooldefK0) AliInfo("Warning:: Event pool not defined properly");
305
a2ad7da1 306 // Several histograms are more conveniently managed in a TList
307 fOutputMass = new TList();
308 fOutputMass->SetOwner();
309 fOutputMass->SetName("listMass");
310
311 fOutputCorr = new TList();
312 fOutputCorr->SetOwner();
313 fOutputCorr->SetName("correlationslist");
314
315 fOutputStudy = new TList();
316 fOutputStudy->SetOwner();
317 fOutputStudy->SetName("MCstudyplots");
318
319 TString nameMass=" ",nameSgn=" ", nameBkg=" ", nameRfl=" ";
320
321 for(Int_t i=0;i<fCutsD0->GetNPtBins();i++){
322
323 nameMass="histMass_";
324 nameMass+=i;
325 nameSgn="histSgn_";
326 nameSgn+=i;
327 nameBkg="histBkg_";
328 nameBkg+=i;
329 nameRfl="histRfl_";
330 nameRfl+=i;
331
332 //histograms of invariant mass distributions
333
334 //MC signal
335 if(fReadMC){
336 TH1F* tmpSt = new TH1F(nameSgn.Data(), "D^{0} invariant mass - MC; M [GeV]; Entries",120,1.5648,2.1648);
337 tmpSt->Sumw2();
338
339 //Reflection: histo filled with D0Mass which pass the cut (also) as D0bar and with D0bar which pass (also) the cut as D0
340 TH1F* tmpRt = new TH1F(nameRfl.Data(), "Reflected signal invariant mass - MC; M [GeV]; Entries",120,1.5648,2.1648);
341 TH1F* tmpBt = new TH1F(nameBkg.Data(), "Background invariant mass - MC; M [GeV]; Entries",120,1.5648,2.1648);
342 tmpBt->Sumw2();
343 tmpRt->Sumw2();
344 fOutputMass->Add(tmpSt);
345 fOutputMass->Add(tmpRt);
346 fOutputMass->Add(tmpBt);
347 }
348
349 //mass
350 TH1F* tmpMt = new TH1F(nameMass.Data(),"D^{0} invariant mass; M [GeV]; Entries",120,1.5648,2.1648);
351 tmpMt->Sumw2();
352 fOutputMass->Add(tmpMt);
353 }
354
355 const char* nameoutput=GetOutputSlot(2)->GetContainer()->GetName();
356
357 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", 18,-0.5,17.5);
358
359 fNentries->GetXaxis()->SetBinLabel(1,"nEventsAnal");
360 fNentries->GetXaxis()->SetBinLabel(2,"nCandSel(Cuts)");
361 fReadMC ? fNentries->GetXaxis()->SetBinLabel(3,"nD0Selected") : fNentries->GetXaxis()->SetBinLabel(3,"Dstar<-D0");
362 fNentries->GetXaxis()->SetBinLabel(4,"nEventsGoodVtxS");
363 fNentries->GetXaxis()->SetBinLabel(5,"ptbin = -1");
364 fNentries->GetXaxis()->SetBinLabel(6,"no daughter");
365 if(fSys==0) fNentries->GetXaxis()->SetBinLabel(7,"nCandSel(Tr)");
366 if(fReadMC && fSys==0){
367 fNentries->GetXaxis()->SetBinLabel(12,"K");
368 fNentries->GetXaxis()->SetBinLabel(13,"Lambda");
369 }
370 fNentries->GetXaxis()->SetBinLabel(14,"Pile-up Rej");
371 fNentries->GetXaxis()->SetBinLabel(15,"N. of 0SMH");
372 if(fSys==1) fNentries->GetXaxis()->SetBinLabel(16,"Nev in centr");
373 if(fIsRejectSDDClusters) fNentries->GetXaxis()->SetBinLabel(17,"SDD-Cls Rej");
374 fNentries->GetXaxis()->SetBinLabel(18,"Phys.Sel.Rej");
375 fNentries->GetXaxis()->SetNdivisions(1,kFALSE);
376
377 fCounter = new AliNormalizationCounter(Form("%s",GetOutputSlot(4)->GetContainer()->GetName()));
378 fCounter->Init();
379
380 CreateCorrelationsObjs(); //creates histos for correlations analysis
381
382 // Post the data
383 PostData(1,fOutputMass);
384 PostData(2,fNentries);
385 PostData(4,fCounter);
386 PostData(5,fOutputCorr);
387 PostData(6,fOutputStudy);
388
389 return;
390}
391
392//________________________________________________________________________
393void AliAnalysisTaskSED0Correlations::UserExec(Option_t */*option*/)
394{
395 // Execute analysis for current event:
396 // heavy flavor candidates association to MC truth
397 //cout<<"I'm in UserExec"<<endl;
398
399
400 //cuts order
401 // printf(" |M-MD0| [GeV] < %f\n",fD0toKpiCuts[0]);
402 // printf(" dca [cm] < %f\n",fD0toKpiCuts[1]);
403 // printf(" cosThetaStar < %f\n",fD0toKpiCuts[2]);
404 // printf(" pTK [GeV/c] > %f\n",fD0toKpiCuts[3]);
405 // printf(" pTpi [GeV/c] > %f\n",fD0toKpiCuts[4]);
406 // printf(" |d0K| [cm] < %f\n",fD0toKpiCuts[5]);
407 // printf(" |d0pi| [cm] < %f\n",fD0toKpiCuts[6]);
408 // printf(" d0d0 [cm^2] < %f\n",fD0toKpiCuts[7]);
409 // printf(" cosThetaPoint > %f\n",fD0toKpiCuts[8]);
410
411
412 AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
413 fEvents++;
414
415 TString bname="D0toKpi";
416
417 TClonesArray *inputArray=0;
418
419 if(!aod && AODEvent() && IsStandardAOD()) {
420 // In case there is an AOD handler writing a standard AOD, use the AOD
421 // event in memory rather than the input (ESD) event.
422 aod = dynamic_cast<AliAODEvent*> (AODEvent());
423 // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
424 // have to taken from the AOD event hold by the AliAODExtension
425 AliAODHandler* aodHandler = (AliAODHandler*)
426 ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
427
428 if(aodHandler->GetExtensions()) {
429 AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
430 AliAODEvent* aodFromExt = ext->GetAOD();
431 inputArray=(TClonesArray*)aodFromExt->GetList()->FindObject(bname.Data());
432 }
433 } else if(aod) {
434 inputArray=(TClonesArray*)aod->GetList()->FindObject(bname.Data());
435 }
436
437 if(!inputArray || !aod) {
438 printf("AliAnalysisTaskSED0Correlations::UserExec: input branch not found!\n");
439 return;
440 }
441
442 // fix for temporary bug in ESDfilter
443 // the AODs with null vertex pointer didn't pass the PhysSel
444 if(!aod->GetPrimaryVertex() || TMath::Abs(aod->GetMagneticField())<0.001) return;
445
446 TClonesArray *mcArray = 0;
447 AliAODMCHeader *mcHeader = 0;
448
449 if(fReadMC) {
450 // load MC particles
451 mcArray = (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName());
452 if(!mcArray) {
453 printf("AliAnalysisTaskSED0Correlations::UserExec: MC particles branch not found!\n");
454 return;
455 }
456
457 // load MC header
458 mcHeader = (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
459 if(!mcHeader) {
460 printf("AliAnalysisTaskSED0Correlations::UserExec: MC header branch not found!\n");
461 return;
462 }
463 }
bce70c96 464
a2ad7da1 465 //histogram filled with 1 for every AOD
466 fNentries->Fill(0);
467 fCounter->StoreEvent(aod,fCutsD0,fReadMC);
468
469 // trigger class for PbPb C0SMH-B-NOPF-ALLNOTRD, C0SMH-B-NOPF-ALL
470 TString trigclass=aod->GetFiredTriggerClasses();
471 if(trigclass.Contains("C0SMH-B-NOPF-ALLNOTRD") || trigclass.Contains("C0SMH-B-NOPF-ALL")) fNentries->Fill(14);
472
473 if(!fCutsD0->IsEventSelected(aod)) {
474 if(fCutsD0->GetWhyRejection()==1) // rejected for pileup
475 fNentries->Fill(13);
476 if(fSys==1 && (fCutsD0->GetWhyRejection()==2 || fCutsD0->GetWhyRejection()==3)) fNentries->Fill(15);
477 if(fCutsD0->GetWhyRejection()==7) fNentries->Fill(17);
478 return;
479 }
480
481 // Check the Nb of SDD clusters
482 if (fIsRejectSDDClusters) {
483 Bool_t skipEvent = kFALSE;
484 Int_t ntracks = 0;
485 if (aod) ntracks = aod->GetNTracks();
486 for(Int_t itrack=0; itrack<ntracks; itrack++) { // loop on tacks
487 // ... get the track
488 AliAODTrack * track = aod->GetTrack(itrack);
489 if(TESTBIT(track->GetITSClusterMap(),2) || TESTBIT(track->GetITSClusterMap(),3) ){
490 skipEvent=kTRUE;
491 fNentries->Fill(16);
492 break;
493 }
494 }
495 if (skipEvent) return;
496 }
497
bce70c96 498 //HFCorrelators initialization (for this event)
499 fCorrelatorTr->SetAODEvent(aod); // set the AOD event from which you are processing
500 fCorrelatorKc->SetAODEvent(aod);
501 fCorrelatorK0->SetAODEvent(aod);
502 Bool_t correlatorONTr = fCorrelatorTr->Initialize(); // initialize the pool for event mixing
503 Bool_t correlatorONKc = fCorrelatorKc->Initialize();
504 Bool_t correlatorONK0 = fCorrelatorK0->Initialize();
505 if(!correlatorONTr) {AliInfo("AliHFCorrelator (tracks) didn't initialize the pool correctly or processed a bad event"); return;}
506 if(!correlatorONKc) {AliInfo("AliHFCorrelator (charged K) didn't initialize the pool correctly or processed a bad event"); return;}
507 if(!correlatorONK0) {AliInfo("AliHFCorrelator (K0) didn't initialize the pool correctly or processed a bad event"); return;}
508
509 if(fReadMC) {
510 fCorrelatorTr->SetMCArray(mcArray); // set the TClonesArray *fmcArray for analysis on monte carlo
511 fCorrelatorKc->SetMCArray(mcArray);
512 fCorrelatorK0->SetMCArray(mcArray);
513 }
514
a2ad7da1 515 // AOD primary vertex
516 AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex();
a2ad7da1 517 Bool_t isGoodVtx=kFALSE;
518
519 //vtx1->Print();
520 TString primTitle = vtx1->GetTitle();
521 if(primTitle.Contains("VertexerTracks") && vtx1->GetNContributors()>0) {
522 isGoodVtx=kTRUE;
523 fNentries->Fill(3);
524 }
525
526 //Reset flag for tracks distributions fill
527 fAlreadyFilled=kFALSE;
528
bce70c96 529 //***** Loop over D0 candidates *****
a2ad7da1 530 Int_t nInD0toKpi = inputArray->GetEntriesFast();
531 if(fDebug>2) printf("Number of D0->Kpi: %d\n",nInD0toKpi);
532
533 if(fFillGlobal) { //loop on V0 and tracks for each event, to fill Pt distr. and InvMass distr.
534
535 TClonesArray *v0array = (TClonesArray*)aod->GetList()->FindObject("v0s");
536 Int_t pdgCodes[2] = {211,211};
537 Int_t idArrayV0[v0array->GetEntriesFast()][2];
538 for(int iV0=0; iV0<v0array->GetEntriesFast(); iV0++) {
539 for(int j=0; j<2; j++) {idArrayV0[iV0][j]=-2;}
540 AliAODv0 *v0 = (AliAODv0*)v0array->UncheckedAt(iV0);
541 if(SelectV0(v0,vtx1,2,idArrayV0)) { //option 2 = for mass inv plots only
bce70c96 542 if(fReadMC && (v0->MatchToMC(310,mcArray,2,pdgCodes)<0)) continue; //310 = K0s, 311 = K0 generico!!
a2ad7da1 543 ((TH2F*)fOutputStudy->FindObject("hK0MassInv"))->Fill(v0->MassK0Short(),v0->Pt()); //invariant mass plot
f80e7bba 544 ((TH1F*)fOutputStudy->FindObject("hist_Pt_K0_AllEv"))->Fill(v0->Pt()); //pT distribution (in all events), K0 case
a2ad7da1 545 }
546 }
547
548 for(Int_t itrack=0; itrack<aod->GetNTracks(); itrack++) { // loop on tacks
549 AliAODTrack * track = aod->GetTrack(itrack);
550 //rejection of tracks
551 if(track->GetID() < 0) continue; //discard negative ID tracks
552 if(track->Pt() < fPtThreshLow.at(0) || track->Pt() > fPtThreshUp.at(0)) continue; //discard tracks outside pt range for hadrons/K
bce70c96 553 if(!fCutsTracks->IsHadronSelected(track) || !fCutsTracks->CheckHadronKinematic(track->Pt(),0.1)) continue; //0.1 = dummy (d0 value, no cut on it for me)
a2ad7da1 554 //pT distribution (in all events), charg and hadr cases
555 ((TH1F*)fOutputStudy->FindObject("hist_Pt_Charg_AllEv"))->Fill(track->Pt());
556 if(fCutsTracks->CheckKaonCompatibility(track,kFALSE,0,2)) ((TH1F*)fOutputStudy->FindObject("hist_Pt_Kcharg_AllEv"))->Fill(track->Pt());
557 }
558
559 } //end of loops for global plot fill
560
561 Int_t nSelectedloose=0,nSelectedtight=0;
562 for (Int_t iD0toKpi = 0; iD0toKpi < nInD0toKpi; iD0toKpi++) {
563 AliAODRecoDecayHF2Prong *d = (AliAODRecoDecayHF2Prong*)inputArray->UncheckedAt(iD0toKpi);
564
565 if(d->GetSelectionMap()) if(!d->HasSelectionBit(AliRDHFCuts::kD0toKpiCuts)){
566 fNentries->Fill(2);
567 continue; //skip the D0 from Dstar
568 }
569
570 if (fCutsD0->IsInFiducialAcceptance(d->Pt(),d->Y(421)) ) {
571 nSelectedloose++;
572 nSelectedtight++;
573 if(fSys==0){
574 if(fCutsD0->IsSelected(d,AliRDHFCuts::kTracks,aod))fNentries->Fill(6);
575 }
576 Int_t ptbin=fCutsD0->PtBin(d->Pt());
577 if(ptbin==-1) {fNentries->Fill(4); continue;} //out of bounds
578
bce70c96 579 fIsSelectedCandidate=fCutsD0->IsSelected(d,AliRDHFCuts::kAll,aod); //D0 selected
a2ad7da1 580 if(!fIsSelectedCandidate) continue;
581
bce70c96 582 //D0 infos
583 Double_t phiD0 = fCorrelatorTr->SetCorrectPhiRange(d->Phi());
584 phiD0 = fCorrelatorKc->SetCorrectPhiRange(d->Phi()); //bad usage, but returns a Double_t...
585 phiD0 = fCorrelatorK0->SetCorrectPhiRange(d->Phi());
586 fCorrelatorTr->SetTriggerParticleProperties(d->Pt(),phiD0,d->Eta()); // sets the parameters of the trigger particles that are needed
587 fCorrelatorKc->SetTriggerParticleProperties(d->Pt(),phiD0,d->Eta());
588 fCorrelatorK0->SetTriggerParticleProperties(d->Pt(),phiD0,d->Eta());
589 fCorrelatorTr->SetD0Properties(d,fIsSelectedCandidate); //sets special properties for D0
590 fCorrelatorKc->SetD0Properties(d,fIsSelectedCandidate);
591 fCorrelatorK0->SetD0Properties(d,fIsSelectedCandidate);
592
593 if(!fReadMC) {
594 CalculateCorrelations(d); //correlations on real data
595 } else { //correlations on MC -> association of selected D0 to MCinfo with MCtruth
a2ad7da1 596 Int_t pdgDgD0toKpi[2]={321,211};
597 Int_t labD0 = d->MatchToMC(421,mcArray,2,pdgDgD0toKpi); //return MC particle label if the array corresponds to a D0, -1 if not
bce70c96 598 if (labD0>-1) {
599 CalculateCorrelations(d,labD0,mcArray);
600 }
a2ad7da1 601 }
602
603 FillMassHists(d,mcArray,fCutsD0,fOutputMass);
604 }
605
606 } //end for prongs
bce70c96 607
608 if(fMixing /* && fAlreadyFilled*/) { // update the pool for Event Mixing, if: enabled, event is ok, at least a SelD0 found! (fAlreadyFilled's role!)
609 Bool_t updatedTr = fCorrelatorTr->PoolUpdate();
610 Bool_t updatedKc = fCorrelatorKc->PoolUpdate();
611 Bool_t updatedK0 = fCorrelatorK0->PoolUpdate();
612 if(!updatedTr || !updatedKc || !updatedK0) AliInfo("Pool was not updated");
613 }
614
a2ad7da1 615 fCounter->StoreCandidates(aod,nSelectedloose,kTRUE);
616 fCounter->StoreCandidates(aod,nSelectedtight,kFALSE);
617
618 // Post the data
619 PostData(1,fOutputMass);
620 PostData(2,fNentries);
621 PostData(4,fCounter);
622 PostData(5,fOutputCorr);
623 PostData(6,fOutputStudy);
624
625 return;
626}
627
628//____________________________________________________________________________
629void AliAnalysisTaskSED0Correlations::FillMassHists(AliAODRecoDecayHF2Prong *part, TClonesArray *arrMC, AliRDHFCutsD0toKpi* cuts, TList *listout){
630 //
631 // function used in UserExec to fill mass histograms:
632 //
633 if (!fIsSelectedCandidate) return;
634
635 if(fDebug>2) cout<<"Candidate selected"<<endl;
636
637 Double_t invmassD0 = part->InvMassD0(), invmassD0bar = part->InvMassD0bar();
638 Int_t ptbin = cuts->PtBin(part->Pt());
639
640 TString fillthis="";
641 Int_t pdgDgD0toKpi[2]={321,211};
642 Int_t labD0=-1;
643 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)
644
645 //count candidates selected by cuts
646 fNentries->Fill(1);
647 //count true D0 selected by cuts
648 if (fReadMC && labD0>=0) fNentries->Fill(2);
649
650 if ((fIsSelectedCandidate==1 || fIsSelectedCandidate==3) && fFillOnlyD0D0bar<2) { //D0
651
652 if(fReadMC){
653 if(labD0>=0) {
654 AliAODMCParticle *partD0 = (AliAODMCParticle*)arrMC->At(labD0);
655 Int_t pdgD0 = partD0->GetPdgCode();
656 if (pdgD0==421){ //D0
657 fillthis="histSgn_";
658 fillthis+=ptbin;
659 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0);
660 } else{ //it was a D0bar
661 fillthis="histRfl_";
662 fillthis+=ptbin;
663 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0);
664 }
665 } else {//background
666 fillthis="histBkg_";
667 fillthis+=ptbin;
668 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0);
669 }
670 }else{
671 fillthis="histMass_";
672 fillthis+=ptbin;
673 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0);
674 }
675
676 }
677 if (fIsSelectedCandidate>1 && (fFillOnlyD0D0bar==0 || fFillOnlyD0D0bar==2)) { //D0bar
678
679 if(fReadMC){
680 if(labD0>=0) {
681 AliAODMCParticle *partD0 = (AliAODMCParticle*)arrMC->At(labD0);
682 Int_t pdgD0 = partD0->GetPdgCode();
683
684 if (pdgD0==-421){ //D0bar
685 fillthis="histSgn_";
686 fillthis+=ptbin;
687 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar);
688 } else{
689 fillthis="histRfl_";
690 fillthis+=ptbin;
691 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar);
692 }
693 } else {//background or LS
694 fillthis="histBkg_";
695 fillthis+=ptbin;
696 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar);
697 }
698 }else{
699 fillthis="histMass_";
700 fillthis+=ptbin;
701 ((TH1F*)listout->FindObject(fillthis))->Fill(invmassD0bar);
702 }
703 }
704
705 return;
706}
707
708//________________________________________________________________________
709void AliAnalysisTaskSED0Correlations::Terminate(Option_t */*option*/)
710{
711 // Terminate analysis
712 //
713 if(fDebug > 1) printf("AnalysisTaskSED0Correlations: Terminate() \n");
714
715 fOutputMass = dynamic_cast<TList*> (GetOutputData(1));
716 if (!fOutputMass) {
717 printf("ERROR: fOutputMass not available\n");
718 return;
719 }
720
721 fNentries = dynamic_cast<TH1F*>(GetOutputData(2));
722
723 if(!fNentries){
724 printf("ERROR: fNEntries not available\n");
725 return;
726 }
727
728 fCutsD0 = dynamic_cast<AliRDHFCutsD0toKpi*>(GetOutputData(3));
729 if(!fCutsD0){
730 printf("ERROR: fCuts not available\n");
731 return;
732 }
733
734 fCounter = dynamic_cast<AliNormalizationCounter*>(GetOutputData(4));
735 if (!fCounter) {
736 printf("ERROR: fCounter not available\n");
737 return;
738 }
739 fOutputCorr = dynamic_cast<TList*> (GetOutputData(5));
740 if (!fOutputCorr) {
741 printf("ERROR: fOutputCorr not available\n");
742 return;
743 }
744 fOutputStudy = dynamic_cast<TList*> (GetOutputData(6));
745 if (!fOutputStudy) {
746 printf("ERROR: fOutputStudy not available\n");
747 return;
748 }
749 fCutsTracks = dynamic_cast<AliHFAssociatedTrackCuts*>(GetOutputData(7));
750 if(!fCutsTracks){
751 printf("ERROR: fCutsTracks not available\n");
752 return;
753 }
754
755 return;
756}
757
758//_________________________________________________________________________________________________
759Int_t AliAnalysisTaskSED0Correlations::CheckD0Origin(TClonesArray* arrayMC, AliAODMCParticle *mcPartCandidate) const {
760 //
761 // checking whether the mother of the particles come from a charm or a bottom quark
762 //
bce70c96 763 printf("AliAnalysisTaskSED0Correlations::CheckD0Origin() \n");
a2ad7da1 764
765 Int_t pdgGranma = 0;
766 Int_t mother = 0;
767 mother = mcPartCandidate->GetMother();
768 Int_t abspdgGranma =0;
769 Bool_t isFromB=kFALSE;
770 Bool_t isQuarkFound=kFALSE;
771
772 while (mother > 0){
773 AliAODMCParticle* mcGranma = dynamic_cast<AliAODMCParticle*>(arrayMC->At(mother));
774 if (mcGranma){
775 pdgGranma = mcGranma->GetPdgCode();
776 abspdgGranma = TMath::Abs(pdgGranma);
777 if ((abspdgGranma > 500 && abspdgGranma < 600) || (abspdgGranma > 5000 && abspdgGranma < 6000)){
778 isFromB=kTRUE;
779 }
780 if(abspdgGranma==4 || abspdgGranma==5) isQuarkFound=kTRUE;
781 mother = mcGranma->GetMother();
782 }else{
783 AliError("Failed casting the mother particle!");
784 break;
785 }
786 }
787
788 if(isQuarkFound) {
789 if(isFromB) return 5;
790 else return 4;
791 }
792 else return 1;
793}
794
795//________________________________________________________________________
796void AliAnalysisTaskSED0Correlations::CreateCorrelationsObjs() {
797//
798
799 TString namePlot = "";
800
801 //These for limits in THnSparse (one per bin, same limits).
bce70c96 802 //Vars: DeltaPhi, InvMass, PtTrack, Displacement, DeltaEta
f80e7bba 803 Int_t nBinsPhi[5] = {32,150,30,3,16};
804 Double_t binMinPhi[5] = {-1.6,1.6,0.,0.,-1.6}; //is the minimum for all the bins
805 Double_t binMaxPhi[5] = {4.8,2.2,3.0,3.,1.6}; //is the maximum for all the bins
a2ad7da1 806
bce70c96 807 //Vars: DeltaPhi, InvMass, DeltaEta
808 Int_t nBinsMix[3] = {32,150,16};
809 Double_t binMinMix[3] = {-1.6,1.6,-1.6}; //is the minimum for all the bins
810 Double_t binMaxMix[3] = {4.8,2.2,1.6}; //is the maximum for all the bins
a2ad7da1 811
bce70c96 812 for(Int_t i=0;i<fNPtBinsCorr;i++){
a2ad7da1 813
bce70c96 814 if(!fMixing) {
815 //THnSparse plots: correlations for various invariant mass (MC and data)
816 namePlot="hPhi_K0_Bin";
a2ad7da1 817 namePlot+=i;
818
bce70c96 819 THnSparseI *hPhiK = new THnSparseI(namePlot.Data(), "Azimuthal correlation; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
820 hPhiK->Sumw2();
821 fOutputCorr->Add(hPhiK);
a2ad7da1 822
bce70c96 823 namePlot="hPhi_Kcharg_Bin";
a2ad7da1 824 namePlot+=i;
825
bce70c96 826 THnSparseI *hPhiH = new THnSparseI(namePlot.Data(), "Azimuthal correlation; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
827 hPhiH->Sumw2();
828 fOutputCorr->Add(hPhiH);
a2ad7da1 829
bce70c96 830 namePlot="hPhi_Charg_Bin";
a2ad7da1 831 namePlot+=i;
832
bce70c96 833 THnSparseI *hPhiC = new THnSparseI(namePlot.Data(), "Azimuthal correlation; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
834 hPhiC->Sumw2();
835 fOutputCorr->Add(hPhiC);
a2ad7da1 836
bce70c96 837 //histos for c/b origin for D0 (MC only)
838 if (fReadMC) {
a2ad7da1 839
bce70c96 840 //generic origin for tracks
841 namePlot="hPhi_K0_From_c_Bin";
842 namePlot+=i;
a2ad7da1 843
bce70c96 844 THnSparseI *hPhiK_c = new THnSparseI(namePlot.Data(), "Azimuthal correlation - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
845 hPhiK_c->Sumw2();
846 fOutputCorr->Add(hPhiK_c);
a2ad7da1 847
bce70c96 848 namePlot="hPhi_Kcharg_From_c_Bin";
849 namePlot+=i;
a2ad7da1 850
bce70c96 851 THnSparseI *hPhiH_c = new THnSparseI(namePlot.Data(), "Azimuthal correlation - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
852 hPhiH_c->Sumw2();
853 fOutputCorr->Add(hPhiH_c);
a2ad7da1 854
bce70c96 855 namePlot="hPhi_Charg_From_c_Bin";
856 namePlot+=i;
a2ad7da1 857
bce70c96 858 THnSparseI *hPhiC_c = new THnSparseI(namePlot.Data(), "Azimuthal correlation - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
859 hPhiC_c->Sumw2();
860 fOutputCorr->Add(hPhiC_c);
861
862 namePlot="hPhi_K0_From_b_Bin";
863 namePlot+=i;
a2ad7da1 864
bce70c96 865 THnSparseI *hPhiK_b = new THnSparseI(namePlot.Data(), "Azimuthal correlation - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
866 hPhiK_b->Sumw2();
867 fOutputCorr->Add(hPhiK_b);
a2ad7da1 868
bce70c96 869 namePlot="hPhi_Kcharg_From_b_Bin";
870 namePlot+=i;
a2ad7da1 871
bce70c96 872 THnSparseI *hPhiH_b = new THnSparseI(namePlot.Data(), "Azimuthal correlation - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
873 hPhiH_b->Sumw2();
874 fOutputCorr->Add(hPhiH_b);
a2ad7da1 875
bce70c96 876 namePlot="hPhi_Charg_From_b_Bin";
877 namePlot+=i;
a2ad7da1 878
bce70c96 879 THnSparseI *hPhiC_b = new THnSparseI(namePlot.Data(), "Azimuthal correlation - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
880 hPhiC_b->Sumw2();
881 fOutputCorr->Add(hPhiC_b);
a2ad7da1 882
bce70c96 883 //HF-only tracks (c for c->D0, b for b->D0)
884 namePlot="hPhi_K0_HF_From_c_Bin";
885 namePlot+=i;
a2ad7da1 886
bce70c96 887 THnSparseI *hPhiK_HF_c = new THnSparseI(namePlot.Data(), "Azimuthal correlation HF - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
888 hPhiK_HF_c->Sumw2();
889 fOutputCorr->Add(hPhiK_HF_c);
a2ad7da1 890
bce70c96 891 namePlot="hPhi_Kcharg_HF_From_c_Bin";
892 namePlot+=i;
a2ad7da1 893
bce70c96 894 THnSparseI *hPhiH_HF_c = new THnSparseI(namePlot.Data(), "Azimuthal correlation HF - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
895 hPhiH_HF_c->Sumw2();
896 fOutputCorr->Add(hPhiH_HF_c);
a2ad7da1 897
bce70c96 898 namePlot="hPhi_Charg_HF_From_c_Bin";
899 namePlot+=i;
a2ad7da1 900
bce70c96 901 THnSparseI *hPhiC_HF_c = new THnSparseI(namePlot.Data(), "Azimuthal correlation HF - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
902 hPhiC_HF_c->Sumw2();
903 fOutputCorr->Add(hPhiC_HF_c);
a2ad7da1 904
bce70c96 905 namePlot="hPhi_K0_HF_From_b_Bin";
906 namePlot+=i;
a2ad7da1 907
bce70c96 908 THnSparseI *hPhiK_HF_b = new THnSparseI(namePlot.Data(), "Azimuthal correlation HF - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
909 hPhiK_HF_b->Sumw2();
910 fOutputCorr->Add(hPhiK_HF_b);
a2ad7da1 911
bce70c96 912 namePlot="hPhi_Kcharg_HF_From_b_Bin";
913 namePlot+=i;
a2ad7da1 914
bce70c96 915 THnSparseI *hPhiH_HF_b = new THnSparseI(namePlot.Data(), "Azimuthal correlation HF - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
916 hPhiH_HF_b->Sumw2();
917 fOutputCorr->Add(hPhiH_HF_b);
a2ad7da1 918
bce70c96 919 namePlot="hPhi_Charg_HF_From_b_Bin";
920 namePlot+=i;
a2ad7da1 921
bce70c96 922 THnSparseI *hPhiC_HF_b = new THnSparseI(namePlot.Data(), "Azimuthal correlation HF - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
923 hPhiC_HF_b->Sumw2();
924 fOutputCorr->Add(hPhiC_HF_b);
925 }
a2ad7da1 926
bce70c96 927 //leading hadron correlations
928 namePlot="hPhi_Lead_Bin";
a2ad7da1 929 namePlot+=i;
930
bce70c96 931 TH2F *hCorrLead = new TH2F(namePlot.Data(), "Leading particle correlation; #Delta#phi",32,-1.6,4.8,300,1.6,2.2);
932 hCorrLead->Sumw2();
933 fOutputCorr->Add(hCorrLead);
a2ad7da1 934
bce70c96 935 if (fReadMC) {
936 namePlot="hPhi_Lead_From_c_Bin";
937 namePlot+=i;
a2ad7da1 938
bce70c96 939 TH2F *hCorrLead_c = new TH2F(namePlot.Data(), "Leading particle correlation - c origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2);
940 hCorrLead_c->Sumw2();
941 fOutputCorr->Add(hCorrLead_c);
942
943 namePlot="hPhi_Lead_From_b_Bin";
944 namePlot+=i;
945
946 TH2F *hCorrLead_b = new TH2F(namePlot.Data(), "Leading particle correlation - b origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2);
947 hCorrLead_b->Sumw2();
948 fOutputCorr->Add(hCorrLead_b);
949
950 namePlot="hPhi_Lead_HF_From_c_Bin";
951 namePlot+=i;
952
953 TH2F *hCorrLead_HF_c = new TH2F(namePlot.Data(), "Leading particle correlation HF - c origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2);
954 hCorrLead_HF_c->Sumw2();
955 fOutputCorr->Add(hCorrLead_HF_c);
956
957 namePlot="hPhi_Lead_HF_From_b_Bin";
958 namePlot+=i;
959
960 TH2F *hCorrLead_HF_b = new TH2F(namePlot.Data(), "Leading particle correlation HF - b origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2);
961 hCorrLead_HF_b->Sumw2();
962 fOutputCorr->Add(hCorrLead_HF_b);
963 }
964
965 //pT weighted correlations
966 namePlot="hPhi_Weig_Bin";
a2ad7da1 967 namePlot+=i;
bce70c96 968
969 TH2F *hCorrWeig = new TH2F(namePlot.Data(), "Charged particle correlation (pT weighted); #Delta#phi",32,-1.6,4.8,300,1.6,2.2);
970 fOutputCorr->Add(hCorrWeig);
971
972 if (fReadMC) {
973 namePlot="hPhi_Weig_From_c_Bin";
974 namePlot+=i;
975
976 TH2F *hCorrWeig_c = new TH2F(namePlot.Data(), "Charged particle correlation (pT weighted) - c origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2);
977 fOutputCorr->Add(hCorrWeig_c);
978
979 namePlot="hPhi_Weig_From_b_Bin";
980 namePlot+=i;
981
982 TH2F *hCorrWeig_b = new TH2F(namePlot.Data(), "Charged particle correlation (pT weighted) - b origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2);
983 fOutputCorr->Add(hCorrWeig_b);
984
985 namePlot="hPhi_Weig_HF_From_c_Bin";
986 namePlot+=i;
987
988 TH2F *hCorrWeig_HF_c = new TH2F(namePlot.Data(), "Charged particle correlation (pT weighted) HF - c origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2);
989 fOutputCorr->Add(hCorrWeig_HF_c);
990
991 namePlot="hPhi_Weig_HF_From_b_Bin";
992 namePlot+=i;
993
994 TH2F *hCorrWeig_HF_b = new TH2F(namePlot.Data(), "Charged particle correlation (pT weighted) HF - b origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2);
995 fOutputCorr->Add(hCorrWeig_HF_b);
996 }
a2ad7da1 997
998 //pT distribution histos
999 namePlot = "hist_Pt_Charg_Bin"; namePlot+=i;
1000 TH1F *hPtC = new TH1F(namePlot.Data(), "Charged track pT (in D0 evs); p_{T} (GeV/c)",240,0.,12.);
1001 hPtC->SetMinimum(0);
1002 fOutputStudy->Add(hPtC);
1003
1004 namePlot = "hist_Pt_Kcharg_Bin"; namePlot+=i;
1005 TH1F *hPtH = new TH1F(namePlot.Data(), "Hadrons pT (in D0 evs); p_{T} (GeV/c)",240,0.,12.);
1006 hPtH->SetMinimum(0);
1007 fOutputStudy->Add(hPtH);
1008
f80e7bba 1009 namePlot = "hist_Pt_K0_Bin"; namePlot+=i;
a2ad7da1 1010 TH1F *hPtK = new TH1F(namePlot.Data(), "Kaons pT (in D0 evs); p_{T} (GeV/c)",240,0.,12.);
1011 hPtK->SetMinimum(0);
1012 fOutputStudy->Add(hPtK);
1013
1014 //D* feeddown pions rejection histos
1015 namePlot = "hDstarPions_Bin"; namePlot+=i;
f80e7bba 1016 TH2F *hDstarPions = new TH2F(namePlot.Data(), "Tracks rejected for D* inv.mass cut; # Tracks",2,0.,2.,300,1.6,2.2);
a2ad7da1 1017 hDstarPions->GetXaxis()->SetBinLabel(1,"Not rejected");
1018 hDstarPions->GetXaxis()->SetBinLabel(2,"Rejected");
1019 hDstarPions->SetMinimum(0);
bce70c96 1020 fOutputStudy->Add(hDstarPions);
1021
1022 }
1023
1024 if(fMixing) {
1025 //THnSparse plots for event mixing!
1026 namePlot="hPhi_K0_Bin";
1027 namePlot+=i;namePlot+="_EvMix";
1028
1029 THnSparseI *hPhiK_EvMix = new THnSparseI(namePlot.Data(), "Az. corr. EvMix; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix);
1030 hPhiK_EvMix->Sumw2();
1031 fOutputCorr->Add(hPhiK_EvMix);
1032
1033 namePlot="hPhi_Kcharg_Bin";
1034 namePlot+=i;namePlot+="_EvMix";
1035
1036 THnSparseI *hPhiH_EvMix = new THnSparseI(namePlot.Data(), "Az. corr. EvMix; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix);
1037 hPhiH_EvMix->Sumw2();
1038 fOutputCorr->Add(hPhiH_EvMix);
1039
1040 namePlot="hPhi_Charg_Bin";
1041 namePlot+=i;namePlot+="_EvMix";
1042
1043 THnSparseI *hPhiC_EvMix = new THnSparseI(namePlot.Data(), "Az. corr. EvMix; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix);
1044 hPhiC_EvMix->Sumw2();
1045 fOutputCorr->Add(hPhiC_EvMix);
1046 }
a2ad7da1 1047 }
1048
1049 //out of bin loop
bce70c96 1050 TH1F *hCountC = new TH1F("hist_Count_Charg", "Charged track counter; # Tracks",100,0.,100.);
1051 hCountC->SetMinimum(0);
1052 fOutputStudy->Add(hCountC);
1053
1054 TH1F *hCountH = new TH1F("hist_Count_Kcharg", "Hadrons counter; # Tracks",20,0.,20.);
1055 hCountH->SetMinimum(0);
1056 fOutputStudy->Add(hCountH);
1057
1058 TH1F *hCountK = new TH1F("hist_Count_K0", "Kaons counter; # Tracks",20,0.,20.);
1059 hCountK->SetMinimum(0);
1060 fOutputStudy->Add(hCountK);
1061
a2ad7da1 1062 if (fFillGlobal) { //all-events plots
1063 //pt distributions
1064 TH1F *hPtCAll = new TH1F("hist_Pt_Charg_AllEv", "Charged track pT (All); p_{T} (GeV/c)",240,0.,12.);
1065 hPtCAll->SetMinimum(0);
1066 fOutputStudy->Add(hPtCAll);
1067
1068 TH1F *hPtHAll = new TH1F("hist_Pt_Kcharg_AllEv", "Hadrons pT (All); p_{T} (GeV/c)",240,0.,12.);
1069 hPtHAll->SetMinimum(0);
1070 fOutputStudy->Add(hPtHAll);
1071
bce70c96 1072 TH1F *hPtKAll = new TH1F("hist_Pt_K0_AllEv", "Kaons pT (All); p_{T} (GeV/c)",240,0.,12.);
a2ad7da1 1073 hPtKAll->SetMinimum(0);
1074 fOutputStudy->Add(hPtKAll);
1075
bce70c96 1076 //phi distributions
1077 TH1F *hPhiDistCAll = new TH1F("hist_PhiDistr_Charg", "Charged track phi distr. (All); p_{T} (GeV/c)",32,-1.6,4.8);
1078 hPhiDistCAll->SetMinimum(0);
1079 fOutputStudy->Add(hPhiDistCAll);
1080
1081 TH1F *hPhiDistHAll = new TH1F("hist_PhiDistr_Kcharg", "Hadrons phi distr. (All); p_{T} (GeV/c)",32,-1.6,4.8);
1082 hPhiDistHAll->SetMinimum(0);
1083 fOutputStudy->Add(hPhiDistHAll);
1084
1085 TH1F *hPhiDistKAll = new TH1F("hist_PhiDistr_K0", "Kaons phi distr. (All); p_{T} (GeV/c)",32,-1.6,4.8);
1086 hPhiDistKAll->SetMinimum(0);
1087 fOutputStudy->Add(hPhiDistKAll);
1088
1089 TH1F *hPhiDistDAll = new TH1F("hist_PhiDistr_D0", "D^{0} phi distr. (All); p_{T} (GeV/c)",32,-1.6,4.8);
1090 hPhiDistDAll->SetMinimum(0);
1091 fOutputStudy->Add(hPhiDistDAll);
1092
a2ad7da1 1093 //K0 Invariant Mass plots
1094 TH2F *hK0MassInv = new TH2F("hK0MassInv", "K0 invariant mass; Invariant mass (MeV/c^{2}); pT (GeV/c)",200,0.4,0.6,100,0.,10.);
1095 hK0MassInv->SetMinimum(0);
1096 fOutputStudy->Add(hK0MassInv);
1097 }
1098
1099 //for MC analysis only
1100 if (fReadMC) {
1101
1102 for(Int_t i=0;i<fNPtBinsCorr;i++){
1103
1104 //displacement histos
f80e7bba 1105 namePlot="histDispl_K0_Bin"; namePlot+=i;
a2ad7da1 1106 TH1F *hDisplK = new TH1F(namePlot.Data(), "Kaons Displacement; DCA",150,0.,0.15);
1107 hDisplK->SetMinimum(0);
1108 fOutputStudy->Add(hDisplK);
1109
f80e7bba 1110 namePlot="histDispl_K0_HF_Bin"; namePlot+=i;
a2ad7da1 1111 TH1F *hDisplK_HF = new TH1F(namePlot.Data(), "Kaons Displacement (from HF decay only); DCA",150,0.,0.15);
1112 hDisplK_HF->SetMinimum(0);
1113 fOutputStudy->Add(hDisplK_HF);
1114
1115 namePlot="histDispl_Kcharg_Bin"; namePlot+=i;
1116 TH1F *hDisplHadr = new TH1F(namePlot.Data(), "Hadrons Displacement; DCA",150,0.,0.15);
1117 hDisplHadr->SetMinimum(0);
1118 fOutputStudy->Add(hDisplHadr);
1119
1120 namePlot="histDispl_Kcharg_HF_Bin"; namePlot+=i;
1121 TH1F *hDisplHadr_HF = new TH1F(namePlot.Data(), "Hadrons Displacement (from HF decay only); DCA",150,0.,0.15);
1122 hDisplHadr_HF->SetMinimum(0);
1123 fOutputStudy->Add(hDisplHadr_HF);
1124
1125 namePlot="histDispl_Charg_Bin"; namePlot+=i;
1126 TH1F *hDisplCharg = new TH1F(namePlot.Data(), "Charged tracks Displacement; DCA",150,0.,0.15);
1127 hDisplCharg->SetMinimum(0);
1128 fOutputStudy->Add(hDisplCharg);
1129
1130 namePlot="histDispl_Charg_HF_Bin"; namePlot+=i;
1131 TH1F *hDisplCharg_HF = new TH1F(namePlot.Data(), "Charged tracks Displacement (from HF decay only); DCA",150,0.,0.15);
1132 hDisplCharg_HF->SetMinimum(0);
1133 fOutputStudy->Add(hDisplCharg_HF);
1134
f80e7bba 1135 namePlot="histDispl_K0_From_c_Bin"; namePlot+=i;
a2ad7da1 1136 TH1F *hDisplK_c = new TH1F(namePlot.Data(), "Kaons Displacement - c origin; DCA",150,0.,0.15);
1137 hDisplK_c->SetMinimum(0);
1138 fOutputStudy->Add(hDisplK_c);
1139
f80e7bba 1140 namePlot="histDispl_K0_HF_From_c_Bin"; namePlot+=i;
a2ad7da1 1141 TH1F *hDisplK_HF_c = new TH1F(namePlot.Data(), "Kaons Displacement (from HF decay only) - c origin; DCA",150,0.,0.15);
1142 hDisplK_HF_c->SetMinimum(0);
1143 fOutputStudy->Add(hDisplK_HF_c);
1144
1145 namePlot="histDispl_Kcharg_From_c_Bin"; namePlot+=i;
1146 TH1F *hDisplHadr_c = new TH1F(namePlot.Data(), "Hadrons Displacement - c origin; DCA",150,0.,0.15);
1147 hDisplHadr_c->SetMinimum(0);
1148 fOutputStudy->Add(hDisplHadr_c);
1149
1150 namePlot="histDispl_Kcharg_HF_From_c_Bin"; namePlot+=i;
1151 TH1F *hDisplHadr_HF_c = new TH1F(namePlot.Data(), "Hadrons Displacement (from HF decay only) - c origin; DCA",150,0.,0.15);
1152 hDisplHadr_HF_c->SetMinimum(0);
1153 fOutputStudy->Add(hDisplHadr_HF_c);
1154
1155 namePlot="histDispl_Charg_From_c_Bin"; namePlot+=i;
1156 TH1F *hDisplCharg_c = new TH1F(namePlot.Data(), "Charged tracks Displacement - c origin; DCA",150,0.,0.15);
1157 hDisplCharg_c->Sumw2();
1158 hDisplCharg_c->SetMinimum(0);
1159 fOutputStudy->Add(hDisplCharg_c);
1160
1161 namePlot="histDispl_Charg_HF_From_c_Bin"; namePlot+=i;
1162 TH1F *hDisplCharg_HF_c = new TH1F(namePlot.Data(), "Charged tracks Displacement (from HF decay only) - c origin; DCA",150,0.,0.15);
1163 hDisplCharg_HF_c->SetMinimum(0);
1164 fOutputStudy->Add(hDisplCharg_HF_c);
1165
f80e7bba 1166 namePlot="histDispl_K0_From_b_Bin"; namePlot+=i;
a2ad7da1 1167 TH1F *hDisplK_b = new TH1F(namePlot.Data(), "Kaons Displacement - b origin; DCA",150,0.,0.15);
1168 hDisplK_b->SetMinimum(0);
1169 fOutputStudy->Add(hDisplK_b);
1170
f80e7bba 1171 namePlot="histDispl_K0_HF_From_b_Bin"; namePlot+=i;
a2ad7da1 1172 TH1F *hDisplK_HF_b = new TH1F(namePlot.Data(), "Kaons Displacement (from HF decay only) - b origin; DCA",150,0.,0.15);
1173 hDisplK_HF_b->SetMinimum(0);
1174 fOutputStudy->Add(hDisplK_HF_b);
1175
1176 namePlot="histDispl_Kcharg_From_b_Bin"; namePlot+=i;
1177 TH1F *hDisplHadr_b = new TH1F(namePlot.Data(), "Hadrons Displacement - b origin; DCA",150,0.,0.15);
1178 hDisplHadr_b->SetMinimum(0);
1179 fOutputStudy->Add(hDisplHadr_b);
1180
1181 namePlot="histDispl_Kcharg_HF_From_b_Bin"; namePlot+=i;
1182 TH1F *hDisplHadr_HF_b = new TH1F(namePlot.Data(), "Hadrons Displacement (from HF decay only) - b origin; DCA",150,0.,0.15);
1183 hDisplHadr_HF_b->SetMinimum(0);
1184 fOutputStudy->Add(hDisplHadr_HF_b);
1185
1186 namePlot="histDispl_Charg_From_b_Bin"; namePlot+=i;
1187 TH1F *hDisplCharg_b = new TH1F(namePlot.Data(), "Charged tracks Displacement - b origin; DCA",150,0.,0.15);
1188 hDisplCharg_b->SetMinimum(0);
1189 fOutputStudy->Add(hDisplCharg_b);
1190
1191 namePlot="histDispl_Charg_HF_From_b_Bin"; namePlot+=i;
1192 TH1F *hDisplCharg_HF_b = new TH1F(namePlot.Data(), "Charged tracks Displacement (from HF decay only) - b origin; DCA",150,0.,0.15);
1193 hDisplCharg_HF_b->SetMinimum(0);
1194 fOutputStudy->Add(hDisplCharg_HF_b);
1195
1196 //origin of tracks histos
1197 namePlot="histOrig_Charg_Bin"; namePlot+=i;
1198 TH1F *hOrigin_Charm = new TH1F(namePlot.Data(), "Origin of charged tracks",9,0.,9.);
1199 hOrigin_Charm->SetMinimum(0);
1200 hOrigin_Charm->GetXaxis()->SetBinLabel(1,"Not HF");
1201 hOrigin_Charm->GetXaxis()->SetBinLabel(2,"D->#");
1202 hOrigin_Charm->GetXaxis()->SetBinLabel(3,"D->X->#");
1203 hOrigin_Charm->GetXaxis()->SetBinLabel(4,"B->#");
1204 hOrigin_Charm->GetXaxis()->SetBinLabel(5,"B->X-># (X!=D)");
1205 hOrigin_Charm->GetXaxis()->SetBinLabel(6,"B->D->#");
1206 hOrigin_Charm->GetXaxis()->SetBinLabel(7,"B->D->X->#");
1207 hOrigin_Charm->GetXaxis()->SetBinLabel(8,"c hadr.");
1208 hOrigin_Charm->GetXaxis()->SetBinLabel(9,"b hadr.");
1209 fOutputStudy->Add(hOrigin_Charm);
1210
1211 namePlot="histOrig_Kcharg_Bin"; namePlot+=i;
1212 TH1F *hOrigin_Kcharg = new TH1F(namePlot.Data(), "Origin of hadrons",9,0.,9.);
1213 hOrigin_Kcharg->SetMinimum(0);
1214 hOrigin_Kcharg->GetXaxis()->SetBinLabel(1,"Not HF");
1215 hOrigin_Kcharg->GetXaxis()->SetBinLabel(2,"D->#");
1216 hOrigin_Kcharg->GetXaxis()->SetBinLabel(3,"D->X->#");
1217 hOrigin_Kcharg->GetXaxis()->SetBinLabel(4,"B->#");
1218 hOrigin_Kcharg->GetXaxis()->SetBinLabel(5,"B->X-># (X!=D)");
1219 hOrigin_Kcharg->GetXaxis()->SetBinLabel(6,"B->D->#");
1220 hOrigin_Kcharg->GetXaxis()->SetBinLabel(7,"B->D->X->#");
1221 hOrigin_Kcharg->GetXaxis()->SetBinLabel(8,"c hadr.");
1222 hOrigin_Kcharg->GetXaxis()->SetBinLabel(9,"b hadr.");
1223 fOutputStudy->Add(hOrigin_Kcharg);
1224
f80e7bba 1225 namePlot="histOrig_K0_Bin"; namePlot+=i;
a2ad7da1 1226 TH1F *hOrigin_K = new TH1F(namePlot.Data(), "Origin of kaons",9,0.,9.);
1227 hOrigin_K->SetMinimum(0);
1228 hOrigin_K->GetXaxis()->SetBinLabel(1,"Not HF");
1229 hOrigin_K->GetXaxis()->SetBinLabel(2,"D->#");
1230 hOrigin_K->GetXaxis()->SetBinLabel(3,"D->X->#");
1231 hOrigin_K->GetXaxis()->SetBinLabel(4,"B->#");
1232 hOrigin_K->GetXaxis()->SetBinLabel(5,"B->X-># (X!=D)");
1233 hOrigin_K->GetXaxis()->SetBinLabel(6,"B->D->#");
1234 hOrigin_K->GetXaxis()->SetBinLabel(7,"B->D->X->#");
1235 hOrigin_K->GetXaxis()->SetBinLabel(8,"c hadr.");
1236 hOrigin_K->GetXaxis()->SetBinLabel(9,"b hadr.");
1237 fOutputStudy->Add(hOrigin_K);
1238
1239 //origin of D0 histos
1240 namePlot="histOrig_D0_Bin"; namePlot+=i;
1241 TH1F *hOrigin_D0 = new TH1F(namePlot.Data(), "Origin of D0",2,0.,2.);
1242 hOrigin_D0->SetMinimum(0);
1243 hOrigin_D0->GetXaxis()->SetBinLabel(1,"From c");
1244 hOrigin_D0->GetXaxis()->SetBinLabel(2,"From b");
1245 fOutputStudy->Add(hOrigin_D0);
1246 }
1247 }
1248
1249}
1250
1251//________________________________________________________________________
bce70c96 1252void AliAnalysisTaskSED0Correlations::CalculateCorrelations(AliAODRecoDecayHF2Prong* d, Int_t labD0, TClonesArray* mcArray) {
a2ad7da1 1253//
1254// Method for correlations D0-hadrons study
1255//
bce70c96 1256 Int_t N_Charg = 0, N_KCharg = 0, N_Kaons = 0;
1257 Double_t mD0, mD0bar;
1258 Int_t origD0 = 0, PDGD0 = 0;
a2ad7da1 1259 d->InvMassD0(mD0,mD0bar);
1260 Int_t ptbin = PtBinCorr(d->Pt());
1261 if(ptbin < 0) return;
a2ad7da1 1262
bce70c96 1263 //Origin of D0
1264 TString orig="";
1265 if(fReadMC) {
1266 origD0=CheckD0Origin(mcArray,(AliAODMCParticle*)mcArray->At(labD0));
1267 PDGD0 = ((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode();
1268 switch (CheckD0Origin(mcArray,(AliAODMCParticle*)mcArray->At(labD0))) {
1269 case 4:
1270 orig = "_From_c";
1271 ((TH1F*)fOutputStudy->FindObject(Form("histOrig_D0_Bin%d",ptbin)))->Fill(0.);
1272 break;
1273 case 5:
1274 orig = "_From_b";
1275 ((TH1F*)fOutputStudy->FindObject(Form("histOrig_D0_Bin%d",ptbin)))->Fill(1.);
1276 break;
1277 default:
1278 return;
1279 }
1280 }
1281
1282 Double_t highPt = 0; Double_t lead[3] = {0,0,0}; //infos for leading particle (pt,deltaphi)
a2ad7da1 1283
bce70c96 1284 //loop over the tracks in the pool
ad8d2dfd 1285 Bool_t execPoolTr = fCorrelatorTr->ProcessEventPool(); //pool is ready? (only in ME, in SE returns kFALSE)
1286 Bool_t execPoolKc = fCorrelatorKc->ProcessEventPool(); //pool is ready? (only in ME, in SE returns kFALSE)
1287 Bool_t execPoolK0 = fCorrelatorK0->ProcessEventPool(); //pool is ready? (only in ME, in SE returns kFALSE)
bce70c96 1288
1289 Int_t NofEventsinPool = 1;
ad8d2dfd 1290 if(fMixing) {
1291 NofEventsinPool = fCorrelatorTr->GetNofEventsInPool();
1292 if(!execPoolTr) {
1293 AliInfo("Mixed event analysis: track pool is not ready");
1294 NofEventsinPool = 0;
1295 }
1296 }
1297
1298 //Charged tracks
bce70c96 1299 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)
1300 Bool_t analyzetracksTr = fCorrelatorTr->ProcessAssociatedTracks(jMix);// process all the tracks in the aodEvent, by applying the selection cuts
ad8d2dfd 1301 if(!analyzetracksTr) {
bce70c96 1302 AliInfo("AliHFCorrelator::Cannot process the track array");
1303 continue;
1304 }
ad8d2dfd 1305
bce70c96 1306 for(Int_t iTrack = 0; iTrack<fCorrelatorTr->GetNofTracks(); iTrack++){ // looping on track candidates
1307
1308 Bool_t runcorrelation = fCorrelatorTr->Correlate(iTrack);
1309 if(!runcorrelation) continue;
1310
1311 AliReducedParticle* track = fCorrelatorTr->GetAssociatedParticle();
1312
1313 if(!fMixing) {
1314 if(!track->CheckSoftPi()) { //removal of soft pions
1315 if (fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(1.,mD0);
1316 if (fIsSelectedCandidate >= 2) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(1.,mD0bar);
1317 continue;
1318 } else { //not a soft pion
1319 if (fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(0.,mD0);
1320 if (fIsSelectedCandidate >= 2) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(0.,mD0bar);
1321 }
bce70c96 1322 Int_t idDaughs[2] = {((AliVTrack*)d->GetDaughter(0))->GetID(),((AliVTrack*)d->GetDaughter(1))->GetID()}; //IDs of daughters to be skipped
1323 if(track->GetID() == idDaughs[0] || track->GetID() == idDaughs[1]) continue; //discards daughters of candidate
ad8d2dfd 1324 }
1325 if(track->Pt() < fPtThreshLow.at(ptbin) || track->Pt() > fPtThreshUp.at(ptbin)) continue; //discard tracks outside pt range for hadrons/K
bce70c96 1326
1327 FillSparsePlots(mcArray,d,origD0,PDGD0,track,kTrack); //fills for charged tracks
1328
1329 if(!fMixing) N_Charg++;
1330
1331 //retrieving leading info...
1332 if(track->Pt() > highPt) {
1333 lead[0] = fCorrelatorTr->GetDeltaPhi();
1334 lead[1] = fCorrelatorTr->GetDeltaEta();
1335 lead[2] = fReadMC ? CheckTrackOrigin(mcArray,(AliAODMCParticle*)mcArray->At(track->GetLabel())) : 0;
1336 highPt = track->Pt();
a2ad7da1 1337 }
bce70c96 1338
1339 } // end of tracks loop
ad8d2dfd 1340 } //end of event loop for fCorrelatorTr
1341
1342 if(fMixing) {
1343 NofEventsinPool = fCorrelatorKc->GetNofEventsInPool();
1344 if(!execPoolKc) {
1345 AliInfo("Mixed event analysis: K+/- pool is not ready");
1346 NofEventsinPool = 0;
1347 }
1348 }
1349
1350 //Charged Kaons loop
1351 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)
1352 Bool_t analyzetracksKc = fCorrelatorKc->ProcessAssociatedTracks(jMix);
1353 if(!analyzetracksKc) {
1354 AliInfo("AliHFCorrelator::Cannot process the K+/- array");
1355 continue;
1356 }
bce70c96 1357
bce70c96 1358 for(Int_t iTrack = 0; iTrack<fCorrelatorKc->GetNofTracks(); iTrack++){ // looping on charged kaons candidates
1359
1360 Bool_t runcorrelation = fCorrelatorKc->Correlate(iTrack);
1361 if(!runcorrelation) continue;
1362
1363 AliReducedParticle* kCharg = fCorrelatorKc->GetAssociatedParticle();
1364
1365 if(!fMixing) {
ad8d2dfd 1366 if(!kCharg->CheckSoftPi()) { //removal of soft pions
bce70c96 1367 if (fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(1.,mD0);
1368 if (fIsSelectedCandidate >= 2) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(1.,mD0bar);
1369 continue;
1370 } else {
1371 if (fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(0.,mD0);
1372 if (fIsSelectedCandidate >= 2) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(0.,mD0bar);
a2ad7da1 1373 }
bce70c96 1374 Int_t idDaughs[2] = {((AliVTrack*)d->GetDaughter(0))->GetID(),((AliVTrack*)d->GetDaughter(1))->GetID()}; //IDs of daughters to be skipped
1375 if(kCharg->GetID() == idDaughs[0] || kCharg->GetID() == idDaughs[1]) continue; //discards daughters of candidate
ad8d2dfd 1376 }
1377 if(kCharg->Pt() < fPtThreshLow.at(ptbin) || kCharg->Pt() > fPtThreshUp.at(ptbin)) continue; //discard tracks outside pt range for hadrons/K
1378
bce70c96 1379 FillSparsePlots(mcArray,d,origD0,PDGD0,kCharg,kKCharg); //fills for charged tracks
a2ad7da1 1380
bce70c96 1381 if(!fMixing) N_KCharg++;
a2ad7da1 1382
bce70c96 1383 } // end of charged kaons loop
ad8d2dfd 1384 } //end of event loop for fCorrelatorKc
1385
1386 if(fMixing) {
1387 NofEventsinPool = fCorrelatorK0->GetNofEventsInPool();
1388 if(!execPoolK0) {
1389 AliInfo("Mixed event analysis: K0 pool is not ready");
1390 NofEventsinPool = 0;
1391 }
1392 }
1393
1394 //K0 loop
1395 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)
1396 Bool_t analyzetracksK0 = fCorrelatorK0->ProcessAssociatedTracks(jMix);
1397 if(!analyzetracksK0) {
1398 AliInfo("AliHFCorrelator::Cannot process the K0 array");
1399 continue;
1400 }
a2ad7da1 1401
bce70c96 1402 for(Int_t iTrack = 0; iTrack<fCorrelatorK0->GetNofTracks(); iTrack++){ // looping on k0 candidates
1403
1404 Bool_t runcorrelation = fCorrelatorK0->Correlate(iTrack);
1405 if(!runcorrelation) continue;
1406
1407 AliReducedParticle* k0 = fCorrelatorK0->GetAssociatedParticle();
a2ad7da1 1408
bce70c96 1409 if(k0->Pt() < fPtThreshLow.at(ptbin) || k0->Pt() > fPtThreshUp.at(ptbin)) continue; //discard tracks outside pt range for hadrons/K
bce70c96 1410
1411 FillSparsePlots(mcArray,d,origD0,PDGD0,k0,kK0); //fills for charged tracks
1412
1413 if(!fMixing) N_Kaons++;
1414
1415 } // end of charged kaons loop
ad8d2dfd 1416 } //end of event loop for fCorrelatorK0
bce70c96 1417
1418 //leading track correlations fill
ad8d2dfd 1419 if(!fMixing) {
1420 if(fReadMC) {
1421 if(((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode()==421) { //D0
1422 ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(lead[1],mD0); //c and b D0
1423 ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead%s_Bin%d",orig.Data(),ptbin)))->Fill(lead[0],mD0); //c or b D0
1424 if(origD0==4&&(int)lead[2]>=1&&(int)lead[2]<=2) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(lead[0],mD0);
1425 if(origD0==5&&(int)lead[2]>=3&&(int)lead[2]<=6) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(lead[0],mD0);
1426 }
1427 if(((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode()==-421) { //D0bar
1428 ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(lead[1],mD0bar);
1429 ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead%s_Bin%d",orig.Data(),ptbin)))->Fill(lead[0],mD0bar); //c or b D0
1430 if(origD0==4&&(int)lead[2]>=1&&(int)lead[2]<=2) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(lead[0],mD0bar);
1431 if(origD0==5&&(int)lead[2]>=3&&(int)lead[2]<=6) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(lead[0],mD0bar);
1432 }
1433 } else {
1434 if(fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(lead[0],mD0); //c and b D0
1435 if(fIsSelectedCandidate == 2 || fIsSelectedCandidate == 3) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(lead[0],mD0bar);
a2ad7da1 1436 }
1437
ad8d2dfd 1438 //Fill of count histograms
1439 if (!fAlreadyFilled) {
1440 ((TH1F*)fOutputStudy->FindObject("hist_Count_Charg"))->Fill(N_Charg);
1441 ((TH1F*)fOutputStudy->FindObject("hist_Count_Kcharg"))->Fill(N_KCharg);
1442 ((TH1F*)fOutputStudy->FindObject("hist_Count_K0"))->Fill(N_Kaons);
1443 }
bce70c96 1444 }
a2ad7da1 1445
bce70c96 1446 fAlreadyFilled=kTRUE; //at least a D0 analyzed in the event; distribution plots already filled
a2ad7da1 1447
bce70c96 1448}
1449
1450//________________________________________________________________________
1451void AliAnalysisTaskSED0Correlations::FillSparsePlots(TClonesArray* mcArray, AliAODRecoDecayHF2Prong *d, Int_t origD0, Int_t PdgD0, AliReducedParticle* track, Int_t type) {
1452 //
1453 //fills the THnSparse for correlations, calculating the variables
1454 //
1455
1456 //Initialization of variables
1457 Int_t ptbin = PtBinCorr(d->Pt());
1458 if(ptbin < 0) return;
1459 Double_t mD0, mD0bar, deltaphi = 0., deltaeta = 0.;
1460 d->InvMassD0(mD0,mD0bar);
1461
1462 Double_t ptTrack = track->Pt();
1463 Double_t d0Track = type!=kK0 ? track->GetImpPar() : 0.;
1464 Double_t phiTr = track->Phi();
1465 Double_t origTr = fReadMC ? CheckTrackOrigin(mcArray,(AliAODMCParticle*)mcArray->At(track->GetLabel())) : 0;
1466
1467 TString part = "", orig = "";
1468
1469 switch (type) {
1470 case(kTrack): {
1471 part = "Charg";
1472 deltaphi = fCorrelatorTr->GetDeltaPhi();
1473 deltaeta = fCorrelatorTr->GetDeltaEta();
1474 break;
1475 }
1476 case(kKCharg): {
1477 part = "Kcharg";
1478 deltaphi = fCorrelatorKc->GetDeltaPhi();
1479 deltaeta = fCorrelatorKc->GetDeltaEta();
1480 break;
1481 }
1482 case(kK0): {
1483 part = "K0";
1484 deltaphi = fCorrelatorK0->GetDeltaPhi();
1485 deltaeta = fCorrelatorK0->GetDeltaEta();
1486 break;
a2ad7da1 1487 }
bce70c96 1488 }
bce70c96 1489
1490 if(fMixing == kSE) {
1491
ad8d2dfd 1492 //Fixes limits; needed to include overflow into THnSparse projections!
1493 Double_t pTorig = track->Pt();
1494 Double_t d0orig = track->GetImpPar();
1495 Double_t ptLim_Sparse = ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->GetAxis(2)->GetXmax(); //all plots have same axes...
1496 Double_t displLim_Sparse = ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->GetAxis(3)->GetXmax();
1497 Double_t EtaLim_Sparse = ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->GetAxis(4)->GetXmax();
1498 if(ptTrack > ptLim_Sparse) ptTrack = ptLim_Sparse-0.01;
1499 if(d0Track > displLim_Sparse) d0Track = (displLim_Sparse-0.001);
1500 if(deltaeta > EtaLim_Sparse) deltaeta = EtaLim_Sparse-0.01;
1501 if(deltaeta < -EtaLim_Sparse) deltaeta = -EtaLim_Sparse+0.01;
1502
1503 //variables for filling histos
1504 Double_t fillSpPhiD0[5] = {deltaphi,mD0,ptTrack,d0Track,deltaeta};
1505 Double_t fillSpPhiD0bar[5] = {deltaphi,mD0bar,ptTrack,d0Track,deltaeta};
a2ad7da1 1506
bce70c96 1507 if(fReadMC == 0) {
1508 //sparse fill for data (tracks, K+-, K0) + weighted
1509 if(fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) { //D0
1510 ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0);
1511 ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_Bin%d",ptbin)))->Fill(deltaphi,mD0,pTorig);
a2ad7da1 1512 }
bce70c96 1513 if(fIsSelectedCandidate == 2 || fIsSelectedCandidate == 3) { //D0bar
1514 ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0bar);
1515 ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_Bin%d",ptbin)))->Fill(deltaphi,mD0bar,pTorig);
a2ad7da1 1516 }
1517 if(!fAlreadyFilled) {
bce70c96 1518 ((TH1F*)fOutputStudy->FindObject(Form("hist_Pt_%s_Bin%d",part.Data(),ptbin)))->Fill(pTorig);
1519 ((TH1F*)fOutputStudy->FindObject(Form("hist_PhiDistr_%s",part.Data())))->Fill(phiTr);
a2ad7da1 1520 }
bce70c96 1521 }
a2ad7da1 1522
bce70c96 1523 if(fReadMC) {
a2ad7da1 1524
bce70c96 1525 if(origD0==4) {orig = "_From_c";} else {orig = "_From_b";}
a2ad7da1 1526
bce70c96 1527 //sparse fill for data (tracks, K+-, K0) + weighted
1528 if(PdgD0==421) { //D0 (from MCTruth)
1529 ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0);
1530 ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0);
1531 if(origD0==4&&origTr>=1&&origTr<=2) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_HF%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0);
1532 if(origD0==5&&origTr>=3&&origTr<=6) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_HF%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0);
1533 ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_Bin%d",ptbin)))->Fill(deltaphi,mD0,pTorig);
1534 ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0,pTorig);
1535 if(origD0==4&&origTr>=1&&origTr<=2) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0,pTorig);
1536 if(origD0==5&&origTr>=3&&origTr<=6) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0,pTorig);
a2ad7da1 1537 }
bce70c96 1538 if(PdgD0==-421) { //D0bar
1539 ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0bar);
1540 ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0bar);
1541 if(origD0==4&&origTr>=1&&origTr<=2) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_HF%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0bar);
1542 if(origD0==5&&origTr>=3&&origTr<=6) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_HF%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0bar);
1543 ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_Bin%d",ptbin)))->Fill(deltaphi,mD0bar,pTorig);
1544 ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0bar,pTorig);
1545 if(origD0==4&&origTr>=1&&origTr<=2) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0bar,pTorig);
1546 if(origD0==5&&origTr>=3&&origTr<=6) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0bar,pTorig);
1547 }
a2ad7da1 1548 if(!fAlreadyFilled) {
bce70c96 1549 ((TH1F*)fOutputStudy->FindObject(Form("histDispl_%s_Bin%d",part.Data(),ptbin)))->Fill(d0orig); //Fills displacement histos
1550 if (origTr>=1&&origTr<=6) ((TH1F*)fOutputStudy->FindObject(Form("histDispl_%s_HF_Bin%d",part.Data(),ptbin)))->Fill(d0orig);
1551 if (origTr>=1&&origTr<=6) ((TH1F*)fOutputStudy->FindObject(Form("histDispl_%s_HF%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(d0orig);
1552 ((TH1F*)fOutputStudy->FindObject(Form("histDispl_%s%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(d0orig); //Fills displacement histos
1553 ((TH1F*)fOutputStudy->FindObject(Form("hist_Pt_%s_Bin%d",part.Data(),ptbin)))->Fill(pTorig);
1554 ((TH1F*)fOutputStudy->FindObject(Form("histOrig_%s_Bin%d",part.Data(),ptbin)))->Fill(origTr);
1555 ((TH1F*)fOutputStudy->FindObject(Form("hist_PhiDistr_%s",part.Data())))->Fill(phiTr);
a2ad7da1 1556 }
bce70c96 1557 }//end MC case
a2ad7da1 1558
bce70c96 1559 } //end of SE fill
a2ad7da1 1560
bce70c96 1561 if(fMixing == kME) {
a2ad7da1 1562
ad8d2dfd 1563 //Fixes limits; needed to include overflow into THnSparse projections!
1564 Double_t EtaLim_Sparse = ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d_EvMix",ptbin)))->GetAxis(2)->GetXmax();
1565 if(deltaeta > EtaLim_Sparse) deltaeta = EtaLim_Sparse-0.01;
1566 if(deltaeta < -EtaLim_Sparse) deltaeta = -EtaLim_Sparse+0.01;
1567
1568 //variables for filling histos
1569 Double_t fillSpPhiD0[3] = {deltaphi,mD0,deltaeta};
1570 Double_t fillSpPhiD0bar[3] = {deltaphi,mD0bar,deltaeta};
a2ad7da1 1571
bce70c96 1572 if(fReadMC == 0) {
1573 //sparse fill for data (tracks, K+-, K0)
1574 if(fIsSelectedCandidate == 1||fIsSelectedCandidate == 3) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d_EvMix",part.Data(),ptbin)))->Fill(fillSpPhiD0);
1575 if(fIsSelectedCandidate == 2||fIsSelectedCandidate == 3) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d_EvMix",part.Data(),ptbin)))->Fill(fillSpPhiD0bar);
1576 }
1577 if(fReadMC == 1) {
1578 //sparse fill for data (tracks, K+-, K0)
1579 if(PdgD0==421) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d_EvMix",part.Data(),ptbin)))->Fill(fillSpPhiD0);
1580 if(PdgD0==-421) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d_EvMix",part.Data(),ptbin)))->Fill(fillSpPhiD0bar);
1581 }//end MC case
a2ad7da1 1582
bce70c96 1583 } //end of ME fill
1584
1585 return;
a2ad7da1 1586}
1587
1588//_________________________________________________________________________________________________
1589Int_t AliAnalysisTaskSED0Correlations::CheckTrackOrigin(TClonesArray* arrayMC, AliAODMCParticle *mcPartCandidate) const {
1590 //
1591 // checks on particle (#) origin:
1592 // 0) Not HF
1593 // 1) D->#
1594 // 2) D->X->#
1595 // 3) B->#
1596 // 4) B->X-># (X!=D)
1597 // 5) B->D->#
1598 // 6) B->D->X->#
1599 // 7) c hadronization
1600 // 8) b hadronization
1601 //
bce70c96 1602 if(fDebug>2) printf("AliAnalysisTaskSED0Correlations::CheckTrkOrigin() \n");
a2ad7da1 1603
1604 Int_t pdgGranma = 0;
1605 Int_t mother = 0;
1606 mother = mcPartCandidate->GetMother();
1607 Int_t istep = 0;
1608 Int_t abspdgGranma =0;
1609 Bool_t isFromB=kFALSE;
1610 Bool_t isDdaugh=kFALSE;
1611 Bool_t isDchaindaugh=kFALSE;
1612 Bool_t isBdaugh=kFALSE;
1613 Bool_t isBchaindaugh=kFALSE;
1614 Bool_t isQuarkFound=kFALSE;
1615
1616 while (mother > 0){
1617 istep++;
1618 AliAODMCParticle* mcGranma = dynamic_cast<AliAODMCParticle*>(arrayMC->At(mother));
1619 if (mcGranma){
1620 pdgGranma = mcGranma->GetPdgCode();
1621 abspdgGranma = TMath::Abs(pdgGranma);
1622 if ((abspdgGranma > 500 && abspdgGranma < 600) || (abspdgGranma > 5000 && abspdgGranma < 6000)){
1623 isBchaindaugh=kTRUE;
1624 if(istep==1) isBdaugh=kTRUE;
1625 }
1626 if ((abspdgGranma > 400 && abspdgGranma < 500) || (abspdgGranma > 4000 && abspdgGranma < 5000)){
1627 isDchaindaugh=kTRUE;
1628 if(istep==1) isDdaugh=kTRUE;
1629 }
1630 if(abspdgGranma==4 || abspdgGranma==5) {isQuarkFound=kTRUE; if(abspdgGranma==5) isFromB = kTRUE;}
1631 mother = mcGranma->GetMother();
1632 }else{
1633 AliError("Failed casting the mother particle!");
1634 break;
1635 }
1636 }
1637
1638 //decides what to return based on the flag status
1639 if(isQuarkFound) {
1640 if(!isFromB) { //charm
1641 if(isDdaugh) return 1; //charm immediate
1642 else if(isDchaindaugh) return 2; //charm chain
1643 else return 7; //charm hadronization
1644 }
1645 else { //beauty
1646 if(isBdaugh) return 3; //b immediate
1647 else if(isBchaindaugh) { //b chain
1648 if(isDchaindaugh) {
1649 if(isDdaugh) return 5; //d immediate
1650 return 6; //d chain
1651 }
1652 else return 4; //b, not d
1653 }
1654 else return 8; //b hadronization
1655 }
1656 }
1657 else return 0; //no HF quark
1658}
1659
a2ad7da1 1660//________________________________________________________________________
1661Int_t AliAnalysisTaskSED0Correlations::PtBinCorr(Double_t pt) const {
1662 //
1663 //give the pt bin where the pt lies.
1664 //
1665 Int_t ptbin=-1;
1666 if(pt<fBinLimsCorr.at(0)) return ptbin; //out of bounds
1667
1668 Int_t i = 0;
1669 while(pt>fBinLimsCorr.at(i)) {ptbin=i; i++;}
1670
1671 return ptbin;
1672}
1673
a2ad7da1 1674//---------------------------------------------------------------------------
1675Bool_t AliAnalysisTaskSED0Correlations::SelectV0(AliAODv0* v0, AliAODVertex *vtx, Int_t opt, Int_t idArrayV0[][2]) const
1676{
1677 //
1678 // Selection for K0 hypotheses
1679 // options: 1 = selects mass invariant about 3 sigma inside the peak + threshold of 0.3 GeV
1680 // 2 = no previous selections
1681
1682 if(!fCutsTracks->IsKZeroSelected(v0,vtx)) return kFALSE;
1683
1684 AliAODTrack *v0Daug1 = (AliAODTrack*)v0->GetDaughter(0);
1685 AliAODTrack *v0Daug2 = (AliAODTrack*)v0->GetDaughter(1);
1686
1687 if(opt==1) { //additional cuts for correlations (V0 has to be closer than 3 sigma from K0 mass)
bce70c96 1688 if(TMath::Abs(v0->MassK0Short()-0.4976) > 3*0.004) return kFALSE;
a2ad7da1 1689 }
1690
1691 //This part removes double counting for swapped tracks!
1692 Int_t i = 0; //while loop (until the last-written entry pair of ID!
1693 while(idArrayV0[i][0]!=-2 && idArrayV0[i][1]!=-2) {
1694 if((v0Daug1->GetID()==idArrayV0[i][0] && v0Daug2->GetID()==idArrayV0[i][1])||
1695 (v0Daug1->GetID()==idArrayV0[i][1] && v0Daug2->GetID()==idArrayV0[i][0])) return kFALSE;
1696 i++;
1697 }
1698 idArrayV0[i][0]=v0Daug1->GetID();
1699 idArrayV0[i][1]=v0Daug2->GetID();
1700
1701 return kTRUE;
1702}
1703
1704//________________________________________________________________________
1705void AliAnalysisTaskSED0Correlations::PrintBinsAndLimits() {
1706
1707 cout << "--------------------------\n";
1708 cout << "PtBins = " << fNPtBinsCorr << "\n";
1709 cout << "PtBin limits--------------\n";
1710 for (int i=0; i<fNPtBinsCorr; i++) {
1711 cout << "Bin "<<i+1<<" = "<<fBinLimsCorr.at(i)<<" to "<<fBinLimsCorr.at(i+1)<<"\n";
1712 }
1713 cout << "\n--------------------------\n";
1714 cout << "PtBin tresh. tracks low---\n";
1715 for (int i=0; i<fNPtBinsCorr; i++) {
1716 cout << fPtThreshLow.at(i) << "; ";
1717 }
1718 cout << "PtBin tresh. tracks up----\n";
1719 for (int i=0; i<fNPtBinsCorr; i++) {
1720 cout << fPtThreshUp.at(i) << "; ";
1721 }
1722 cout << "\n--------------------------\n";
1723 cout << "MC Truth = "<<fReadMC<<"\n";
1724 cout << "--------------------------\n";
bce70c96 1725 cout << "Ev Mixing = "<<fMixing<<"\n";
1726 cout << "--------------------------\n";
a2ad7da1 1727}
1728