]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGHF/correlationHF/AliAnalysisTaskSED0Correlations.cxx
New class for steering D-hadron correlation analysis (Sandro)
[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
1050 TH1F *hRejTracks = new TH1F("hRejTracks", "Tracks accepted/rejected; # Tracks",2,0.,2.);
1051 hRejTracks->SetMinimum(0);
1052 hRejTracks->GetXaxis()->SetBinLabel(1,"Accepted");
1053 hRejTracks->GetXaxis()->SetBinLabel(2,"Rejected");
1054 fOutputStudy->Add(hRejTracks);
1055
bce70c96 1056 TH1F *hCountC = new TH1F("hist_Count_Charg", "Charged track counter; # Tracks",100,0.,100.);
1057 hCountC->SetMinimum(0);
1058 fOutputStudy->Add(hCountC);
1059
1060 TH1F *hCountH = new TH1F("hist_Count_Kcharg", "Hadrons counter; # Tracks",20,0.,20.);
1061 hCountH->SetMinimum(0);
1062 fOutputStudy->Add(hCountH);
1063
1064 TH1F *hCountK = new TH1F("hist_Count_K0", "Kaons counter; # Tracks",20,0.,20.);
1065 hCountK->SetMinimum(0);
1066 fOutputStudy->Add(hCountK);
1067
a2ad7da1 1068 if (fFillGlobal) { //all-events plots
1069 //pt distributions
1070 TH1F *hPtCAll = new TH1F("hist_Pt_Charg_AllEv", "Charged track pT (All); p_{T} (GeV/c)",240,0.,12.);
1071 hPtCAll->SetMinimum(0);
1072 fOutputStudy->Add(hPtCAll);
1073
1074 TH1F *hPtHAll = new TH1F("hist_Pt_Kcharg_AllEv", "Hadrons pT (All); p_{T} (GeV/c)",240,0.,12.);
1075 hPtHAll->SetMinimum(0);
1076 fOutputStudy->Add(hPtHAll);
1077
bce70c96 1078 TH1F *hPtKAll = new TH1F("hist_Pt_K0_AllEv", "Kaons pT (All); p_{T} (GeV/c)",240,0.,12.);
a2ad7da1 1079 hPtKAll->SetMinimum(0);
1080 fOutputStudy->Add(hPtKAll);
1081
bce70c96 1082 //phi distributions
1083 TH1F *hPhiDistCAll = new TH1F("hist_PhiDistr_Charg", "Charged track phi distr. (All); p_{T} (GeV/c)",32,-1.6,4.8);
1084 hPhiDistCAll->SetMinimum(0);
1085 fOutputStudy->Add(hPhiDistCAll);
1086
1087 TH1F *hPhiDistHAll = new TH1F("hist_PhiDistr_Kcharg", "Hadrons phi distr. (All); p_{T} (GeV/c)",32,-1.6,4.8);
1088 hPhiDistHAll->SetMinimum(0);
1089 fOutputStudy->Add(hPhiDistHAll);
1090
1091 TH1F *hPhiDistKAll = new TH1F("hist_PhiDistr_K0", "Kaons phi distr. (All); p_{T} (GeV/c)",32,-1.6,4.8);
1092 hPhiDistKAll->SetMinimum(0);
1093 fOutputStudy->Add(hPhiDistKAll);
1094
1095 TH1F *hPhiDistDAll = new TH1F("hist_PhiDistr_D0", "D^{0} phi distr. (All); p_{T} (GeV/c)",32,-1.6,4.8);
1096 hPhiDistDAll->SetMinimum(0);
1097 fOutputStudy->Add(hPhiDistDAll);
1098
a2ad7da1 1099 //K0 Invariant Mass plots
1100 TH2F *hK0MassInv = new TH2F("hK0MassInv", "K0 invariant mass; Invariant mass (MeV/c^{2}); pT (GeV/c)",200,0.4,0.6,100,0.,10.);
1101 hK0MassInv->SetMinimum(0);
1102 fOutputStudy->Add(hK0MassInv);
1103 }
1104
1105 //for MC analysis only
1106 if (fReadMC) {
1107
1108 for(Int_t i=0;i<fNPtBinsCorr;i++){
1109
1110 //displacement histos
f80e7bba 1111 namePlot="histDispl_K0_Bin"; namePlot+=i;
a2ad7da1 1112 TH1F *hDisplK = new TH1F(namePlot.Data(), "Kaons Displacement; DCA",150,0.,0.15);
1113 hDisplK->SetMinimum(0);
1114 fOutputStudy->Add(hDisplK);
1115
f80e7bba 1116 namePlot="histDispl_K0_HF_Bin"; namePlot+=i;
a2ad7da1 1117 TH1F *hDisplK_HF = new TH1F(namePlot.Data(), "Kaons Displacement (from HF decay only); DCA",150,0.,0.15);
1118 hDisplK_HF->SetMinimum(0);
1119 fOutputStudy->Add(hDisplK_HF);
1120
1121 namePlot="histDispl_Kcharg_Bin"; namePlot+=i;
1122 TH1F *hDisplHadr = new TH1F(namePlot.Data(), "Hadrons Displacement; DCA",150,0.,0.15);
1123 hDisplHadr->SetMinimum(0);
1124 fOutputStudy->Add(hDisplHadr);
1125
1126 namePlot="histDispl_Kcharg_HF_Bin"; namePlot+=i;
1127 TH1F *hDisplHadr_HF = new TH1F(namePlot.Data(), "Hadrons Displacement (from HF decay only); DCA",150,0.,0.15);
1128 hDisplHadr_HF->SetMinimum(0);
1129 fOutputStudy->Add(hDisplHadr_HF);
1130
1131 namePlot="histDispl_Charg_Bin"; namePlot+=i;
1132 TH1F *hDisplCharg = new TH1F(namePlot.Data(), "Charged tracks Displacement; DCA",150,0.,0.15);
1133 hDisplCharg->SetMinimum(0);
1134 fOutputStudy->Add(hDisplCharg);
1135
1136 namePlot="histDispl_Charg_HF_Bin"; namePlot+=i;
1137 TH1F *hDisplCharg_HF = new TH1F(namePlot.Data(), "Charged tracks Displacement (from HF decay only); DCA",150,0.,0.15);
1138 hDisplCharg_HF->SetMinimum(0);
1139 fOutputStudy->Add(hDisplCharg_HF);
1140
f80e7bba 1141 namePlot="histDispl_K0_From_c_Bin"; namePlot+=i;
a2ad7da1 1142 TH1F *hDisplK_c = new TH1F(namePlot.Data(), "Kaons Displacement - c origin; DCA",150,0.,0.15);
1143 hDisplK_c->SetMinimum(0);
1144 fOutputStudy->Add(hDisplK_c);
1145
f80e7bba 1146 namePlot="histDispl_K0_HF_From_c_Bin"; namePlot+=i;
a2ad7da1 1147 TH1F *hDisplK_HF_c = new TH1F(namePlot.Data(), "Kaons Displacement (from HF decay only) - c origin; DCA",150,0.,0.15);
1148 hDisplK_HF_c->SetMinimum(0);
1149 fOutputStudy->Add(hDisplK_HF_c);
1150
1151 namePlot="histDispl_Kcharg_From_c_Bin"; namePlot+=i;
1152 TH1F *hDisplHadr_c = new TH1F(namePlot.Data(), "Hadrons Displacement - c origin; DCA",150,0.,0.15);
1153 hDisplHadr_c->SetMinimum(0);
1154 fOutputStudy->Add(hDisplHadr_c);
1155
1156 namePlot="histDispl_Kcharg_HF_From_c_Bin"; namePlot+=i;
1157 TH1F *hDisplHadr_HF_c = new TH1F(namePlot.Data(), "Hadrons Displacement (from HF decay only) - c origin; DCA",150,0.,0.15);
1158 hDisplHadr_HF_c->SetMinimum(0);
1159 fOutputStudy->Add(hDisplHadr_HF_c);
1160
1161 namePlot="histDispl_Charg_From_c_Bin"; namePlot+=i;
1162 TH1F *hDisplCharg_c = new TH1F(namePlot.Data(), "Charged tracks Displacement - c origin; DCA",150,0.,0.15);
1163 hDisplCharg_c->Sumw2();
1164 hDisplCharg_c->SetMinimum(0);
1165 fOutputStudy->Add(hDisplCharg_c);
1166
1167 namePlot="histDispl_Charg_HF_From_c_Bin"; namePlot+=i;
1168 TH1F *hDisplCharg_HF_c = new TH1F(namePlot.Data(), "Charged tracks Displacement (from HF decay only) - c origin; DCA",150,0.,0.15);
1169 hDisplCharg_HF_c->SetMinimum(0);
1170 fOutputStudy->Add(hDisplCharg_HF_c);
1171
f80e7bba 1172 namePlot="histDispl_K0_From_b_Bin"; namePlot+=i;
a2ad7da1 1173 TH1F *hDisplK_b = new TH1F(namePlot.Data(), "Kaons Displacement - b origin; DCA",150,0.,0.15);
1174 hDisplK_b->SetMinimum(0);
1175 fOutputStudy->Add(hDisplK_b);
1176
f80e7bba 1177 namePlot="histDispl_K0_HF_From_b_Bin"; namePlot+=i;
a2ad7da1 1178 TH1F *hDisplK_HF_b = new TH1F(namePlot.Data(), "Kaons Displacement (from HF decay only) - b origin; DCA",150,0.,0.15);
1179 hDisplK_HF_b->SetMinimum(0);
1180 fOutputStudy->Add(hDisplK_HF_b);
1181
1182 namePlot="histDispl_Kcharg_From_b_Bin"; namePlot+=i;
1183 TH1F *hDisplHadr_b = new TH1F(namePlot.Data(), "Hadrons Displacement - b origin; DCA",150,0.,0.15);
1184 hDisplHadr_b->SetMinimum(0);
1185 fOutputStudy->Add(hDisplHadr_b);
1186
1187 namePlot="histDispl_Kcharg_HF_From_b_Bin"; namePlot+=i;
1188 TH1F *hDisplHadr_HF_b = new TH1F(namePlot.Data(), "Hadrons Displacement (from HF decay only) - b origin; DCA",150,0.,0.15);
1189 hDisplHadr_HF_b->SetMinimum(0);
1190 fOutputStudy->Add(hDisplHadr_HF_b);
1191
1192 namePlot="histDispl_Charg_From_b_Bin"; namePlot+=i;
1193 TH1F *hDisplCharg_b = new TH1F(namePlot.Data(), "Charged tracks Displacement - b origin; DCA",150,0.,0.15);
1194 hDisplCharg_b->SetMinimum(0);
1195 fOutputStudy->Add(hDisplCharg_b);
1196
1197 namePlot="histDispl_Charg_HF_From_b_Bin"; namePlot+=i;
1198 TH1F *hDisplCharg_HF_b = new TH1F(namePlot.Data(), "Charged tracks Displacement (from HF decay only) - b origin; DCA",150,0.,0.15);
1199 hDisplCharg_HF_b->SetMinimum(0);
1200 fOutputStudy->Add(hDisplCharg_HF_b);
1201
1202 //origin of tracks histos
1203 namePlot="histOrig_Charg_Bin"; namePlot+=i;
1204 TH1F *hOrigin_Charm = new TH1F(namePlot.Data(), "Origin of charged tracks",9,0.,9.);
1205 hOrigin_Charm->SetMinimum(0);
1206 hOrigin_Charm->GetXaxis()->SetBinLabel(1,"Not HF");
1207 hOrigin_Charm->GetXaxis()->SetBinLabel(2,"D->#");
1208 hOrigin_Charm->GetXaxis()->SetBinLabel(3,"D->X->#");
1209 hOrigin_Charm->GetXaxis()->SetBinLabel(4,"B->#");
1210 hOrigin_Charm->GetXaxis()->SetBinLabel(5,"B->X-># (X!=D)");
1211 hOrigin_Charm->GetXaxis()->SetBinLabel(6,"B->D->#");
1212 hOrigin_Charm->GetXaxis()->SetBinLabel(7,"B->D->X->#");
1213 hOrigin_Charm->GetXaxis()->SetBinLabel(8,"c hadr.");
1214 hOrigin_Charm->GetXaxis()->SetBinLabel(9,"b hadr.");
1215 fOutputStudy->Add(hOrigin_Charm);
1216
1217 namePlot="histOrig_Kcharg_Bin"; namePlot+=i;
1218 TH1F *hOrigin_Kcharg = new TH1F(namePlot.Data(), "Origin of hadrons",9,0.,9.);
1219 hOrigin_Kcharg->SetMinimum(0);
1220 hOrigin_Kcharg->GetXaxis()->SetBinLabel(1,"Not HF");
1221 hOrigin_Kcharg->GetXaxis()->SetBinLabel(2,"D->#");
1222 hOrigin_Kcharg->GetXaxis()->SetBinLabel(3,"D->X->#");
1223 hOrigin_Kcharg->GetXaxis()->SetBinLabel(4,"B->#");
1224 hOrigin_Kcharg->GetXaxis()->SetBinLabel(5,"B->X-># (X!=D)");
1225 hOrigin_Kcharg->GetXaxis()->SetBinLabel(6,"B->D->#");
1226 hOrigin_Kcharg->GetXaxis()->SetBinLabel(7,"B->D->X->#");
1227 hOrigin_Kcharg->GetXaxis()->SetBinLabel(8,"c hadr.");
1228 hOrigin_Kcharg->GetXaxis()->SetBinLabel(9,"b hadr.");
1229 fOutputStudy->Add(hOrigin_Kcharg);
1230
f80e7bba 1231 namePlot="histOrig_K0_Bin"; namePlot+=i;
a2ad7da1 1232 TH1F *hOrigin_K = new TH1F(namePlot.Data(), "Origin of kaons",9,0.,9.);
1233 hOrigin_K->SetMinimum(0);
1234 hOrigin_K->GetXaxis()->SetBinLabel(1,"Not HF");
1235 hOrigin_K->GetXaxis()->SetBinLabel(2,"D->#");
1236 hOrigin_K->GetXaxis()->SetBinLabel(3,"D->X->#");
1237 hOrigin_K->GetXaxis()->SetBinLabel(4,"B->#");
1238 hOrigin_K->GetXaxis()->SetBinLabel(5,"B->X-># (X!=D)");
1239 hOrigin_K->GetXaxis()->SetBinLabel(6,"B->D->#");
1240 hOrigin_K->GetXaxis()->SetBinLabel(7,"B->D->X->#");
1241 hOrigin_K->GetXaxis()->SetBinLabel(8,"c hadr.");
1242 hOrigin_K->GetXaxis()->SetBinLabel(9,"b hadr.");
1243 fOutputStudy->Add(hOrigin_K);
1244
1245 //origin of D0 histos
1246 namePlot="histOrig_D0_Bin"; namePlot+=i;
1247 TH1F *hOrigin_D0 = new TH1F(namePlot.Data(), "Origin of D0",2,0.,2.);
1248 hOrigin_D0->SetMinimum(0);
1249 hOrigin_D0->GetXaxis()->SetBinLabel(1,"From c");
1250 hOrigin_D0->GetXaxis()->SetBinLabel(2,"From b");
1251 fOutputStudy->Add(hOrigin_D0);
1252 }
1253 }
1254
bce70c96 1255 if(fMixing) { //check for event mixing (events in bins)
1256 TH2I* EventMixingCheckTr = new TH2I("EventMixingCheckTr","EventMixingCheckTr",5,-0.5,4.5,7,-0.5,6.5);
1257 fOutputStudy->Add(EventMixingCheckTr);
1258
1259 TH2I* EventMixingCheckK0 = new TH2I("EventMixingCheckK0","EventMixingCheckK0",5,-0.5,4.5,7,-0.5,6.5);
1260 fOutputStudy->Add(EventMixingCheckK0);
1261 }
1262
a2ad7da1 1263}
1264
1265//________________________________________________________________________
bce70c96 1266void AliAnalysisTaskSED0Correlations::CalculateCorrelations(AliAODRecoDecayHF2Prong* d, Int_t labD0, TClonesArray* mcArray) {
a2ad7da1 1267//
1268// Method for correlations D0-hadrons study
1269//
bce70c96 1270 Int_t N_Charg = 0, N_KCharg = 0, N_Kaons = 0;
1271 Double_t mD0, mD0bar;
1272 Int_t origD0 = 0, PDGD0 = 0;
a2ad7da1 1273 d->InvMassD0(mD0,mD0bar);
1274 Int_t ptbin = PtBinCorr(d->Pt());
1275 if(ptbin < 0) return;
a2ad7da1 1276
bce70c96 1277 //Origin of D0
1278 TString orig="";
1279 if(fReadMC) {
1280 origD0=CheckD0Origin(mcArray,(AliAODMCParticle*)mcArray->At(labD0));
1281 PDGD0 = ((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode();
1282 switch (CheckD0Origin(mcArray,(AliAODMCParticle*)mcArray->At(labD0))) {
1283 case 4:
1284 orig = "_From_c";
1285 ((TH1F*)fOutputStudy->FindObject(Form("histOrig_D0_Bin%d",ptbin)))->Fill(0.);
1286 break;
1287 case 5:
1288 orig = "_From_b";
1289 ((TH1F*)fOutputStudy->FindObject(Form("histOrig_D0_Bin%d",ptbin)))->Fill(1.);
1290 break;
1291 default:
1292 return;
1293 }
1294 }
1295
1296 Double_t highPt = 0; Double_t lead[3] = {0,0,0}; //infos for leading particle (pt,deltaphi)
a2ad7da1 1297
bce70c96 1298 //loop over the tracks in the pool
1299 Bool_t execPool = fCorrelatorTr->ProcessEventPool(); //pool is ready? (only in ME, in SE returns kFALSE)
1300
1301 if(fMixing && !execPool) {
1302 AliInfo("Mixed event analysis: pool is not ready");
1303 return;
1304 }
1305
1306 Int_t NofEventsinPool = 1;
1307 if(fMixing) NofEventsinPool = fCorrelatorTr->GetNofEventsInPool();
1308
1309 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)
1310 Bool_t analyzetracksTr = fCorrelatorTr->ProcessAssociatedTracks(jMix);// process all the tracks in the aodEvent, by applying the selection cuts
1311 Bool_t analyzetracksKc = fCorrelatorKc->ProcessAssociatedTracks(jMix);
1312 Bool_t analyzetracksK0 = fCorrelatorK0->ProcessAssociatedTracks(jMix);
1313 if(!analyzetracksTr || !analyzetracksKc || !analyzetracksK0) {
1314 AliInfo("AliHFCorrelator::Cannot process the track array");
1315 continue;
1316 }
a2ad7da1 1317
bce70c96 1318 //Charged tracks
1319 for(Int_t iTrack = 0; iTrack<fCorrelatorTr->GetNofTracks(); iTrack++){ // looping on track candidates
1320
1321 Bool_t runcorrelation = fCorrelatorTr->Correlate(iTrack);
1322 if(!runcorrelation) continue;
1323
1324 AliReducedParticle* track = fCorrelatorTr->GetAssociatedParticle();
1325
1326 if(!fMixing) {
1327 if(!track->CheckSoftPi()) { //removal of soft pions
1328 if (fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(1.,mD0);
1329 if (fIsSelectedCandidate >= 2) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(1.,mD0bar);
1330 continue;
1331 } else { //not a soft pion
1332 if (fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(0.,mD0);
1333 if (fIsSelectedCandidate >= 2) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(0.,mD0bar);
1334 }
a2ad7da1 1335 }
bce70c96 1336 if(track->Pt() < fPtThreshLow.at(ptbin) || track->Pt() > fPtThreshUp.at(ptbin)) continue; //discard tracks outside pt range for hadrons/K
1337
1338 Int_t idDaughs[2] = {((AliVTrack*)d->GetDaughter(0))->GetID(),((AliVTrack*)d->GetDaughter(1))->GetID()}; //IDs of daughters to be skipped
1339 if(track->GetID() == idDaughs[0] || track->GetID() == idDaughs[1]) continue; //discards daughters of candidate
1340
1341 FillSparsePlots(mcArray,d,origD0,PDGD0,track,kTrack); //fills for charged tracks
1342
1343 if(!fMixing) N_Charg++;
1344
1345 //retrieving leading info...
1346 if(track->Pt() > highPt) {
1347 lead[0] = fCorrelatorTr->GetDeltaPhi();
1348 lead[1] = fCorrelatorTr->GetDeltaEta();
1349 lead[2] = fReadMC ? CheckTrackOrigin(mcArray,(AliAODMCParticle*)mcArray->At(track->GetLabel())) : 0;
1350 highPt = track->Pt();
a2ad7da1 1351 }
bce70c96 1352
1353 } // end of tracks loop
1354
1355 //Charged Kaons loop
1356 for(Int_t iTrack = 0; iTrack<fCorrelatorKc->GetNofTracks(); iTrack++){ // looping on charged kaons candidates
1357
1358 Bool_t runcorrelation = fCorrelatorKc->Correlate(iTrack);
1359 if(!runcorrelation) continue;
1360
1361 AliReducedParticle* kCharg = fCorrelatorKc->GetAssociatedParticle();
1362
1363 if(!fMixing) {
1364 if(kCharg->CheckSoftPi()) { //removal of soft pions
1365 if (fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(1.,mD0);
1366 if (fIsSelectedCandidate >= 2) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(1.,mD0bar);
1367 continue;
1368 } else {
1369 if (fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(0.,mD0);
1370 if (fIsSelectedCandidate >= 2) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(0.,mD0bar);
a2ad7da1 1371 }
bce70c96 1372 }
1373 if(kCharg->Pt() < fPtThreshLow.at(ptbin) || kCharg->Pt() > fPtThreshUp.at(ptbin)) continue; //discard tracks outside pt range for hadrons/K
a2ad7da1 1374
bce70c96 1375 Int_t idDaughs[2] = {((AliVTrack*)d->GetDaughter(0))->GetID(),((AliVTrack*)d->GetDaughter(1))->GetID()}; //IDs of daughters to be skipped
1376 if(kCharg->GetID() == idDaughs[0] || kCharg->GetID() == idDaughs[1]) continue; //discards daughters of candidate
1377
1378 FillSparsePlots(mcArray,d,origD0,PDGD0,kCharg,kKCharg); //fills for charged tracks
a2ad7da1 1379
bce70c96 1380 if(!fMixing) N_KCharg++;
a2ad7da1 1381
bce70c96 1382 } // end of charged kaons loop
a2ad7da1 1383
bce70c96 1384 //K0 loop
1385 for(Int_t iTrack = 0; iTrack<fCorrelatorK0->GetNofTracks(); iTrack++){ // looping on k0 candidates
1386
1387 Bool_t runcorrelation = fCorrelatorK0->Correlate(iTrack);
1388 if(!runcorrelation) continue;
1389
1390 AliReducedParticle* k0 = fCorrelatorK0->GetAssociatedParticle();
a2ad7da1 1391
bce70c96 1392 if(k0->Pt() < fPtThreshLow.at(ptbin) || k0->Pt() > fPtThreshUp.at(ptbin)) continue; //discard tracks outside pt range for hadrons/K
1393
1394 Int_t idDaughs[2] = {((AliVTrack*)d->GetDaughter(0))->GetID(),((AliVTrack*)d->GetDaughter(1))->GetID()}; //IDs of daughters to be skipped
1395 if(k0->GetID() == idDaughs[0] || k0->GetID() == idDaughs[1]) continue; //discards daughters of candidate
1396
1397 FillSparsePlots(mcArray,d,origD0,PDGD0,k0,kK0); //fills for charged tracks
1398
1399 if(!fMixing) N_Kaons++;
1400
1401 } // end of charged kaons loop
1402
1403 } //end of events loop
1404
1405 //leading track correlations fill
1406 if(fReadMC) {
1407 if(((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode()==421) { //D0
1408 ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(lead[1],mD0); //c and b D0
1409 ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead%s_Bin%d",orig.Data(),ptbin)))->Fill(lead[0],mD0); //c or b D0
1410 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);
1411 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);
a2ad7da1 1412 }
bce70c96 1413 if(((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode()==-421) { //D0bar
1414 ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(lead[1],mD0bar);
1415 ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead%s_Bin%d",orig.Data(),ptbin)))->Fill(lead[0],mD0bar); //c or b D0
1416 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);
1417 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);
a2ad7da1 1418 }
bce70c96 1419 } else {
1420 if(fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(lead[0],mD0); //c and b D0
1421 if(fIsSelectedCandidate == 2 || fIsSelectedCandidate == 3) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(lead[0],mD0bar);
1422 }
a2ad7da1 1423
bce70c96 1424 //Fill of count histograms
1425 if (!fAlreadyFilled) {
1426 ((TH1F*)fOutputStudy->FindObject("hist_Count_Charg"))->Fill(N_Charg);
1427 ((TH1F*)fOutputStudy->FindObject("hist_Count_Kcharg"))->Fill(N_KCharg);
1428 ((TH1F*)fOutputStudy->FindObject("hist_Count_K0"))->Fill(N_Kaons);
1429 }
a2ad7da1 1430
bce70c96 1431 fAlreadyFilled=kTRUE; //at least a D0 analyzed in the event; distribution plots already filled
a2ad7da1 1432
bce70c96 1433}
1434
1435//________________________________________________________________________
1436void AliAnalysisTaskSED0Correlations::FillSparsePlots(TClonesArray* mcArray, AliAODRecoDecayHF2Prong *d, Int_t origD0, Int_t PdgD0, AliReducedParticle* track, Int_t type) {
1437 //
1438 //fills the THnSparse for correlations, calculating the variables
1439 //
1440
1441 //Initialization of variables
1442 Int_t ptbin = PtBinCorr(d->Pt());
1443 if(ptbin < 0) return;
1444 Double_t mD0, mD0bar, deltaphi = 0., deltaeta = 0.;
1445 d->InvMassD0(mD0,mD0bar);
1446
1447 Double_t ptTrack = track->Pt();
1448 Double_t d0Track = type!=kK0 ? track->GetImpPar() : 0.;
1449 Double_t phiTr = track->Phi();
1450 Double_t origTr = fReadMC ? CheckTrackOrigin(mcArray,(AliAODMCParticle*)mcArray->At(track->GetLabel())) : 0;
1451
1452 TString part = "", orig = "";
1453
1454 switch (type) {
1455 case(kTrack): {
1456 part = "Charg";
1457 deltaphi = fCorrelatorTr->GetDeltaPhi();
1458 deltaeta = fCorrelatorTr->GetDeltaEta();
1459 break;
1460 }
1461 case(kKCharg): {
1462 part = "Kcharg";
1463 deltaphi = fCorrelatorKc->GetDeltaPhi();
1464 deltaeta = fCorrelatorKc->GetDeltaEta();
1465 break;
1466 }
1467 case(kK0): {
1468 part = "K0";
1469 deltaphi = fCorrelatorK0->GetDeltaPhi();
1470 deltaeta = fCorrelatorK0->GetDeltaEta();
1471 break;
a2ad7da1 1472 }
bce70c96 1473 }
a2ad7da1 1474
bce70c96 1475 //Fixes limits; needed to include overflow into THnSparse projections!
1476 Double_t pTorig = track->Pt();
1477 Double_t d0orig = track->GetImpPar();
1478 Double_t ptLim_Sparse = ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->GetAxis(2)->GetXmax(); //all plots have same axes...
1479 Double_t displLim_Sparse = ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->GetAxis(3)->GetXmax();
1480 Double_t EtaLim_Sparse = ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->GetAxis(4)->GetXmax();
1481 if(ptTrack > ptLim_Sparse) ptTrack = ptLim_Sparse-0.01;
1482 if(d0Track > displLim_Sparse) d0Track = (displLim_Sparse-0.001);
1483 if(deltaeta > EtaLim_Sparse) deltaeta = EtaLim_Sparse-0.01;
1484 if(deltaeta < -EtaLim_Sparse) deltaeta = -EtaLim_Sparse+0.01;
1485
1486 if(fMixing == kSE) {
1487
1488 //variables for filling histos
1489 Double_t fillSpPhiD0[5] = {deltaphi,mD0,ptTrack,d0Track,deltaeta};
1490 Double_t fillSpPhiD0bar[5] = {deltaphi,mD0bar,ptTrack,d0Track,deltaeta};
a2ad7da1 1491
bce70c96 1492 if(fReadMC == 0) {
1493 //sparse fill for data (tracks, K+-, K0) + weighted
1494 if(fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) { //D0
1495 ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0);
1496 ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_Bin%d",ptbin)))->Fill(deltaphi,mD0,pTorig);
a2ad7da1 1497 }
bce70c96 1498 if(fIsSelectedCandidate == 2 || fIsSelectedCandidate == 3) { //D0bar
1499 ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0bar);
1500 ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_Bin%d",ptbin)))->Fill(deltaphi,mD0bar,pTorig);
a2ad7da1 1501 }
1502 if(!fAlreadyFilled) {
bce70c96 1503 ((TH1F*)fOutputStudy->FindObject(Form("hist_Pt_%s_Bin%d",part.Data(),ptbin)))->Fill(pTorig);
1504 ((TH1F*)fOutputStudy->FindObject(Form("hist_PhiDistr_%s",part.Data())))->Fill(phiTr);
a2ad7da1 1505 }
bce70c96 1506 }
a2ad7da1 1507
bce70c96 1508 if(fReadMC) {
a2ad7da1 1509
bce70c96 1510 if(origD0==4) {orig = "_From_c";} else {orig = "_From_b";}
a2ad7da1 1511
bce70c96 1512 //sparse fill for data (tracks, K+-, K0) + weighted
1513 if(PdgD0==421) { //D0 (from MCTruth)
1514 ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0);
1515 ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0);
1516 if(origD0==4&&origTr>=1&&origTr<=2) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_HF%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0);
1517 if(origD0==5&&origTr>=3&&origTr<=6) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_HF%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0);
1518 ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_Bin%d",ptbin)))->Fill(deltaphi,mD0,pTorig);
1519 ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0,pTorig);
1520 if(origD0==4&&origTr>=1&&origTr<=2) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0,pTorig);
1521 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 1522 }
bce70c96 1523 if(PdgD0==-421) { //D0bar
1524 ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0bar);
1525 ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0bar);
1526 if(origD0==4&&origTr>=1&&origTr<=2) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_HF%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0bar);
1527 if(origD0==5&&origTr>=3&&origTr<=6) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_HF%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0bar);
1528 ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_Bin%d",ptbin)))->Fill(deltaphi,mD0bar,pTorig);
1529 ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0bar,pTorig);
1530 if(origD0==4&&origTr>=1&&origTr<=2) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0bar,pTorig);
1531 if(origD0==5&&origTr>=3&&origTr<=6) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0bar,pTorig);
1532 }
a2ad7da1 1533 if(!fAlreadyFilled) {
bce70c96 1534 ((TH1F*)fOutputStudy->FindObject(Form("histDispl_%s_Bin%d",part.Data(),ptbin)))->Fill(d0orig); //Fills displacement histos
1535 if (origTr>=1&&origTr<=6) ((TH1F*)fOutputStudy->FindObject(Form("histDispl_%s_HF_Bin%d",part.Data(),ptbin)))->Fill(d0orig);
1536 if (origTr>=1&&origTr<=6) ((TH1F*)fOutputStudy->FindObject(Form("histDispl_%s_HF%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(d0orig);
1537 ((TH1F*)fOutputStudy->FindObject(Form("histDispl_%s%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(d0orig); //Fills displacement histos
1538 ((TH1F*)fOutputStudy->FindObject(Form("hist_Pt_%s_Bin%d",part.Data(),ptbin)))->Fill(pTorig);
1539 ((TH1F*)fOutputStudy->FindObject(Form("histOrig_%s_Bin%d",part.Data(),ptbin)))->Fill(origTr);
1540 ((TH1F*)fOutputStudy->FindObject(Form("hist_PhiDistr_%s",part.Data())))->Fill(phiTr);
a2ad7da1 1541 }
bce70c96 1542 }//end MC case
a2ad7da1 1543
bce70c96 1544 } //end of SE fill
a2ad7da1 1545
bce70c96 1546 if(fMixing == kME) {
a2ad7da1 1547
bce70c96 1548 //variables for filling histos
1549 Double_t fillSpPhiD0[3] = {deltaphi,mD0,deltaeta};
1550 Double_t fillSpPhiD0bar[3] = {deltaphi,mD0bar,deltaeta};
a2ad7da1 1551
bce70c96 1552 if(fReadMC == 0) {
1553 //sparse fill for data (tracks, K+-, K0)
1554 if(fIsSelectedCandidate == 1||fIsSelectedCandidate == 3) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d_EvMix",part.Data(),ptbin)))->Fill(fillSpPhiD0);
1555 if(fIsSelectedCandidate == 2||fIsSelectedCandidate == 3) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d_EvMix",part.Data(),ptbin)))->Fill(fillSpPhiD0bar);
1556 }
1557 if(fReadMC == 1) {
1558 //sparse fill for data (tracks, K+-, K0)
1559 if(PdgD0==421) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d_EvMix",part.Data(),ptbin)))->Fill(fillSpPhiD0);
1560 if(PdgD0==-421) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d_EvMix",part.Data(),ptbin)))->Fill(fillSpPhiD0bar);
1561 }//end MC case
a2ad7da1 1562
bce70c96 1563 } //end of ME fill
1564
1565 return;
a2ad7da1 1566}
1567
1568//_________________________________________________________________________________________________
1569Int_t AliAnalysisTaskSED0Correlations::CheckTrackOrigin(TClonesArray* arrayMC, AliAODMCParticle *mcPartCandidate) const {
1570 //
1571 // checks on particle (#) origin:
1572 // 0) Not HF
1573 // 1) D->#
1574 // 2) D->X->#
1575 // 3) B->#
1576 // 4) B->X-># (X!=D)
1577 // 5) B->D->#
1578 // 6) B->D->X->#
1579 // 7) c hadronization
1580 // 8) b hadronization
1581 //
bce70c96 1582 if(fDebug>2) printf("AliAnalysisTaskSED0Correlations::CheckTrkOrigin() \n");
a2ad7da1 1583
1584 Int_t pdgGranma = 0;
1585 Int_t mother = 0;
1586 mother = mcPartCandidate->GetMother();
1587 Int_t istep = 0;
1588 Int_t abspdgGranma =0;
1589 Bool_t isFromB=kFALSE;
1590 Bool_t isDdaugh=kFALSE;
1591 Bool_t isDchaindaugh=kFALSE;
1592 Bool_t isBdaugh=kFALSE;
1593 Bool_t isBchaindaugh=kFALSE;
1594 Bool_t isQuarkFound=kFALSE;
1595
1596 while (mother > 0){
1597 istep++;
1598 AliAODMCParticle* mcGranma = dynamic_cast<AliAODMCParticle*>(arrayMC->At(mother));
1599 if (mcGranma){
1600 pdgGranma = mcGranma->GetPdgCode();
1601 abspdgGranma = TMath::Abs(pdgGranma);
1602 if ((abspdgGranma > 500 && abspdgGranma < 600) || (abspdgGranma > 5000 && abspdgGranma < 6000)){
1603 isBchaindaugh=kTRUE;
1604 if(istep==1) isBdaugh=kTRUE;
1605 }
1606 if ((abspdgGranma > 400 && abspdgGranma < 500) || (abspdgGranma > 4000 && abspdgGranma < 5000)){
1607 isDchaindaugh=kTRUE;
1608 if(istep==1) isDdaugh=kTRUE;
1609 }
1610 if(abspdgGranma==4 || abspdgGranma==5) {isQuarkFound=kTRUE; if(abspdgGranma==5) isFromB = kTRUE;}
1611 mother = mcGranma->GetMother();
1612 }else{
1613 AliError("Failed casting the mother particle!");
1614 break;
1615 }
1616 }
1617
1618 //decides what to return based on the flag status
1619 if(isQuarkFound) {
1620 if(!isFromB) { //charm
1621 if(isDdaugh) return 1; //charm immediate
1622 else if(isDchaindaugh) return 2; //charm chain
1623 else return 7; //charm hadronization
1624 }
1625 else { //beauty
1626 if(isBdaugh) return 3; //b immediate
1627 else if(isBchaindaugh) { //b chain
1628 if(isDchaindaugh) {
1629 if(isDdaugh) return 5; //d immediate
1630 return 6; //d chain
1631 }
1632 else return 4; //b, not d
1633 }
1634 else return 8; //b hadronization
1635 }
1636 }
1637 else return 0; //no HF quark
1638}
1639
a2ad7da1 1640//________________________________________________________________________
1641Int_t AliAnalysisTaskSED0Correlations::PtBinCorr(Double_t pt) const {
1642 //
1643 //give the pt bin where the pt lies.
1644 //
1645 Int_t ptbin=-1;
1646 if(pt<fBinLimsCorr.at(0)) return ptbin; //out of bounds
1647
1648 Int_t i = 0;
1649 while(pt>fBinLimsCorr.at(i)) {ptbin=i; i++;}
1650
1651 return ptbin;
1652}
1653
a2ad7da1 1654//---------------------------------------------------------------------------
1655Bool_t AliAnalysisTaskSED0Correlations::SelectV0(AliAODv0* v0, AliAODVertex *vtx, Int_t opt, Int_t idArrayV0[][2]) const
1656{
1657 //
1658 // Selection for K0 hypotheses
1659 // options: 1 = selects mass invariant about 3 sigma inside the peak + threshold of 0.3 GeV
1660 // 2 = no previous selections
1661
1662 if(!fCutsTracks->IsKZeroSelected(v0,vtx)) return kFALSE;
1663
1664 AliAODTrack *v0Daug1 = (AliAODTrack*)v0->GetDaughter(0);
1665 AliAODTrack *v0Daug2 = (AliAODTrack*)v0->GetDaughter(1);
1666
1667 if(opt==1) { //additional cuts for correlations (V0 has to be closer than 3 sigma from K0 mass)
bce70c96 1668 if(TMath::Abs(v0->MassK0Short()-0.4976) > 3*0.004) return kFALSE;
a2ad7da1 1669 }
1670
1671 //This part removes double counting for swapped tracks!
1672 Int_t i = 0; //while loop (until the last-written entry pair of ID!
1673 while(idArrayV0[i][0]!=-2 && idArrayV0[i][1]!=-2) {
1674 if((v0Daug1->GetID()==idArrayV0[i][0] && v0Daug2->GetID()==idArrayV0[i][1])||
1675 (v0Daug1->GetID()==idArrayV0[i][1] && v0Daug2->GetID()==idArrayV0[i][0])) return kFALSE;
1676 i++;
1677 }
1678 idArrayV0[i][0]=v0Daug1->GetID();
1679 idArrayV0[i][1]=v0Daug2->GetID();
1680
1681 return kTRUE;
1682}
1683
1684//________________________________________________________________________
1685void AliAnalysisTaskSED0Correlations::PrintBinsAndLimits() {
1686
1687 cout << "--------------------------\n";
1688 cout << "PtBins = " << fNPtBinsCorr << "\n";
1689 cout << "PtBin limits--------------\n";
1690 for (int i=0; i<fNPtBinsCorr; i++) {
1691 cout << "Bin "<<i+1<<" = "<<fBinLimsCorr.at(i)<<" to "<<fBinLimsCorr.at(i+1)<<"\n";
1692 }
1693 cout << "\n--------------------------\n";
1694 cout << "PtBin tresh. tracks low---\n";
1695 for (int i=0; i<fNPtBinsCorr; i++) {
1696 cout << fPtThreshLow.at(i) << "; ";
1697 }
1698 cout << "PtBin tresh. tracks up----\n";
1699 for (int i=0; i<fNPtBinsCorr; i++) {
1700 cout << fPtThreshUp.at(i) << "; ";
1701 }
1702 cout << "\n--------------------------\n";
1703 cout << "MC Truth = "<<fReadMC<<"\n";
1704 cout << "--------------------------\n";
bce70c96 1705 cout << "Ev Mixing = "<<fMixing<<"\n";
1706 cout << "--------------------------\n";
a2ad7da1 1707}
1708