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