]> 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
Merge branch 'feature-movesplit'
[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;
4640c275 610 if (aod) ntracks = aod->GetNumberOfTracks();
a2ad7da1 611 for(Int_t itrack=0; itrack<ntracks; itrack++) { // loop on tacks
612 // ... get the track
f15c1f69 613 AliAODTrack * track = dynamic_cast<AliAODTrack*>(aod->GetTrack(itrack));
ad718085 614 if(!track){
615 AliWarning("Error in casting to AOD track. Not a standard AOD?");
616 continue;
617 }
a2ad7da1 618 if(TESTBIT(track->GetITSClusterMap(),2) || TESTBIT(track->GetITSClusterMap(),3) ){
619 skipEvent=kTRUE;
620 fNentries->Fill(16);
621 break;
622 }
623 }
624 if (skipEvent) return;
625 }
626
bce70c96 627 //HFCorrelators initialization (for this event)
628 fCorrelatorTr->SetAODEvent(aod); // set the AOD event from which you are processing
629 fCorrelatorKc->SetAODEvent(aod);
630 fCorrelatorK0->SetAODEvent(aod);
631 Bool_t correlatorONTr = fCorrelatorTr->Initialize(); // initialize the pool for event mixing
632 Bool_t correlatorONKc = fCorrelatorKc->Initialize();
633 Bool_t correlatorONK0 = fCorrelatorK0->Initialize();
634 if(!correlatorONTr) {AliInfo("AliHFCorrelator (tracks) didn't initialize the pool correctly or processed a bad event"); return;}
635 if(!correlatorONKc) {AliInfo("AliHFCorrelator (charged K) didn't initialize the pool correctly or processed a bad event"); return;}
636 if(!correlatorONK0) {AliInfo("AliHFCorrelator (K0) didn't initialize the pool correctly or processed a bad event"); return;}
637
638 if(fReadMC) {
639 fCorrelatorTr->SetMCArray(mcArray); // set the TClonesArray *fmcArray for analysis on monte carlo
640 fCorrelatorKc->SetMCArray(mcArray);
641 fCorrelatorK0->SetMCArray(mcArray);
642 }
643
a2ad7da1 644 // AOD primary vertex
645 AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex();
a2ad7da1 646
647 //vtx1->Print();
648 TString primTitle = vtx1->GetTitle();
649 if(primTitle.Contains("VertexerTracks") && vtx1->GetNContributors()>0) {
a2ad7da1 650 fNentries->Fill(3);
651 }
652
653 //Reset flag for tracks distributions fill
654 fAlreadyFilled=kFALSE;
655
bce70c96 656 //***** Loop over D0 candidates *****
a2ad7da1 657 Int_t nInD0toKpi = inputArray->GetEntriesFast();
658 if(fDebug>2) printf("Number of D0->Kpi: %d\n",nInD0toKpi);
659
660 if(fFillGlobal) { //loop on V0 and tracks for each event, to fill Pt distr. and InvMass distr.
661
662 TClonesArray *v0array = (TClonesArray*)aod->GetList()->FindObject("v0s");
663 Int_t pdgCodes[2] = {211,211};
664 Int_t idArrayV0[v0array->GetEntriesFast()][2];
665 for(int iV0=0; iV0<v0array->GetEntriesFast(); iV0++) {
666 for(int j=0; j<2; j++) {idArrayV0[iV0][j]=-2;}
667 AliAODv0 *v0 = (AliAODv0*)v0array->UncheckedAt(iV0);
668 if(SelectV0(v0,vtx1,2,idArrayV0)) { //option 2 = for mass inv plots only
4d978a4a 669 if(fReadMC && fRecoTr && (v0->MatchToMC(310,mcArray,2,pdgCodes)<0)) continue; //310 = K0s, 311 = K0 generico!!
a2ad7da1 670 ((TH2F*)fOutputStudy->FindObject("hK0MassInv"))->Fill(v0->MassK0Short(),v0->Pt()); //invariant mass plot
f80e7bba 671 ((TH1F*)fOutputStudy->FindObject("hist_Pt_K0_AllEv"))->Fill(v0->Pt()); //pT distribution (in all events), K0 case
a2ad7da1 672 }
673 }
674
4640c275 675 for(Int_t itrack=0; itrack<aod->GetNumberOfTracks(); itrack++) { // loop on tacks
f15c1f69 676 AliAODTrack * track = dynamic_cast<AliAODTrack*>(aod->GetTrack(itrack));
ad718085 677 if(!track){
678 AliWarning("Error in casting to AOD track. Not a standard AOD?");
679 continue;
680 }
a2ad7da1 681 //rejection of tracks
682 if(track->GetID() < 0) continue; //discard negative ID tracks
683 if(track->Pt() < fPtThreshLow.at(0) || track->Pt() > fPtThreshUp.at(0)) continue; //discard tracks outside pt range for hadrons/K
bce70c96 684 if(!fCutsTracks->IsHadronSelected(track) || !fCutsTracks->CheckHadronKinematic(track->Pt(),0.1)) continue; //0.1 = dummy (d0 value, no cut on it for me)
a2ad7da1 685 //pT distribution (in all events), charg and hadr cases
686 ((TH1F*)fOutputStudy->FindObject("hist_Pt_Charg_AllEv"))->Fill(track->Pt());
687 if(fCutsTracks->CheckKaonCompatibility(track,kFALSE,0,2)) ((TH1F*)fOutputStudy->FindObject("hist_Pt_Kcharg_AllEv"))->Fill(track->Pt());
688 }
689
690 } //end of loops for global plot fill
691
692 Int_t nSelectedloose=0,nSelectedtight=0;
4d978a4a 693
241994ab 694 //Fill Event Multiplicity (needed only in Reco)
695 fMultEv = (Double_t)(AliVertexingHFUtils::GetNumberOfTrackletsInEtaRange(aod,-1.,1.));
696
4d978a4a 697 //RecoD0 case ************************************************
698 if(fRecoD0) {
699
700 for (Int_t iD0toKpi = 0; iD0toKpi < nInD0toKpi; iD0toKpi++) {
701 AliAODRecoDecayHF2Prong *d = (AliAODRecoDecayHF2Prong*)inputArray->UncheckedAt(iD0toKpi);
241994ab 702
7f221b36 703 if(d->Pt()<2.) continue; //to save time and merging memory...
704
4d978a4a 705 if(d->GetSelectionMap()) if(!d->HasSelectionBit(AliRDHFCuts::kD0toKpiCuts)){
706 fNentries->Fill(2);
707 continue; //skip the D0 from Dstar
a2ad7da1 708 }
709
4d978a4a 710 if (fCutsD0->IsInFiducialAcceptance(d->Pt(),d->Y(421)) ) {
711 nSelectedloose++;
712 nSelectedtight++;
713 if(fSys==0){
714 if(fCutsD0->IsSelected(d,AliRDHFCuts::kTracks,aod))fNentries->Fill(6);
715 }
716 Int_t ptbin=fCutsD0->PtBin(d->Pt());
717 if(ptbin==-1) {fNentries->Fill(4); continue;} //out of bounds
718
719 fIsSelectedCandidate=fCutsD0->IsSelected(d,AliRDHFCuts::kAll,aod); //D0 selected
720 if(!fIsSelectedCandidate) continue;
721
722 //D0 infos
723 Double_t phiD0 = fCorrelatorTr->SetCorrectPhiRange(d->Phi());
724 phiD0 = fCorrelatorKc->SetCorrectPhiRange(d->Phi()); //bad usage, but returns a Double_t...
725 phiD0 = fCorrelatorK0->SetCorrectPhiRange(d->Phi());
726 fCorrelatorTr->SetTriggerParticleProperties(d->Pt(),phiD0,d->Eta()); // sets the parameters of the trigger particles that are needed
727 fCorrelatorKc->SetTriggerParticleProperties(d->Pt(),phiD0,d->Eta());
728 fCorrelatorK0->SetTriggerParticleProperties(d->Pt(),phiD0,d->Eta());
729 fCorrelatorTr->SetD0Properties(d,fIsSelectedCandidate); //sets special properties for D0
730 fCorrelatorKc->SetD0Properties(d,fIsSelectedCandidate);
731 fCorrelatorK0->SetD0Properties(d,fIsSelectedCandidate);
732
733 if(!fReadMC) {
241994ab 734 if (TMath::Abs(d->Eta())<fEtaForCorrel) {
735 CalculateCorrelations(d); //correlations on real data
736 if(!fMixing) ((TH1F*)fOutputStudy->FindObject(Form("hMultiplEvt_Bin%d",ptbin)))->Fill(fMultEv);
737 }
4d978a4a 738 } else { //correlations on MC -> association of selected D0 to MCinfo with MCtruth
739 if (TMath::Abs(d->Eta())<fEtaForCorrel) {
740 Int_t pdgDgD0toKpi[2]={321,211};
741 Int_t labD0 = d->MatchToMC(421,mcArray,2,pdgDgD0toKpi); //return MC particle label if the array corresponds to a D0, -1 if not
241994ab 742 if (labD0>-1) {
743 CalculateCorrelations(d,labD0,mcArray);
744 if(!fMixing) ((TH1F*)fOutputStudy->FindObject(Form("hMultiplEvt_Bin%d",ptbin)))->Fill(fMultEv); //Fill multiplicity histo
745 }
4d978a4a 746 }
747 }
748
749 FillMassHists(d,mcArray,fCutsD0,fOutputMass);
a2ad7da1 750 }
4d978a4a 751 }
752 }
753 //End RecoD0 case ************************************************
754
755 //MCKineD0 case ************************************************
756 if(fReadMC && !fRecoD0) {
757
758 for (Int_t iPart=0; iPart<mcArray->GetEntriesFast(); iPart++) { //Loop over all the tracks of MCArray
759 AliAODMCParticle* mcPart = dynamic_cast<AliAODMCParticle*>(mcArray->At(iPart));
760 if (!mcPart) {
761 AliWarning("Particle not found in tree, skipping");
762 continue;
763 }
764
765 if(TMath::Abs(mcPart->GetPdgCode()) == 421){ // THIS IS A D0
766 if (fCutsD0->IsInFiducialAcceptance(mcPart->Pt(),mcPart->Y()) ) {
767 nSelectedloose++;
768 nSelectedtight++;
dea90a65 769
770 //Removal of cases in which D0 decay is not in Kpi!
771 if(mcPart->GetNDaughters()!=2) continue;
772 AliAODMCParticle* mcDau1 = dynamic_cast<AliAODMCParticle*>(mcArray->At(mcPart->GetDaughter(0)));
773 AliAODMCParticle* mcDau2 = dynamic_cast<AliAODMCParticle*>(mcArray->At(mcPart->GetDaughter(1)));
774 if(!mcDau1 || !mcDau2) continue;
775 Int_t pdg1 = TMath::Abs(mcDau1->GetPdgCode());
776 Int_t pdg2 = TMath::Abs(mcDau2->GetPdgCode());
777 if(!((pdg1 == 211 && pdg2 == 321) || (pdg2 == 211 && pdg1 == 321))) continue;
778 if(TMath::Abs(mcDau1->Eta())>0.8||TMath::Abs(mcDau2->Eta())>0.8) continue;
779 //Check momentum conservation (to exclude 4-prong decays with tracks outside y=1.5)
780 Double_t p1[3] = {mcDau1->Px(),mcDau1->Py(),mcDau1->Pz()};
781 Double_t p2[3] = {mcDau2->Px(),mcDau2->Py(),mcDau2->Pz()};
782 Double_t pD0[3] = {mcPart->Px(),mcPart->Py(),mcPart->Pz()};
783 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;
784
4d978a4a 785 if(fSys==0) fNentries->Fill(6);
786 Int_t ptbin=fCutsD0->PtBin(mcPart->Pt());
787 if(ptbin==-1) {fNentries->Fill(4); continue;} //out of bounds
788
789 //D0 infos
790 Double_t phiD0 = fCorrelatorTr->SetCorrectPhiRange(mcPart->Phi());
791 phiD0 = fCorrelatorKc->SetCorrectPhiRange(mcPart->Phi()); //bad usage, but returns a Double_t...
792 phiD0 = fCorrelatorK0->SetCorrectPhiRange(mcPart->Phi());
793 fCorrelatorTr->SetTriggerParticleProperties(mcPart->Pt(),phiD0,mcPart->Eta()); // sets the parameters of the trigger particles that are needed
794 fCorrelatorKc->SetTriggerParticleProperties(mcPart->Pt(),phiD0,mcPart->Eta());
795 fCorrelatorK0->SetTriggerParticleProperties(mcPart->Pt(),phiD0,mcPart->Eta());
796 //fCorrelatorTr->SetD0Properties(mcPart,fIsSelectedCandidate); //needed for D* soft pions rejection, useless in MCKine
797 //fCorrelatorKc->SetD0Properties(mcPart,fIsSelectedCandidate);
798 //fCorrelatorK0->SetD0Properties(mcPart,fIsSelectedCandidate);
799
800 if (TMath::Abs(mcPart->Eta())<fEtaForCorrel) {
801
802 //Removal of D0 from D* feeddown! This solves also the problem of soft pions, now excluded
803 /* Int_t mother = mcPart->GetMother();
804 AliAODMCParticle* mcMoth = dynamic_cast<AliAODMCParticle*>(mcArray->At(mother));
805 if(!mcMoth) continue;
241994ab 806 if(TMath::Abs(mcMoth->GetPdgCode())==413) continue;*/
807
4d978a4a 808 if (mcPart->GetPdgCode()==421) fIsSelectedCandidate = 1;
809 else fIsSelectedCandidate = 2;
810
241994ab 811 TString fillthis="histSgn_";
812 if(CheckD0Origin(mcArray,mcPart)==4) fillthis+="c_";
813 else if(CheckD0Origin(mcArray,mcPart)==5) fillthis+="b_";
814 else continue;
815 fillthis+=ptbin;
4d978a4a 816 ((TH1F*)(fOutputMass->FindObject(fillthis)))->Fill(1.864);
817
818 CalculateCorrelationsMCKine(mcPart,mcArray);
241994ab 819 if(!fMixing) ((TH1F*)fOutputStudy->FindObject(Form("hMultiplEvt_Bin%d",ptbin)))->Fill(fMultEv); //Fill multiplicity histo
4d978a4a 820 }
bce70c96 821 }
a2ad7da1 822 }
a2ad7da1 823 }
4d978a4a 824 }
825 //End MCKineD0 case ************************************************
bce70c96 826
827 if(fMixing /* && fAlreadyFilled*/) { // update the pool for Event Mixing, if: enabled, event is ok, at least a SelD0 found! (fAlreadyFilled's role!)
828 Bool_t updatedTr = fCorrelatorTr->PoolUpdate();
829 Bool_t updatedKc = fCorrelatorKc->PoolUpdate();
830 Bool_t updatedK0 = fCorrelatorK0->PoolUpdate();
831 if(!updatedTr || !updatedKc || !updatedK0) AliInfo("Pool was not updated");
832 }
833
a2ad7da1 834 fCounter->StoreCandidates(aod,nSelectedloose,kTRUE);
835 fCounter->StoreCandidates(aod,nSelectedtight,kFALSE);
836
837 // Post the data
838 PostData(1,fOutputMass);
839 PostData(2,fNentries);
840 PostData(4,fCounter);
841 PostData(5,fOutputCorr);
842 PostData(6,fOutputStudy);
843
844 return;
845}
846
847//____________________________________________________________________________
848void AliAnalysisTaskSED0Correlations::FillMassHists(AliAODRecoDecayHF2Prong *part, TClonesArray *arrMC, AliRDHFCutsD0toKpi* cuts, TList *listout){
849 //
850 // function used in UserExec to fill mass histograms:
851 //
852 if (!fIsSelectedCandidate) return;
853
854 if(fDebug>2) cout<<"Candidate selected"<<endl;
855
856 Double_t invmassD0 = part->InvMassD0(), invmassD0bar = part->InvMassD0bar();
857 Int_t ptbin = cuts->PtBin(part->Pt());
858
859 TString fillthis="";
860 Int_t pdgDgD0toKpi[2]={321,211};
861 Int_t labD0=-1;
862 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)
863
864 //count candidates selected by cuts
865 fNentries->Fill(1);
866 //count true D0 selected by cuts
867 if (fReadMC && labD0>=0) fNentries->Fill(2);
868
869 if ((fIsSelectedCandidate==1 || fIsSelectedCandidate==3) && fFillOnlyD0D0bar<2) { //D0
870
871 if(fReadMC){
241994ab 872 if(labD0>=0 && CheckD0Origin(arrMC,(AliAODMCParticle*)arrMC->At(labD0))==4) {
873 AliAODMCParticle *partD0 = (AliAODMCParticle*)arrMC->At(labD0);
a2ad7da1 874 Int_t pdgD0 = partD0->GetPdgCode();
875 if (pdgD0==421){ //D0
241994ab 876 fillthis="histSgn_c_";
a2ad7da1 877 fillthis+=ptbin;
878 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0);
241994ab 879 fillthis="histSgn_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 } else{ //it was a D0bar
884 fillthis="histRfl_c_";
a2ad7da1 885 fillthis+=ptbin;
241994ab 886 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0);
887 fillthis="histRfl_WeigD0Eff_c_";
7f221b36 888 fillthis+=ptbin;
241994ab 889 Double_t effD0 = fCutsTracks->GetTrigWeight(part->Pt(),fMultEv);
7f221b36 890 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0,1./effD0);
241994ab 891 }
892 } else if(labD0>=0 && CheckD0Origin(arrMC,(AliAODMCParticle*)arrMC->At(labD0))==5) {
893 AliAODMCParticle *partD0 = (AliAODMCParticle*)arrMC->At(labD0);
894 Int_t pdgD0 = partD0->GetPdgCode();
895 if (pdgD0==421){ //D0
896 fillthis="histSgn_b_";
897 fillthis+=ptbin;
898 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0);
899 fillthis="histSgn_WeigD0Eff_b_";
900 fillthis+=ptbin;
901 Double_t effD0 = fCutsTracks->GetTrigWeightB(part->Pt(),fMultEv);
902 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0,1./effD0);
903 } else{ //it was a D0bar
904 fillthis="histRfl_b_";
905 fillthis+=ptbin;
906 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0);
907 fillthis="histRfl_WeigD0Eff_b_";
908 fillthis+=ptbin;
909 Double_t effD0 = fCutsTracks->GetTrigWeightB(part->Pt(),fMultEv);
910 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0,1./effD0);
911 }
a2ad7da1 912 } else {//background
241994ab 913 fillthis="histBkg_c_";
a2ad7da1 914 fillthis+=ptbin;
915 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0);
241994ab 916 fillthis="histBkg_WeigD0Eff_c_";
7f221b36 917 fillthis+=ptbin;
241994ab 918 Double_t effD0 = fCutsTracks->GetTrigWeight(part->Pt(),fMultEv);
7f221b36 919 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0,1./effD0);
a2ad7da1 920 }
921 }else{
922 fillthis="histMass_";
923 fillthis+=ptbin;
924 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0);
7f221b36 925 fillthis="histMass_WeigD0Eff_";
926 fillthis+=ptbin;
241994ab 927 Double_t effD0 = fCutsTracks->GetTrigWeight(part->Pt(),fMultEv);
7f221b36 928 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0,1./effD0);
a2ad7da1 929 }
930
931 }
932 if (fIsSelectedCandidate>1 && (fFillOnlyD0D0bar==0 || fFillOnlyD0D0bar==2)) { //D0bar
933
934 if(fReadMC){
241994ab 935 if(labD0>=0 && CheckD0Origin(arrMC,(AliAODMCParticle*)arrMC->At(labD0))==4) {
936 AliAODMCParticle *partD0 = (AliAODMCParticle*)arrMC->At(labD0);
a2ad7da1 937 Int_t pdgD0 = partD0->GetPdgCode();
241994ab 938 if (pdgD0==-421){ //D0
939 fillthis="histSgn_c_";
a2ad7da1 940 fillthis+=ptbin;
941 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar);
241994ab 942 fillthis="histSgn_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 } else{ //it was a D0bar
947 fillthis="histRfl_c_";
a2ad7da1 948 fillthis+=ptbin;
949 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar);
241994ab 950 fillthis="histRfl_WeigD0Eff_c_";
7f221b36 951 fillthis+=ptbin;
241994ab 952 Double_t effD0 = fCutsTracks->GetTrigWeight(part->Pt(),fMultEv);
7f221b36 953 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar,1./effD0);
241994ab 954 }
955 } else if(labD0>=0 && CheckD0Origin(arrMC,(AliAODMCParticle*)arrMC->At(labD0))==5) {
956 AliAODMCParticle *partD0 = (AliAODMCParticle*)arrMC->At(labD0);
957 Int_t pdgD0 = partD0->GetPdgCode();
958 if (pdgD0==-421){ //D0
959 fillthis="histSgn_b_";
960 fillthis+=ptbin;
961 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar);
962 fillthis="histSgn_WeigD0Eff_b_";
963 fillthis+=ptbin;
964 Double_t effD0 = fCutsTracks->GetTrigWeightB(part->Pt(),fMultEv);
965 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar,1./effD0);
966 } else{ //it was a D0bar
967 fillthis="histRfl_b_";
968 fillthis+=ptbin;
969 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar);
970 fillthis="histRfl_WeigD0Eff_b_";
971 fillthis+=ptbin;
972 Double_t effD0 = fCutsTracks->GetTrigWeightB(part->Pt(),fMultEv);
973 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar,1./effD0);
974 }
975 } else {//background
976 fillthis="histBkg_c_";
a2ad7da1 977 fillthis+=ptbin;
978 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar);
241994ab 979 fillthis="histBkg_WeigD0Eff_c_";
7f221b36 980 fillthis+=ptbin;
241994ab 981 Double_t effD0 = fCutsTracks->GetTrigWeight(part->Pt(),fMultEv);
7f221b36 982 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar,1./effD0);
a2ad7da1 983 }
984 }else{
985 fillthis="histMass_";
986 fillthis+=ptbin;
241994ab 987 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar);
7f221b36 988 fillthis="histMass_WeigD0Eff_";
989 fillthis+=ptbin;
241994ab 990 Double_t effD0 = fCutsTracks->GetTrigWeight(part->Pt(),fMultEv);
7f221b36 991 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar,1./effD0);
a2ad7da1 992 }
241994ab 993
a2ad7da1 994 }
995
996 return;
997}
998
999//________________________________________________________________________
1000void AliAnalysisTaskSED0Correlations::Terminate(Option_t */*option*/)
1001{
1002 // Terminate analysis
1003 //
1004 if(fDebug > 1) printf("AnalysisTaskSED0Correlations: Terminate() \n");
1005
1006 fOutputMass = dynamic_cast<TList*> (GetOutputData(1));
1007 if (!fOutputMass) {
1008 printf("ERROR: fOutputMass not available\n");
1009 return;
1010 }
1011
1012 fNentries = dynamic_cast<TH1F*>(GetOutputData(2));
1013
1014 if(!fNentries){
1015 printf("ERROR: fNEntries not available\n");
1016 return;
1017 }
1018
1019 fCutsD0 = dynamic_cast<AliRDHFCutsD0toKpi*>(GetOutputData(3));
1020 if(!fCutsD0){
1021 printf("ERROR: fCuts not available\n");
1022 return;
1023 }
1024
1025 fCounter = dynamic_cast<AliNormalizationCounter*>(GetOutputData(4));
1026 if (!fCounter) {
1027 printf("ERROR: fCounter not available\n");
1028 return;
1029 }
1030 fOutputCorr = dynamic_cast<TList*> (GetOutputData(5));
1031 if (!fOutputCorr) {
1032 printf("ERROR: fOutputCorr not available\n");
1033 return;
1034 }
1035 fOutputStudy = dynamic_cast<TList*> (GetOutputData(6));
1036 if (!fOutputStudy) {
1037 printf("ERROR: fOutputStudy not available\n");
1038 return;
1039 }
1040 fCutsTracks = dynamic_cast<AliHFAssociatedTrackCuts*>(GetOutputData(7));
1041 if(!fCutsTracks){
1042 printf("ERROR: fCutsTracks not available\n");
1043 return;
1044 }
1045
1046 return;
1047}
1048
1049//_________________________________________________________________________________________________
1050Int_t AliAnalysisTaskSED0Correlations::CheckD0Origin(TClonesArray* arrayMC, AliAODMCParticle *mcPartCandidate) const {
1051 //
1052 // checking whether the mother of the particles come from a charm or a bottom quark
1053 //
bce70c96 1054 printf("AliAnalysisTaskSED0Correlations::CheckD0Origin() \n");
a2ad7da1 1055
1056 Int_t pdgGranma = 0;
1057 Int_t mother = 0;
1058 mother = mcPartCandidate->GetMother();
1059 Int_t abspdgGranma =0;
1060 Bool_t isFromB=kFALSE;
1061 Bool_t isQuarkFound=kFALSE;
1062
1063 while (mother > 0){
4d978a4a 1064 AliAODMCParticle* mcMoth = dynamic_cast<AliAODMCParticle*>(arrayMC->At(mother));
1065 if (mcMoth){
1066 pdgGranma = mcMoth->GetPdgCode();
a2ad7da1 1067 abspdgGranma = TMath::Abs(pdgGranma);
1068 if ((abspdgGranma > 500 && abspdgGranma < 600) || (abspdgGranma > 5000 && abspdgGranma < 6000)){
1069 isFromB=kTRUE;
1070 }
1071 if(abspdgGranma==4 || abspdgGranma==5) isQuarkFound=kTRUE;
4d978a4a 1072 mother = mcMoth->GetMother();
a2ad7da1 1073 }else{
1074 AliError("Failed casting the mother particle!");
1075 break;
1076 }
1077 }
1078
1079 if(isQuarkFound) {
1080 if(isFromB) return 5;
1081 else return 4;
1082 }
1083 else return 1;
1084}
1085
1086//________________________________________________________________________
1087void AliAnalysisTaskSED0Correlations::CreateCorrelationsObjs() {
1088//
1089
1090 TString namePlot = "";
1091
1092 //These for limits in THnSparse (one per bin, same limits).
bce70c96 1093 //Vars: DeltaPhi, InvMass, PtTrack, Displacement, DeltaEta
bec72d8c 1094 Int_t nBinsPhi[5] = {32,150,6,3,16};
241994ab 1095 Double_t binMinPhi[5] = {-TMath::Pi()/2.,1.6,0.,0.,-1.6}; //is the minimum for all the bins
1096 Double_t binMaxPhi[5] = {3.*TMath::Pi()/2.,2.2,3.0,3.,1.6}; //is the maximum for all the bins
a2ad7da1 1097
bce70c96 1098 //Vars: DeltaPhi, InvMass, DeltaEta
241994ab 1099 Int_t nBinsMix[4] = {32,150,16,6};
1100 Double_t binMinMix[4] = {-TMath::Pi()/2.,1.6,-1.6,0.}; //is the minimum for all the bins
1101 Double_t binMaxMix[4] = {3.*TMath::Pi()/2.,2.2,1.6,3.}; //is the maximum for all the bins
a2ad7da1 1102
bce70c96 1103 for(Int_t i=0;i<fNPtBinsCorr;i++){
a2ad7da1 1104
bce70c96 1105 if(!fMixing) {
1106 //THnSparse plots: correlations for various invariant mass (MC and data)
1107 namePlot="hPhi_K0_Bin";
a2ad7da1 1108 namePlot+=i;
1109
4d978a4a 1110 THnSparseF *hPhiK = new THnSparseF(namePlot.Data(), "Azimuthal correlation; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
bce70c96 1111 hPhiK->Sumw2();
1112 fOutputCorr->Add(hPhiK);
a2ad7da1 1113
bce70c96 1114 namePlot="hPhi_Kcharg_Bin";
a2ad7da1 1115 namePlot+=i;
1116
4d978a4a 1117 THnSparseF *hPhiH = new THnSparseF(namePlot.Data(), "Azimuthal correlation; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
bce70c96 1118 hPhiH->Sumw2();
1119 fOutputCorr->Add(hPhiH);
a2ad7da1 1120
bce70c96 1121 namePlot="hPhi_Charg_Bin";
a2ad7da1 1122 namePlot+=i;
1123
4d978a4a 1124 THnSparseF *hPhiC = new THnSparseF(namePlot.Data(), "Azimuthal correlation; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi);
bce70c96 1125 hPhiC->Sumw2();
1126 fOutputCorr->Add(hPhiC);
a2ad7da1 1127
bce70c96 1128 //histos for c/b origin for D0 (MC only)
1129 if (fReadMC) {
a2ad7da1 1130
bce70c96 1131 //generic origin for tracks
1132 namePlot="hPhi_K0_From_c_Bin";
1133 namePlot+=i;
a2ad7da1 1134
4d978a4a 1135 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 1136 hPhiK_c->Sumw2();
1137 fOutputCorr->Add(hPhiK_c);
a2ad7da1 1138
bce70c96 1139 namePlot="hPhi_Kcharg_From_c_Bin";
1140 namePlot+=i;
a2ad7da1 1141
4d978a4a 1142 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 1143 hPhiH_c->Sumw2();
1144 fOutputCorr->Add(hPhiH_c);
a2ad7da1 1145
bce70c96 1146 namePlot="hPhi_Charg_From_c_Bin";
1147 namePlot+=i;
a2ad7da1 1148
4d978a4a 1149 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 1150 hPhiC_c->Sumw2();
1151 fOutputCorr->Add(hPhiC_c);
1152
1153 namePlot="hPhi_K0_From_b_Bin";
1154 namePlot+=i;
a2ad7da1 1155
4d978a4a 1156 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 1157 hPhiK_b->Sumw2();
1158 fOutputCorr->Add(hPhiK_b);
a2ad7da1 1159
bce70c96 1160 namePlot="hPhi_Kcharg_From_b_Bin";
1161 namePlot+=i;
a2ad7da1 1162
4d978a4a 1163 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 1164 hPhiH_b->Sumw2();
1165 fOutputCorr->Add(hPhiH_b);
a2ad7da1 1166
bce70c96 1167 namePlot="hPhi_Charg_From_b_Bin";
1168 namePlot+=i;
a2ad7da1 1169
4d978a4a 1170 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 1171 hPhiC_b->Sumw2();
1172 fOutputCorr->Add(hPhiC_b);
a2ad7da1 1173
bce70c96 1174 //HF-only tracks (c for c->D0, b for b->D0)
1175 namePlot="hPhi_K0_HF_From_c_Bin";
1176 namePlot+=i;
a2ad7da1 1177
4d978a4a 1178 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 1179 hPhiK_HF_c->Sumw2();
1180 fOutputCorr->Add(hPhiK_HF_c);
a2ad7da1 1181
bce70c96 1182 namePlot="hPhi_Kcharg_HF_From_c_Bin";
1183 namePlot+=i;
a2ad7da1 1184
4d978a4a 1185 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 1186 hPhiH_HF_c->Sumw2();
1187 fOutputCorr->Add(hPhiH_HF_c);
a2ad7da1 1188
bce70c96 1189 namePlot="hPhi_Charg_HF_From_c_Bin";
1190 namePlot+=i;
a2ad7da1 1191
4d978a4a 1192 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 1193 hPhiC_HF_c->Sumw2();
1194 fOutputCorr->Add(hPhiC_HF_c);
a2ad7da1 1195
bce70c96 1196 namePlot="hPhi_K0_HF_From_b_Bin";
1197 namePlot+=i;
a2ad7da1 1198
4d978a4a 1199 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 1200 hPhiK_HF_b->Sumw2();
1201 fOutputCorr->Add(hPhiK_HF_b);
a2ad7da1 1202
bce70c96 1203 namePlot="hPhi_Kcharg_HF_From_b_Bin";
1204 namePlot+=i;
a2ad7da1 1205
4d978a4a 1206 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 1207 hPhiH_HF_b->Sumw2();
1208 fOutputCorr->Add(hPhiH_HF_b);
a2ad7da1 1209
bce70c96 1210 namePlot="hPhi_Charg_HF_From_b_Bin";
1211 namePlot+=i;
a2ad7da1 1212
4d978a4a 1213 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 1214 hPhiC_HF_b->Sumw2();
1215 fOutputCorr->Add(hPhiC_HF_b);
4d978a4a 1216
1217 namePlot="hPhi_K0_NonHF_Bin";
1218 namePlot+=i;
1219
1220 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);
1221 hPhiK_Non->Sumw2();
1222 fOutputCorr->Add(hPhiK_Non);
1223
1224 namePlot="hPhi_Kcharg_NonHF_Bin";
1225 namePlot+=i;
1226
1227 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);
1228 hPhiH_Non->Sumw2();
1229 fOutputCorr->Add(hPhiH_Non);
1230
1231 namePlot="hPhi_Charg_NonHF_Bin";
1232 namePlot+=i;
1233
1234 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);
1235 hPhiC_Non->Sumw2();
1236 fOutputCorr->Add(hPhiC_Non);
bce70c96 1237 }
a2ad7da1 1238
bce70c96 1239 //leading hadron correlations
1240 namePlot="hPhi_Lead_Bin";
a2ad7da1 1241 namePlot+=i;
1242
241994ab 1243 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 1244 hCorrLead->Sumw2();
1245 fOutputCorr->Add(hCorrLead);
a2ad7da1 1246
bce70c96 1247 if (fReadMC) {
1248 namePlot="hPhi_Lead_From_c_Bin";
1249 namePlot+=i;
a2ad7da1 1250
241994ab 1251 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 1252 hCorrLead_c->Sumw2();
1253 fOutputCorr->Add(hCorrLead_c);
1254
1255 namePlot="hPhi_Lead_From_b_Bin";
1256 namePlot+=i;
1257
241994ab 1258 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 1259 hCorrLead_b->Sumw2();
1260 fOutputCorr->Add(hCorrLead_b);
1261
1262 namePlot="hPhi_Lead_HF_From_c_Bin";
1263 namePlot+=i;
1264
241994ab 1265 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 1266 hCorrLead_HF_c->Sumw2();
1267 fOutputCorr->Add(hCorrLead_HF_c);
1268
1269 namePlot="hPhi_Lead_HF_From_b_Bin";
1270 namePlot+=i;
1271
241994ab 1272 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 1273 hCorrLead_HF_b->Sumw2();
1274 fOutputCorr->Add(hCorrLead_HF_b);
4d978a4a 1275
1276 namePlot="hPhi_Lead_NonHF_Bin";
1277 namePlot+=i;
1278
241994ab 1279 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 1280 hCorrLead_Non->Sumw2();
1281 fOutputCorr->Add(hCorrLead_Non);
bce70c96 1282 }
1283
1284 //pT weighted correlations
1285 namePlot="hPhi_Weig_Bin";
a2ad7da1 1286 namePlot+=i;
bce70c96 1287
241994ab 1288 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 1289 fOutputCorr->Add(hCorrWeig);
1290
1291 if (fReadMC) {
1292 namePlot="hPhi_Weig_From_c_Bin";
1293 namePlot+=i;
1294
241994ab 1295 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 1296 fOutputCorr->Add(hCorrWeig_c);
1297
1298 namePlot="hPhi_Weig_From_b_Bin";
1299 namePlot+=i;
1300
241994ab 1301 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 1302 fOutputCorr->Add(hCorrWeig_b);
1303
1304 namePlot="hPhi_Weig_HF_From_c_Bin";
1305 namePlot+=i;
1306
241994ab 1307 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 1308 fOutputCorr->Add(hCorrWeig_HF_c);
1309
1310 namePlot="hPhi_Weig_HF_From_b_Bin";
1311 namePlot+=i;
1312
241994ab 1313 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 1314 fOutputCorr->Add(hCorrWeig_HF_b);
4d978a4a 1315
1316 namePlot="hPhi_Weig_NonHF_Bin";
1317 namePlot+=i;
1318
241994ab 1319 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 1320 fOutputCorr->Add(hCorrWeig_Non);
bce70c96 1321 }
a2ad7da1 1322
1323 //pT distribution histos
1324 namePlot = "hist_Pt_Charg_Bin"; namePlot+=i;
1325 TH1F *hPtC = new TH1F(namePlot.Data(), "Charged track pT (in D0 evs); p_{T} (GeV/c)",240,0.,12.);
1326 hPtC->SetMinimum(0);
1327 fOutputStudy->Add(hPtC);
1328
1329 namePlot = "hist_Pt_Kcharg_Bin"; namePlot+=i;
1330 TH1F *hPtH = new TH1F(namePlot.Data(), "Hadrons pT (in D0 evs); p_{T} (GeV/c)",240,0.,12.);
1331 hPtH->SetMinimum(0);
1332 fOutputStudy->Add(hPtH);
1333
f80e7bba 1334 namePlot = "hist_Pt_K0_Bin"; namePlot+=i;
a2ad7da1 1335 TH1F *hPtK = new TH1F(namePlot.Data(), "Kaons pT (in D0 evs); p_{T} (GeV/c)",240,0.,12.);
1336 hPtK->SetMinimum(0);
1337 fOutputStudy->Add(hPtK);
1338
1339 //D* feeddown pions rejection histos
1340 namePlot = "hDstarPions_Bin"; namePlot+=i;
4d978a4a 1341 TH2F *hDstarPions = new TH2F(namePlot.Data(), "Tracks rejected for D* inv.mass cut; # Tracks",2,0.,2.,150,1.6,2.2);
a2ad7da1 1342 hDstarPions->GetXaxis()->SetBinLabel(1,"Not rejected");
1343 hDstarPions->GetXaxis()->SetBinLabel(2,"Rejected");
1344 hDstarPions->SetMinimum(0);
bce70c96 1345 fOutputStudy->Add(hDstarPions);
241994ab 1346
1347 //Events multiplicity
1348 Double_t yAxisMult[13] = {0, 4, 8, 12, 16, 20, 28, 36, 44, 100};
1349 namePlot = "hMultiplEvt_Bin"; namePlot+=i;
1350 TH1F *hMultEv = new TH1F(namePlot.Data(), "Event multiplicity",9,yAxisMult);
1351 hMultEv->SetMinimum(0);
1352 fOutputStudy->Add(hMultEv);
bce70c96 1353
1354 }
1355
1356 if(fMixing) {
1357 //THnSparse plots for event mixing!
1358 namePlot="hPhi_K0_Bin";
1359 namePlot+=i;namePlot+="_EvMix";
1360
241994ab 1361 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 1362 hPhiK_EvMix->Sumw2();
1363 fOutputCorr->Add(hPhiK_EvMix);
1364
1365 namePlot="hPhi_Kcharg_Bin";
1366 namePlot+=i;namePlot+="_EvMix";
1367
241994ab 1368 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 1369 hPhiH_EvMix->Sumw2();
1370 fOutputCorr->Add(hPhiH_EvMix);
1371
1372 namePlot="hPhi_Charg_Bin";
1373 namePlot+=i;namePlot+="_EvMix";
1374
241994ab 1375 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 1376 hPhiC_EvMix->Sumw2();
1377 fOutputCorr->Add(hPhiC_EvMix);
1378 }
a2ad7da1 1379 }
1380
1381 //out of bin loop
da975030 1382 if(!fMixing) {
1383 TH1F *hCountC = new TH1F("hist_Count_Charg", "Charged track counter; # Tracks",100,0.,100.);
1384 hCountC->SetMinimum(0);
1385 fOutputStudy->Add(hCountC);
bce70c96 1386
da975030 1387 TH1F *hCountH = new TH1F("hist_Count_Kcharg", "Hadrons counter; # Tracks",20,0.,20.);
1388 hCountH->SetMinimum(0);
1389 fOutputStudy->Add(hCountH);
bce70c96 1390
da975030 1391 TH1F *hCountK = new TH1F("hist_Count_K0", "Kaons counter; # Tracks",20,0.,20.);
1392 hCountK->SetMinimum(0);
1393 fOutputStudy->Add(hCountK);
1394 }
bce70c96 1395
4d978a4a 1396 if (fReadMC) {
1397 TH1D *hEventTypeMC = new TH1D("EventTypeMC","EventTypeMC",100,-0.5,99.5);
1398 fOutputStudy->Add(hEventTypeMC);
1399 }
1400
a2ad7da1 1401 if (fFillGlobal) { //all-events plots
1402 //pt distributions
1403 TH1F *hPtCAll = new TH1F("hist_Pt_Charg_AllEv", "Charged track pT (All); p_{T} (GeV/c)",240,0.,12.);
1404 hPtCAll->SetMinimum(0);
1405 fOutputStudy->Add(hPtCAll);
1406
1407 TH1F *hPtHAll = new TH1F("hist_Pt_Kcharg_AllEv", "Hadrons pT (All); p_{T} (GeV/c)",240,0.,12.);
1408 hPtHAll->SetMinimum(0);
1409 fOutputStudy->Add(hPtHAll);
1410
bce70c96 1411 TH1F *hPtKAll = new TH1F("hist_Pt_K0_AllEv", "Kaons pT (All); p_{T} (GeV/c)",240,0.,12.);
a2ad7da1 1412 hPtKAll->SetMinimum(0);
1413 fOutputStudy->Add(hPtKAll);
1414
1415 //K0 Invariant Mass plots
1416 TH2F *hK0MassInv = new TH2F("hK0MassInv", "K0 invariant mass; Invariant mass (MeV/c^{2}); pT (GeV/c)",200,0.4,0.6,100,0.,10.);
1417 hK0MassInv->SetMinimum(0);
1418 fOutputStudy->Add(hK0MassInv);
1419 }
1420
dea90a65 1421 if(!fMixing) {
1422 //phi distributions
1423 TH1F *hPhiDistCAll = new TH1F("hist_PhiDistr_Charg", "Charged track phi distr. (All); p_{T} (GeV/c)",64,0,6.283);
1424 hPhiDistCAll->SetMinimum(0);
1425 fOutputStudy->Add(hPhiDistCAll);
1426
1427 TH1F *hPhiDistHAll = new TH1F("hist_PhiDistr_Kcharg", "Hadrons phi distr. (All); p_{T} (GeV/c)",64,0,6.283);
1428 hPhiDistHAll->SetMinimum(0);
1429 fOutputStudy->Add(hPhiDistHAll);
1430
1431 TH1F *hPhiDistKAll = new TH1F("hist_PhiDistr_K0", "Kaons phi distr. (All); p_{T} (GeV/c)",64,0,6.283);
1432 hPhiDistKAll->SetMinimum(0);
1433 fOutputStudy->Add(hPhiDistKAll);
1434
1435 TH1F *hPhiDistDAll = new TH1F("hist_PhiDistr_D0", "D^{0} phi distr. (All); p_{T} (GeV/c)",64,0,6.283);
1436 hPhiDistDAll->SetMinimum(0);
1437 fOutputStudy->Add(hPhiDistDAll);
1438 }
1439
a2ad7da1 1440 //for MC analysis only
da975030 1441 for(Int_t i=0;i<fNPtBinsCorr;i++) {
a2ad7da1 1442
da975030 1443 if (fReadMC && !fMixing) {
a2ad7da1 1444
1445 //displacement histos
f80e7bba 1446 namePlot="histDispl_K0_Bin"; namePlot+=i;
a2ad7da1 1447 TH1F *hDisplK = new TH1F(namePlot.Data(), "Kaons Displacement; DCA",150,0.,0.15);
1448 hDisplK->SetMinimum(0);
1449 fOutputStudy->Add(hDisplK);
1450
f80e7bba 1451 namePlot="histDispl_K0_HF_Bin"; namePlot+=i;
a2ad7da1 1452 TH1F *hDisplK_HF = new TH1F(namePlot.Data(), "Kaons Displacement (from HF decay only); DCA",150,0.,0.15);
1453 hDisplK_HF->SetMinimum(0);
1454 fOutputStudy->Add(hDisplK_HF);
1455
1456 namePlot="histDispl_Kcharg_Bin"; namePlot+=i;
1457 TH1F *hDisplHadr = new TH1F(namePlot.Data(), "Hadrons Displacement; DCA",150,0.,0.15);
1458 hDisplHadr->SetMinimum(0);
1459 fOutputStudy->Add(hDisplHadr);
1460
1461 namePlot="histDispl_Kcharg_HF_Bin"; namePlot+=i;
1462 TH1F *hDisplHadr_HF = new TH1F(namePlot.Data(), "Hadrons Displacement (from HF decay only); DCA",150,0.,0.15);
1463 hDisplHadr_HF->SetMinimum(0);
1464 fOutputStudy->Add(hDisplHadr_HF);
1465
1466 namePlot="histDispl_Charg_Bin"; namePlot+=i;
1467 TH1F *hDisplCharg = new TH1F(namePlot.Data(), "Charged tracks Displacement; DCA",150,0.,0.15);
1468 hDisplCharg->SetMinimum(0);
1469 fOutputStudy->Add(hDisplCharg);
1470
1471 namePlot="histDispl_Charg_HF_Bin"; namePlot+=i;
1472 TH1F *hDisplCharg_HF = new TH1F(namePlot.Data(), "Charged tracks Displacement (from HF decay only); DCA",150,0.,0.15);
1473 hDisplCharg_HF->SetMinimum(0);
1474 fOutputStudy->Add(hDisplCharg_HF);
1475
f80e7bba 1476 namePlot="histDispl_K0_From_c_Bin"; namePlot+=i;
a2ad7da1 1477 TH1F *hDisplK_c = new TH1F(namePlot.Data(), "Kaons Displacement - c origin; DCA",150,0.,0.15);
1478 hDisplK_c->SetMinimum(0);
1479 fOutputStudy->Add(hDisplK_c);
1480
f80e7bba 1481 namePlot="histDispl_K0_HF_From_c_Bin"; namePlot+=i;
a2ad7da1 1482 TH1F *hDisplK_HF_c = new TH1F(namePlot.Data(), "Kaons Displacement (from HF decay only) - c origin; DCA",150,0.,0.15);
1483 hDisplK_HF_c->SetMinimum(0);
1484 fOutputStudy->Add(hDisplK_HF_c);
1485
1486 namePlot="histDispl_Kcharg_From_c_Bin"; namePlot+=i;
1487 TH1F *hDisplHadr_c = new TH1F(namePlot.Data(), "Hadrons Displacement - c origin; DCA",150,0.,0.15);
1488 hDisplHadr_c->SetMinimum(0);
1489 fOutputStudy->Add(hDisplHadr_c);
1490
1491 namePlot="histDispl_Kcharg_HF_From_c_Bin"; namePlot+=i;
1492 TH1F *hDisplHadr_HF_c = new TH1F(namePlot.Data(), "Hadrons Displacement (from HF decay only) - c origin; DCA",150,0.,0.15);
1493 hDisplHadr_HF_c->SetMinimum(0);
1494 fOutputStudy->Add(hDisplHadr_HF_c);
1495
1496 namePlot="histDispl_Charg_From_c_Bin"; namePlot+=i;
1497 TH1F *hDisplCharg_c = new TH1F(namePlot.Data(), "Charged tracks Displacement - c origin; DCA",150,0.,0.15);
1498 hDisplCharg_c->Sumw2();
1499 hDisplCharg_c->SetMinimum(0);
1500 fOutputStudy->Add(hDisplCharg_c);
1501
1502 namePlot="histDispl_Charg_HF_From_c_Bin"; namePlot+=i;
1503 TH1F *hDisplCharg_HF_c = new TH1F(namePlot.Data(), "Charged tracks Displacement (from HF decay only) - c origin; DCA",150,0.,0.15);
1504 hDisplCharg_HF_c->SetMinimum(0);
1505 fOutputStudy->Add(hDisplCharg_HF_c);
1506
f80e7bba 1507 namePlot="histDispl_K0_From_b_Bin"; namePlot+=i;
a2ad7da1 1508 TH1F *hDisplK_b = new TH1F(namePlot.Data(), "Kaons Displacement - b origin; DCA",150,0.,0.15);
1509 hDisplK_b->SetMinimum(0);
1510 fOutputStudy->Add(hDisplK_b);
1511
f80e7bba 1512 namePlot="histDispl_K0_HF_From_b_Bin"; namePlot+=i;
a2ad7da1 1513 TH1F *hDisplK_HF_b = new TH1F(namePlot.Data(), "Kaons Displacement (from HF decay only) - b origin; DCA",150,0.,0.15);
1514 hDisplK_HF_b->SetMinimum(0);
1515 fOutputStudy->Add(hDisplK_HF_b);
1516
1517 namePlot="histDispl_Kcharg_From_b_Bin"; namePlot+=i;
1518 TH1F *hDisplHadr_b = new TH1F(namePlot.Data(), "Hadrons Displacement - b origin; DCA",150,0.,0.15);
1519 hDisplHadr_b->SetMinimum(0);
1520 fOutputStudy->Add(hDisplHadr_b);
1521
1522 namePlot="histDispl_Kcharg_HF_From_b_Bin"; namePlot+=i;
1523 TH1F *hDisplHadr_HF_b = new TH1F(namePlot.Data(), "Hadrons Displacement (from HF decay only) - b origin; DCA",150,0.,0.15);
1524 hDisplHadr_HF_b->SetMinimum(0);
1525 fOutputStudy->Add(hDisplHadr_HF_b);
1526
1527 namePlot="histDispl_Charg_From_b_Bin"; namePlot+=i;
1528 TH1F *hDisplCharg_b = new TH1F(namePlot.Data(), "Charged tracks Displacement - b origin; DCA",150,0.,0.15);
1529 hDisplCharg_b->SetMinimum(0);
1530 fOutputStudy->Add(hDisplCharg_b);
1531
1532 namePlot="histDispl_Charg_HF_From_b_Bin"; namePlot+=i;
1533 TH1F *hDisplCharg_HF_b = new TH1F(namePlot.Data(), "Charged tracks Displacement (from HF decay only) - b origin; DCA",150,0.,0.15);
1534 hDisplCharg_HF_b->SetMinimum(0);
1535 fOutputStudy->Add(hDisplCharg_HF_b);
1536
1537 //origin of tracks histos
1538 namePlot="histOrig_Charg_Bin"; namePlot+=i;
1539 TH1F *hOrigin_Charm = new TH1F(namePlot.Data(), "Origin of charged tracks",9,0.,9.);
1540 hOrigin_Charm->SetMinimum(0);
1541 hOrigin_Charm->GetXaxis()->SetBinLabel(1,"Not HF");
1542 hOrigin_Charm->GetXaxis()->SetBinLabel(2,"D->#");
1543 hOrigin_Charm->GetXaxis()->SetBinLabel(3,"D->X->#");
4d978a4a 1544 hOrigin_Charm->GetXaxis()->SetBinLabel(4,"c hadr.");
1545 hOrigin_Charm->GetXaxis()->SetBinLabel(5,"B->#");
1546 hOrigin_Charm->GetXaxis()->SetBinLabel(6,"B->X-># (X!=D)");
1547 hOrigin_Charm->GetXaxis()->SetBinLabel(7,"B->D->#");
1548 hOrigin_Charm->GetXaxis()->SetBinLabel(8,"B->D->X->#");
a2ad7da1 1549 hOrigin_Charm->GetXaxis()->SetBinLabel(9,"b hadr.");
1550 fOutputStudy->Add(hOrigin_Charm);
1551
1552 namePlot="histOrig_Kcharg_Bin"; namePlot+=i;
1553 TH1F *hOrigin_Kcharg = new TH1F(namePlot.Data(), "Origin of hadrons",9,0.,9.);
1554 hOrigin_Kcharg->SetMinimum(0);
1555 hOrigin_Kcharg->GetXaxis()->SetBinLabel(1,"Not HF");
1556 hOrigin_Kcharg->GetXaxis()->SetBinLabel(2,"D->#");
1557 hOrigin_Kcharg->GetXaxis()->SetBinLabel(3,"D->X->#");
4d978a4a 1558 hOrigin_Kcharg->GetXaxis()->SetBinLabel(4,"c hadr.");
1559 hOrigin_Kcharg->GetXaxis()->SetBinLabel(5,"B->#");
1560 hOrigin_Kcharg->GetXaxis()->SetBinLabel(6,"B->X-># (X!=D)");
1561 hOrigin_Kcharg->GetXaxis()->SetBinLabel(7,"B->D->#");
1562 hOrigin_Kcharg->GetXaxis()->SetBinLabel(8,"B->D->X->#");
a2ad7da1 1563 hOrigin_Kcharg->GetXaxis()->SetBinLabel(9,"b hadr.");
1564 fOutputStudy->Add(hOrigin_Kcharg);
1565
f80e7bba 1566 namePlot="histOrig_K0_Bin"; namePlot+=i;
a2ad7da1 1567 TH1F *hOrigin_K = new TH1F(namePlot.Data(), "Origin of kaons",9,0.,9.);
1568 hOrigin_K->SetMinimum(0);
1569 hOrigin_K->GetXaxis()->SetBinLabel(1,"Not HF");
1570 hOrigin_K->GetXaxis()->SetBinLabel(2,"D->#");
1571 hOrigin_K->GetXaxis()->SetBinLabel(3,"D->X->#");
4d978a4a 1572 hOrigin_K->GetXaxis()->SetBinLabel(4,"c hadr.");
1573 hOrigin_K->GetXaxis()->SetBinLabel(5,"B->#");
1574 hOrigin_K->GetXaxis()->SetBinLabel(6,"B->X-># (X!=D)");
1575 hOrigin_K->GetXaxis()->SetBinLabel(7,"B->D->#");
1576 hOrigin_K->GetXaxis()->SetBinLabel(8,"B->D->X->#");
a2ad7da1 1577 hOrigin_K->GetXaxis()->SetBinLabel(9,"b hadr.");
1578 fOutputStudy->Add(hOrigin_K);
da975030 1579 }
a2ad7da1 1580
da975030 1581 if (fReadMC) {
a2ad7da1 1582 //origin of D0 histos
1583 namePlot="histOrig_D0_Bin"; namePlot+=i;
1584 TH1F *hOrigin_D0 = new TH1F(namePlot.Data(), "Origin of D0",2,0.,2.);
1585 hOrigin_D0->SetMinimum(0);
1586 hOrigin_D0->GetXaxis()->SetBinLabel(1,"From c");
1587 hOrigin_D0->GetXaxis()->SetBinLabel(2,"From b");
1588 fOutputStudy->Add(hOrigin_D0);
241994ab 1589
1590 //primary tracks (Kine & Reco)
1591 namePlot="hPhysPrim_Bin"; namePlot+=i;
1592 TH1F *hPhysPrim = new TH1F(namePlot.Data(), "Origin of hadrons",2,0.,2.);
1593 hPhysPrim->SetMinimum(0);
1594 hPhysPrim->GetXaxis()->SetBinLabel(1,"OK");
1595 hPhysPrim->GetXaxis()->SetBinLabel(2,"NO");
1596 fOutputStudy->Add(hPhysPrim);
a2ad7da1 1597 }
1598 }
a2ad7da1 1599}
1600
1601//________________________________________________________________________
bce70c96 1602void AliAnalysisTaskSED0Correlations::CalculateCorrelations(AliAODRecoDecayHF2Prong* d, Int_t labD0, TClonesArray* mcArray) {
a2ad7da1 1603//
1604// Method for correlations D0-hadrons study
1605//
bce70c96 1606 Int_t N_Charg = 0, N_KCharg = 0, N_Kaons = 0;
1607 Double_t mD0, mD0bar;
4d978a4a 1608 Int_t origD0 = 0, PDGD0 = 0, ptbin = 0;
a2ad7da1 1609 d->InvMassD0(mD0,mD0bar);
4d978a4a 1610 Double_t mInv[2] = {mD0, mD0bar};
1611 ptbin = PtBinCorr(d->Pt());
1612
a2ad7da1 1613 if(ptbin < 0) return;
a2ad7da1 1614
bec72d8c 1615 //Fill of D0 phi distribution
4d978a4a 1616 if (!fMixing) ((TH1F*)fOutputStudy->FindObject("hist_PhiDistr_D0"))->Fill(d->Phi());
bec72d8c 1617
bce70c96 1618 //Origin of D0
1619 TString orig="";
1620 if(fReadMC) {
1621 origD0=CheckD0Origin(mcArray,(AliAODMCParticle*)mcArray->At(labD0));
1622 PDGD0 = ((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode();
1623 switch (CheckD0Origin(mcArray,(AliAODMCParticle*)mcArray->At(labD0))) {
1624 case 4:
1625 orig = "_From_c";
1626 ((TH1F*)fOutputStudy->FindObject(Form("histOrig_D0_Bin%d",ptbin)))->Fill(0.);
1627 break;
1628 case 5:
1629 orig = "_From_b";
1630 ((TH1F*)fOutputStudy->FindObject(Form("histOrig_D0_Bin%d",ptbin)))->Fill(1.);
1631 break;
1632 default:
1633 return;
1634 }
1635 }
1636
4d978a4a 1637 Double_t highPt = 0; Double_t lead[4] = {0,0,0,1}; //infos for leading particle (pt,deltaphi)
a2ad7da1 1638
bce70c96 1639 //loop over the tracks in the pool
ad8d2dfd 1640 Bool_t execPoolTr = fCorrelatorTr->ProcessEventPool(); //pool is ready? (only in ME, in SE returns kFALSE)
1641 Bool_t execPoolKc = fCorrelatorKc->ProcessEventPool(); //pool is ready? (only in ME, in SE returns kFALSE)
1642 Bool_t execPoolK0 = fCorrelatorK0->ProcessEventPool(); //pool is ready? (only in ME, in SE returns kFALSE)
bce70c96 1643
1644 Int_t NofEventsinPool = 1;
ad8d2dfd 1645 if(fMixing) {
1646 NofEventsinPool = fCorrelatorTr->GetNofEventsInPool();
1647 if(!execPoolTr) {
1648 AliInfo("Mixed event analysis: track pool is not ready");
1649 NofEventsinPool = 0;
1650 }
1651 }
1652
1653 //Charged tracks
bce70c96 1654 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)
1655 Bool_t analyzetracksTr = fCorrelatorTr->ProcessAssociatedTracks(jMix);// process all the tracks in the aodEvent, by applying the selection cuts
ad8d2dfd 1656 if(!analyzetracksTr) {
bce70c96 1657 AliInfo("AliHFCorrelator::Cannot process the track array");
1658 continue;
1659 }
ad8d2dfd 1660
bce70c96 1661 for(Int_t iTrack = 0; iTrack<fCorrelatorTr->GetNofTracks(); iTrack++){ // looping on track candidates
1662
1663 Bool_t runcorrelation = fCorrelatorTr->Correlate(iTrack);
1664 if(!runcorrelation) continue;
1665
1666 AliReducedParticle* track = fCorrelatorTr->GetAssociatedParticle();
1667
1668 if(!fMixing) {
241994ab 1669 if(fSoftPiCut && !track->CheckSoftPi()) { //removal of soft pions
bce70c96 1670 if (fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(1.,mD0);
1671 if (fIsSelectedCandidate >= 2) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(1.,mD0bar);
1672 continue;
1673 } else { //not a soft pion
1674 if (fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(0.,mD0);
1675 if (fIsSelectedCandidate >= 2) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(0.,mD0bar);
7f221b36 1676 }
4d978a4a 1677 Int_t idDaughs[2] = {((AliVTrack*)d->GetDaughter(0))->GetID(),((AliVTrack*)d->GetDaughter(1))->GetID()}; //IDs of daughters to be skipped
1678 if(track->GetID() == idDaughs[0] || track->GetID() == idDaughs[1]) continue; //discards daughters of candidate
ad8d2dfd 1679 }
1680 if(track->Pt() < fPtThreshLow.at(ptbin) || track->Pt() > fPtThreshUp.at(ptbin)) continue; //discard tracks outside pt range for hadrons/K
4d978a4a 1681
241994ab 1682 if(fReadMC) {
1683 AliAODMCParticle* trkKine = (AliAODMCParticle*)mcArray->At(track->GetLabel());
1684 if (!trkKine) continue;
1685 if (!trkKine->IsPhysicalPrimary()) {
1686 ((TH1F*)fOutputStudy->FindObject(Form("hPhysPrim_Bin%d",ptbin)))->Fill(1.);
1687 continue; //reject the Reco track if correspondent Kine track is not primary
1688 } else ((TH1F*)fOutputStudy->FindObject(Form("hPhysPrim_Bin%d",ptbin)))->Fill(0.);
1689 }
1690
7f221b36 1691 Double_t effTr = track->GetWeight(); //extract track efficiency
241994ab 1692 Double_t effD0 = 1.;
1693 if(fReadMC) {
1694 if(origD0==4) effD0 = fCutsTracks->GetTrigWeight(d->Pt(),fMultEv);
1695 if(origD0==5) effD0 = fCutsTracks->GetTrigWeightB(d->Pt(),fMultEv);
1696 } else effD0 = fCutsTracks->GetTrigWeight(d->Pt(),fMultEv);
7f221b36 1697 Double_t eff = effTr*effD0;
241994ab 1698
1699 FillSparsePlots(mcArray,mInv,origD0,PDGD0,track,ptbin,kTrack,1./eff); //fills for charged tracks
bce70c96 1700
1701 if(!fMixing) N_Charg++;
1702
1703 //retrieving leading info...
1704 if(track->Pt() > highPt) {
da975030 1705 if(fReadMC && track->GetLabel()<1) continue;
8c2d7467 1706 if(fReadMC && !(AliAODMCParticle*)mcArray->At(track->GetLabel())) continue;
bce70c96 1707 lead[0] = fCorrelatorTr->GetDeltaPhi();
1708 lead[1] = fCorrelatorTr->GetDeltaEta();
1709 lead[2] = fReadMC ? CheckTrackOrigin(mcArray,(AliAODMCParticle*)mcArray->At(track->GetLabel())) : 0;
241994ab 1710 if(fReadMC) {
1711 if(origD0==4) lead[3] = 1./(track->GetWeight()*fCutsTracks->GetTrigWeight(d->Pt(),fMultEv)); //weight is 1./efficiency
1712 if(origD0==5) lead[3] = 1./(track->GetWeight()*fCutsTracks->GetTrigWeightB(d->Pt(),fMultEv)); //weight is 1./efficiency
1713 } else lead[3] = 1./(track->GetWeight()*fCutsTracks->GetTrigWeight(d->Pt(),fMultEv));
bce70c96 1714 highPt = track->Pt();
a2ad7da1 1715 }
bce70c96 1716
1717 } // end of tracks loop
ad8d2dfd 1718 } //end of event loop for fCorrelatorTr
1719
dea90a65 1720 if(fKaonCorr) { //loops for Kcharg and K0
1721
ad8d2dfd 1722 if(fMixing) {
1723 NofEventsinPool = fCorrelatorKc->GetNofEventsInPool();
1724 if(!execPoolKc) {
1725 AliInfo("Mixed event analysis: K+/- pool is not ready");
1726 NofEventsinPool = 0;
1727 }
1728 }
1729
1730 //Charged Kaons loop
4d978a4a 1731 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 1732 Bool_t analyzetracksKc = fCorrelatorKc->ProcessAssociatedTracks(jMix);
1733 if(!analyzetracksKc) {
1734 AliInfo("AliHFCorrelator::Cannot process the K+/- array");
1735 continue;
1736 }
bce70c96 1737
bce70c96 1738 for(Int_t iTrack = 0; iTrack<fCorrelatorKc->GetNofTracks(); iTrack++){ // looping on charged kaons candidates
1739
1740 Bool_t runcorrelation = fCorrelatorKc->Correlate(iTrack);
1741 if(!runcorrelation) continue;
1742
1743 AliReducedParticle* kCharg = fCorrelatorKc->GetAssociatedParticle();
1744
1745 if(!fMixing) {
241994ab 1746 if(fSoftPiCut && !kCharg->CheckSoftPi()) { //removal of soft pions
bce70c96 1747 if (fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(1.,mD0);
1748 if (fIsSelectedCandidate >= 2) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(1.,mD0bar);
1749 continue;
1750 } else {
1751 if (fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(0.,mD0);
1752 if (fIsSelectedCandidate >= 2) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(0.,mD0bar);
a2ad7da1 1753 }
4d978a4a 1754 Int_t idDaughs[2] = {((AliVTrack*)d->GetDaughter(0))->GetID(),((AliVTrack*)d->GetDaughter(1))->GetID()}; //IDs of daughters to be skipped
1755 if(kCharg->GetID() == idDaughs[0] || kCharg->GetID() == idDaughs[1]) continue; //discards daughters of candidate
ad8d2dfd 1756 }
1757 if(kCharg->Pt() < fPtThreshLow.at(ptbin) || kCharg->Pt() > fPtThreshUp.at(ptbin)) continue; //discard tracks outside pt range for hadrons/K
4d978a4a 1758
1759 FillSparsePlots(mcArray,mInv,origD0,PDGD0,kCharg,ptbin,kKCharg); //fills for charged tracks
a2ad7da1 1760
bce70c96 1761 if(!fMixing) N_KCharg++;
a2ad7da1 1762
bce70c96 1763 } // end of charged kaons loop
ad8d2dfd 1764 } //end of event loop for fCorrelatorKc
1765
1766 if(fMixing) {
1767 NofEventsinPool = fCorrelatorK0->GetNofEventsInPool();
1768 if(!execPoolK0) {
1769 AliInfo("Mixed event analysis: K0 pool is not ready");
1770 NofEventsinPool = 0;
1771 }
1772 }
1773
1774 //K0 loop
1775 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)
1776 Bool_t analyzetracksK0 = fCorrelatorK0->ProcessAssociatedTracks(jMix);
1777 if(!analyzetracksK0) {
1778 AliInfo("AliHFCorrelator::Cannot process the K0 array");
1779 continue;
1780 }
a2ad7da1 1781
bce70c96 1782 for(Int_t iTrack = 0; iTrack<fCorrelatorK0->GetNofTracks(); iTrack++){ // looping on k0 candidates
1783
1784 Bool_t runcorrelation = fCorrelatorK0->Correlate(iTrack);
1785 if(!runcorrelation) continue;
1786
1787 AliReducedParticle* k0 = fCorrelatorK0->GetAssociatedParticle();
a2ad7da1 1788
bce70c96 1789 if(k0->Pt() < fPtThreshLow.at(ptbin) || k0->Pt() > fPtThreshUp.at(ptbin)) continue; //discard tracks outside pt range for hadrons/K
4d978a4a 1790
1791 FillSparsePlots(mcArray,mInv,origD0,PDGD0,k0,ptbin,kK0); //fills for charged tracks
bce70c96 1792
1793 if(!fMixing) N_Kaons++;
1794
1795 } // end of charged kaons loop
ad8d2dfd 1796 } //end of event loop for fCorrelatorK0
bce70c96 1797
dea90a65 1798 } //end of 'if(fKaonCorr)'
1799
241994ab 1800 Double_t fillSpLeadD0[4] = {lead[0],mD0,lead[1],0.4}; //dummy value for threshold of leading!
1801 Double_t fillSpLeadD0bar[4] = {lead[0],mD0bar,lead[1],0.4};
4d978a4a 1802
bce70c96 1803 //leading track correlations fill
ad8d2dfd 1804 if(!fMixing) {
1805 if(fReadMC) {
241994ab 1806 if(((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode()==421 && (fIsSelectedCandidate==1||fIsSelectedCandidate==3)) { //D0
4d978a4a 1807 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(fillSpLeadD0,lead[3]); //c and b D0
1808 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadD0,lead[3]); //c or b D0
1809 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]);
1810 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]);
1811 if((int)lead[2]==0) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_NonHF_Bin%d",ptbin)))->Fill(fillSpLeadD0,lead[3]); //non HF
ad8d2dfd 1812 }
241994ab 1813 if(((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode()==-421 && fIsSelectedCandidate>1) { //D0bar
4d978a4a 1814 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(fillSpLeadD0bar,lead[3]);
1815 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadD0bar,lead[3]); //c or b D0
1816 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]);
1817 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]);
1818 if((int)lead[2]==0) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_NonHF_Bin%d",ptbin)))->Fill(fillSpLeadD0bar,lead[3]); //non HF
ad8d2dfd 1819 }
1820 } else {
4d978a4a 1821 if(fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(fillSpLeadD0,lead[3]);
1822 if(fIsSelectedCandidate == 2 || fIsSelectedCandidate == 3) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(fillSpLeadD0bar,lead[3]);
a2ad7da1 1823 }
1824
ad8d2dfd 1825 //Fill of count histograms
1826 if (!fAlreadyFilled) {
1827 ((TH1F*)fOutputStudy->FindObject("hist_Count_Charg"))->Fill(N_Charg);
1828 ((TH1F*)fOutputStudy->FindObject("hist_Count_Kcharg"))->Fill(N_KCharg);
1829 ((TH1F*)fOutputStudy->FindObject("hist_Count_K0"))->Fill(N_Kaons);
1830 }
bce70c96 1831 }
a2ad7da1 1832
bce70c96 1833 fAlreadyFilled=kTRUE; //at least a D0 analyzed in the event; distribution plots already filled
a2ad7da1 1834
bce70c96 1835}
1836
1837//________________________________________________________________________
4d978a4a 1838void AliAnalysisTaskSED0Correlations::CalculateCorrelationsMCKine(AliAODMCParticle* d, TClonesArray* mcArray) {
1839//
1840// Method for correlations D0-hadrons study
1841//
1842 Int_t N_Charg = 0, N_KCharg = 0, N_Kaons = 0;
1843 Double_t mD0 = 1.864, mD0bar = 1.864;
1844 Double_t mInv[2] = {mD0, mD0bar};
1845 Int_t origD0 = 0, PDGD0 = 0;
1846 Int_t ptbin = PtBinCorr(d->Pt());
1847
1848 if(ptbin < 0) return;
1849
1850 //Fill of D0 phi distribution
1851 if (!fMixing) ((TH1F*)fOutputStudy->FindObject("hist_PhiDistr_D0"))->Fill(d->Phi());
1852
1853 //Origin of D0
1854 TString orig="";
1855 origD0=CheckD0Origin(mcArray,d);
1856 PDGD0 = d->GetPdgCode();
1857 switch (CheckD0Origin(mcArray,d)) {
1858 case 4:
1859 orig = "_From_c";
1860 ((TH1F*)fOutputStudy->FindObject(Form("histOrig_D0_Bin%d",ptbin)))->Fill(0.);
1861 break;
1862 case 5:
1863 orig = "_From_b";
1864 ((TH1F*)fOutputStudy->FindObject(Form("histOrig_D0_Bin%d",ptbin)))->Fill(1.);
1865 break;
1866 default:
1867 return;
1868 }
1869
1870 Double_t highPt = 0; Double_t lead[3] = {0,0,0}; //infos for leading particle (pt,deltaphi)
1871
1872 //loop over the tracks in the pool
1873 Bool_t execPoolTr = fCorrelatorTr->ProcessEventPool(); //pool is ready? (only in ME, in SE returns kFALSE)
1874 Bool_t execPoolKc = fCorrelatorKc->ProcessEventPool(); //pool is ready? (only in ME, in SE returns kFALSE)
1875 Bool_t execPoolK0 = fCorrelatorK0->ProcessEventPool(); //pool is ready? (only in ME, in SE returns kFALSE)
1876
1877 Int_t NofEventsinPool = 1;
1878 if(fMixing) {
1879 NofEventsinPool = fCorrelatorTr->GetNofEventsInPool();
1880 if(!execPoolTr) {
1881 AliInfo("Mixed event analysis: track pool is not ready");
1882 NofEventsinPool = 0;
1883 }
1884 }
1885
1886 //Charged tracks
1887 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)
1888
1889 Bool_t analyzetracksTr = fCorrelatorTr->ProcessAssociatedTracks(jMix);// process all the tracks in the aodEvent, by applying the selection cuts
1890 if(!analyzetracksTr) {
1891 AliInfo("AliHFCorrelator::Cannot process the track array");
1892 continue;
1893 }
1894
1895 for(Int_t iTrack = 0; iTrack<fCorrelatorTr->GetNofTracks(); iTrack++){ // looping on track candidates
1896
1897 Bool_t runcorrelation = fCorrelatorTr->Correlate(iTrack);
1898 if(!runcorrelation) continue;
1899
1900 AliReducedParticle* track = fCorrelatorTr->GetAssociatedParticle();
1901 if(track->GetLabel()<0) continue;
1902 if(track->Pt() < fPtThreshLow.at(ptbin) || track->Pt() > fPtThreshUp.at(ptbin)) continue; //discard tracks outside pt range for hadrons/K
1903 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
1904 if(!fMixing) N_Charg++;
1905
1906 AliAODMCParticle *trkMC = (AliAODMCParticle*)mcArray->At(track->GetLabel());
1907 if(!trkMC) continue;
1908
241994ab 1909 if (!trkMC->IsPhysicalPrimary()) { //reject material budget, or other fake tracks
1910 ((TH1F*)fOutputStudy->FindObject(Form("hPhysPrim_Bin%d",ptbin)))->Fill(1.);
1911 continue;
1912 } else ((TH1F*)fOutputStudy->FindObject(Form("hPhysPrim_Bin%d",ptbin)))->Fill(0.);
1913
4d978a4a 1914 if (IsDDaughter(d,trkMC,mcArray)) continue;
241994ab 1915 if (fSoftPiCut && IsSoftPion_MCKine(d,trkMC,mcArray)) continue; //remove soft pions (if requestes, e.g. for templates)
1916
4d978a4a 1917 FillSparsePlots(mcArray,mInv,origD0,PDGD0,track,ptbin,kTrack); //fills for charged tracks
1918
1919 //retrieving leading info...
1920 if(track->Pt() > highPt) {
1921 lead[0] = fCorrelatorTr->GetDeltaPhi();
1922 lead[1] = fCorrelatorTr->GetDeltaEta();
1923 lead[2] = fReadMC ? CheckTrackOrigin(mcArray,trkMC) : 0;
1924 highPt = track->Pt();
1925 }
1926
1927 } // end of tracks loop
1928 } //end of event loop for fCorrelatorTr
1929
dea90a65 1930 if(fKaonCorr) { //loops for Kcharg and K0
1931
4d978a4a 1932 if(fMixing) {
1933 NofEventsinPool = fCorrelatorKc->GetNofEventsInPool();
1934 if(!execPoolKc) {
1935 AliInfo("Mixed event analysis: K+/- pool is not ready");
1936 NofEventsinPool = 0;
1937 }
1938 }
1939
1940 //Charged Kaons loop
1941 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)
1942 Bool_t analyzetracksKc = fCorrelatorKc->ProcessAssociatedTracks(jMix);
1943 if(!analyzetracksKc) {
1944 AliInfo("AliHFCorrelator::Cannot process the K+/- array");
1945 continue;
1946 }
1947
1948 for(Int_t iTrack = 0; iTrack<fCorrelatorKc->GetNofTracks(); iTrack++){ // looping on charged kaons candidates
1949
1950 Bool_t runcorrelation = fCorrelatorKc->Correlate(iTrack);
1951 if(!runcorrelation) continue;
1952
1953 AliReducedParticle* kCharg = fCorrelatorKc->GetAssociatedParticle();
1954 if(kCharg->GetLabel()<1) continue;
1955 if(kCharg->Pt() < fPtThreshLow.at(ptbin) || kCharg->Pt() > fPtThreshUp.at(ptbin)) continue; //discard tracks outside pt range for hadrons/K
1956 if(TMath::Abs(kCharg->Eta())>0.8) continue; //discard tracks outside barrel (since it's kinematic MC and produces tracks all over rapidity region
1957 if(!fMixing) N_KCharg++;
1958
1959 AliAODMCParticle *kChargMC = (AliAODMCParticle*)mcArray->At(kCharg->GetLabel());
1960 if(!kChargMC) continue;
1961
1962 if (IsDDaughter(d,kChargMC,mcArray)) continue;
1963 FillSparsePlots(mcArray,mInv,origD0,PDGD0,kCharg,ptbin,kKCharg); //fills for charged tracks
1964
1965 } // end of charged kaons loop
1966 } //end of event loop for fCorrelatorKc
1967
1968 if(fMixing) {
1969 NofEventsinPool = fCorrelatorK0->GetNofEventsInPool();
1970 if(!execPoolK0) {
1971 AliInfo("Mixed event analysis: K0 pool is not ready");
1972 NofEventsinPool = 0;
1973 }
1974 }
1975
1976 //K0 loop
1977 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)
1978 Bool_t analyzetracksK0 = fCorrelatorK0->ProcessAssociatedTracks(jMix);
1979 if(!analyzetracksK0) {
1980 AliInfo("AliHFCorrelator::Cannot process the K0 array");
1981 continue;
1982 }
1983
1984 for(Int_t iTrack = 0; iTrack<fCorrelatorK0->GetNofTracks(); iTrack++){ // looping on k0 candidates
1985
1986 Bool_t runcorrelation = fCorrelatorK0->Correlate(iTrack);
1987 if(!runcorrelation) continue;
1988
1989 AliReducedParticle* k0 = fCorrelatorK0->GetAssociatedParticle();
1990 if(k0->GetLabel()<1) continue;
1991 if(k0->Pt() < fPtThreshLow.at(ptbin) || k0->Pt() > fPtThreshUp.at(ptbin)) continue; //discard tracks outside pt range for hadrons/K
1992 if(TMath::Abs(k0->Eta())>0.8) continue; //discard tracks outside barrel (since it's kinematic MC and produces tracks all over rapidity region
1993
1994 AliAODMCParticle *k0MC = (AliAODMCParticle*)mcArray->At(k0->GetLabel());
1995 if(!k0MC) continue;
1996
1997 if (IsDDaughter(d,k0MC,mcArray)) continue;
1998 FillSparsePlots(mcArray,mInv,origD0,PDGD0,k0,ptbin,kK0); //fills for charged tracks
1999
2000 if(!fMixing) N_Kaons++;
2001
2002 } // end of charged kaons loop
2003 } //end of event loop for fCorrelatorK0
2004
dea90a65 2005 } //end of 'if(fKaonCorr)'
2006
241994ab 2007 Double_t fillSpLeadMC[4] = {lead[0],mD0,lead[1],0.4}; //mD0 = mD0bar = 1.864
4d978a4a 2008
2009 //leading track correlations fill
2010 if(!fMixing) {
241994ab 2011 if(d->GetPdgCode()==421 && (fIsSelectedCandidate==1||fIsSelectedCandidate==3)) { //D0
4d978a4a 2012 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(fillSpLeadMC); //c and b D0
2013 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadMC); //c or b D0
2014 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);
2015 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);
2016 if((int)lead[2]==0) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_NonHF_Bin%d",ptbin)))->Fill(fillSpLeadMC); //non HF
2017 }
241994ab 2018 if(d->GetPdgCode()==-421 && fIsSelectedCandidate>1) { //D0bar
4d978a4a 2019 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(fillSpLeadMC);
2020 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadMC); //c or b D0
2021 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);
2022 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);
2023 if((int)lead[2]==0) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_NonHF_Bin%d",ptbin)))->Fill(fillSpLeadMC); //non HF
2024 }
2025
2026 //Fill of count histograms
2027 if (!fAlreadyFilled) {
2028 ((TH1F*)fOutputStudy->FindObject("hist_Count_Charg"))->Fill(N_Charg);
2029 ((TH1F*)fOutputStudy->FindObject("hist_Count_Kcharg"))->Fill(N_KCharg);
2030 ((TH1F*)fOutputStudy->FindObject("hist_Count_K0"))->Fill(N_Kaons);
2031 }
2032
2033 fAlreadyFilled=kTRUE; //at least a D0 analyzed in the event; distribution plots already filled
2034 }
2035
2036}
2037
2038//________________________________________________________________________
2039void 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 2040 //
2041 //fills the THnSparse for correlations, calculating the variables
2042 //
2043
2044 //Initialization of variables
bce70c96 2045 Double_t mD0, mD0bar, deltaphi = 0., deltaeta = 0.;
4d978a4a 2046 mD0 = mInv[0];
2047 mD0bar = mInv[1];
bce70c96 2048
da975030 2049 if (fReadMC && track->GetLabel()<1) return;
2050 if (fReadMC && !(AliAODMCParticle*)mcArray->At(track->GetLabel())) return;
bce70c96 2051 Double_t ptTrack = track->Pt();
2052 Double_t d0Track = type!=kK0 ? track->GetImpPar() : 0.;
2053 Double_t phiTr = track->Phi();
2054 Double_t origTr = fReadMC ? CheckTrackOrigin(mcArray,(AliAODMCParticle*)mcArray->At(track->GetLabel())) : 0;
2055
2056 TString part = "", orig = "";
2057
2058 switch (type) {
2059 case(kTrack): {
2060 part = "Charg";
2061 deltaphi = fCorrelatorTr->GetDeltaPhi();
2062 deltaeta = fCorrelatorTr->GetDeltaEta();
2063 break;
2064 }
2065 case(kKCharg): {
2066 part = "Kcharg";
2067 deltaphi = fCorrelatorKc->GetDeltaPhi();
2068 deltaeta = fCorrelatorKc->GetDeltaEta();
2069 break;
2070 }
2071 case(kK0): {
2072 part = "K0";
2073 deltaphi = fCorrelatorK0->GetDeltaPhi();
2074 deltaeta = fCorrelatorK0->GetDeltaEta();
2075 break;
a2ad7da1 2076 }
bce70c96 2077 }
bce70c96 2078
2079 if(fMixing == kSE) {
2080
ad8d2dfd 2081 //Fixes limits; needed to include overflow into THnSparse projections!
2082 Double_t pTorig = track->Pt();
2083 Double_t d0orig = track->GetImpPar();
4d978a4a 2084 Double_t ptLim_Sparse = ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->GetAxis(2)->GetXmax(); //all plots have same axes...
2085 Double_t displLim_Sparse = ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->GetAxis(3)->GetXmax();
2086 Double_t EtaLim_Sparse = ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->GetAxis(4)->GetXmax();
ad8d2dfd 2087 if(ptTrack > ptLim_Sparse) ptTrack = ptLim_Sparse-0.01;
2088 if(d0Track > displLim_Sparse) d0Track = (displLim_Sparse-0.001);
2089 if(deltaeta > EtaLim_Sparse) deltaeta = EtaLim_Sparse-0.01;
2090 if(deltaeta < -EtaLim_Sparse) deltaeta = -EtaLim_Sparse+0.01;
2091
2092 //variables for filling histos
2093 Double_t fillSpPhiD0[5] = {deltaphi,mD0,ptTrack,d0Track,deltaeta};
2094 Double_t fillSpPhiD0bar[5] = {deltaphi,mD0bar,ptTrack,d0Track,deltaeta};
241994ab 2095 Double_t fillSpWeigD0[4] = {deltaphi,mD0,deltaeta,ptTrack};
2096 Double_t fillSpWeigD0bar[4] = {deltaphi,mD0bar,deltaeta,ptTrack};
a2ad7da1 2097
bce70c96 2098 if(fReadMC == 0) {
2099 //sparse fill for data (tracks, K+-, K0) + weighted
2100 if(fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) { //D0
4d978a4a 2101 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0,wg);
2102 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig_Bin%d",ptbin)))->Fill(fillSpWeigD0,pTorig*wg);
a2ad7da1 2103 }
bce70c96 2104 if(fIsSelectedCandidate == 2 || fIsSelectedCandidate == 3) { //D0bar
4d978a4a 2105 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0bar,wg);
2106 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig_Bin%d",ptbin)))->Fill(fillSpWeigD0bar,pTorig*wg);
a2ad7da1 2107 }
2108 if(!fAlreadyFilled) {
bce70c96 2109 ((TH1F*)fOutputStudy->FindObject(Form("hist_Pt_%s_Bin%d",part.Data(),ptbin)))->Fill(pTorig);
2110 ((TH1F*)fOutputStudy->FindObject(Form("hist_PhiDistr_%s",part.Data())))->Fill(phiTr);
a2ad7da1 2111 }
bce70c96 2112 }
a2ad7da1 2113
bce70c96 2114 if(fReadMC) {
a2ad7da1 2115
bce70c96 2116 if(origD0==4) {orig = "_From_c";} else {orig = "_From_b";}
a2ad7da1 2117
bce70c96 2118 //sparse fill for data (tracks, K+-, K0) + weighted
241994ab 2119 if(PdgD0==421 && (fIsSelectedCandidate==1||fIsSelectedCandidate==3)) { //D0 (from MCTruth)
4d978a4a 2120 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0,wg);
2121 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0,wg);
2122 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);
2123 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);
2124 if(origTr==0) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_NonHF_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0,wg);
2125 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig_Bin%d",ptbin)))->Fill(fillSpWeigD0,pTorig*wg);
2126 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpWeigD0,pTorig*wg);
2127 if(origD0==4&&origTr>=1&&origTr<=3) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpWeigD0,pTorig*wg);
2128 if(origD0==5&&origTr>=4&&origTr<=8) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpWeigD0,pTorig*wg);
2129 if(origTr==0) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig_NonHF_Bin%d",ptbin)))->Fill(fillSpWeigD0,pTorig*wg);
a2ad7da1 2130 }
241994ab 2131 if(PdgD0==-421 && fIsSelectedCandidate>1) { //D0bar
4d978a4a 2132 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0bar,wg);
2133 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0bar,wg);
2134 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);
2135 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);
2136 if(origTr==0) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_NonHF_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0bar,wg);
2137 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig_Bin%d",ptbin)))->Fill(fillSpWeigD0bar,pTorig*wg);
2138 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpWeigD0bar,pTorig*wg);
2139 if(origD0==4&&origTr>=1&&origTr<=3) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpWeigD0bar,pTorig*wg);
2140 if(origD0==5&&origTr>=4&&origTr<=8) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpWeigD0bar,pTorig*wg);
2141 if(origTr==0) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig_NonHF_Bin%d",ptbin)))->Fill(fillSpWeigD0bar,pTorig*wg);
bce70c96 2142 }
a2ad7da1 2143 if(!fAlreadyFilled) {
bce70c96 2144 ((TH1F*)fOutputStudy->FindObject(Form("histDispl_%s_Bin%d",part.Data(),ptbin)))->Fill(d0orig); //Fills displacement histos
4d978a4a 2145 if (origTr>=1&&origTr<=8) ((TH1F*)fOutputStudy->FindObject(Form("histDispl_%s_HF_Bin%d",part.Data(),ptbin)))->Fill(d0orig);
2146 if (origTr>=1&&origTr<=8) ((TH1F*)fOutputStudy->FindObject(Form("histDispl_%s_HF%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(d0orig);
bce70c96 2147 ((TH1F*)fOutputStudy->FindObject(Form("histDispl_%s%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(d0orig); //Fills displacement histos
2148 ((TH1F*)fOutputStudy->FindObject(Form("hist_Pt_%s_Bin%d",part.Data(),ptbin)))->Fill(pTorig);
2149 ((TH1F*)fOutputStudy->FindObject(Form("histOrig_%s_Bin%d",part.Data(),ptbin)))->Fill(origTr);
2150 ((TH1F*)fOutputStudy->FindObject(Form("hist_PhiDistr_%s",part.Data())))->Fill(phiTr);
a2ad7da1 2151 }
bce70c96 2152 }//end MC case
a2ad7da1 2153
bce70c96 2154 } //end of SE fill
a2ad7da1 2155
bce70c96 2156 if(fMixing == kME) {
a2ad7da1 2157
ad8d2dfd 2158 //Fixes limits; needed to include overflow into THnSparse projections!
4d978a4a 2159 Double_t EtaLim_Sparse = ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d_EvMix",ptbin)))->GetAxis(2)->GetXmax();
ad8d2dfd 2160 if(deltaeta > EtaLim_Sparse) deltaeta = EtaLim_Sparse-0.01;
2161 if(deltaeta < -EtaLim_Sparse) deltaeta = -EtaLim_Sparse+0.01;
241994ab 2162 Double_t ptLim_Sparse = ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d_EvMix",ptbin)))->GetAxis(3)->GetXmax(); //all plots have same axes...
2163 if(ptTrack > ptLim_Sparse) ptTrack = ptLim_Sparse-0.01;
ad8d2dfd 2164
2165 //variables for filling histos
241994ab 2166 Double_t fillSpPhiD0[4] = {deltaphi,mD0,deltaeta,0.4}; //dummy for ME threshold! unless explicitly set by flag...
2167 Double_t fillSpPhiD0bar[4] = {deltaphi,mD0bar,deltaeta,0.4};
2168 if(fMEAxisThresh) {
2169 fillSpPhiD0[3] = ptTrack;
2170 fillSpPhiD0bar[3] = ptTrack;
2171 }
a2ad7da1 2172
bce70c96 2173 if(fReadMC == 0) {
2174 //sparse fill for data (tracks, K+-, K0)
4d978a4a 2175 if(fIsSelectedCandidate == 1||fIsSelectedCandidate == 3) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d_EvMix",part.Data(),ptbin)))->Fill(fillSpPhiD0,wg);
2176 if(fIsSelectedCandidate == 2||fIsSelectedCandidate == 3) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d_EvMix",part.Data(),ptbin)))->Fill(fillSpPhiD0bar,wg);
bce70c96 2177 }
2178 if(fReadMC == 1) {
2179 //sparse fill for data (tracks, K+-, K0)
241994ab 2180 if(PdgD0==421 && (fIsSelectedCandidate==1||fIsSelectedCandidate==3)) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d_EvMix",part.Data(),ptbin)))->Fill(fillSpPhiD0,wg);
2181 if(PdgD0==-421 && fIsSelectedCandidate>1) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d_EvMix",part.Data(),ptbin)))->Fill(fillSpPhiD0bar,wg);
bce70c96 2182 }//end MC case
a2ad7da1 2183
bce70c96 2184 } //end of ME fill
2185
2186 return;
a2ad7da1 2187}
2188
2189//_________________________________________________________________________________________________
2190Int_t AliAnalysisTaskSED0Correlations::CheckTrackOrigin(TClonesArray* arrayMC, AliAODMCParticle *mcPartCandidate) const {
2191 //
2192 // checks on particle (#) origin:
2193 // 0) Not HF
2194 // 1) D->#
2195 // 2) D->X->#
4d978a4a 2196 // 3) c hadronization
2197 // 4) B->#
2198 // 5) B->X-># (X!=D)
2199 // 6) B->D->#
2200 // 7) B->D->X->#
a2ad7da1 2201 // 8) b hadronization
2202 //
bce70c96 2203 if(fDebug>2) printf("AliAnalysisTaskSED0Correlations::CheckTrkOrigin() \n");
a2ad7da1 2204
2205 Int_t pdgGranma = 0;
2206 Int_t mother = 0;
2207 mother = mcPartCandidate->GetMother();
2208 Int_t istep = 0;
2209 Int_t abspdgGranma =0;
2210 Bool_t isFromB=kFALSE;
2211 Bool_t isDdaugh=kFALSE;
2212 Bool_t isDchaindaugh=kFALSE;
2213 Bool_t isBdaugh=kFALSE;
2214 Bool_t isBchaindaugh=kFALSE;
2215 Bool_t isQuarkFound=kFALSE;
2216
4d978a4a 2217 if (mother<0) return -1;
2218 while (mother >= 0){
a2ad7da1 2219 istep++;
4d978a4a 2220 AliAODMCParticle* mcMoth = dynamic_cast<AliAODMCParticle*>(arrayMC->At(mother));
2221 if (mcMoth){
2222 pdgGranma = mcMoth->GetPdgCode();
a2ad7da1 2223 abspdgGranma = TMath::Abs(pdgGranma);
2224 if ((abspdgGranma > 500 && abspdgGranma < 600) || (abspdgGranma > 5000 && abspdgGranma < 6000)){
2225 isBchaindaugh=kTRUE;
2226 if(istep==1) isBdaugh=kTRUE;
2227 }
2228 if ((abspdgGranma > 400 && abspdgGranma < 500) || (abspdgGranma > 4000 && abspdgGranma < 5000)){
2229 isDchaindaugh=kTRUE;
2230 if(istep==1) isDdaugh=kTRUE;
2231 }
2232 if(abspdgGranma==4 || abspdgGranma==5) {isQuarkFound=kTRUE; if(abspdgGranma==5) isFromB = kTRUE;}
4d978a4a 2233 mother = mcMoth->GetMother();
a2ad7da1 2234 }else{
2235 AliError("Failed casting the mother particle!");
4d978a4a 2236 return -1;
a2ad7da1 2237 }
2238 }
2239
2240 //decides what to return based on the flag status
2241 if(isQuarkFound) {
2242 if(!isFromB) { //charm
2243 if(isDdaugh) return 1; //charm immediate
2244 else if(isDchaindaugh) return 2; //charm chain
4d978a4a 2245 else return 3; //charm hadronization
a2ad7da1 2246 }
2247 else { //beauty
4d978a4a 2248 if(isBdaugh) return 4; //b immediate
a2ad7da1 2249 else if(isBchaindaugh) { //b chain
2250 if(isDchaindaugh) {
4d978a4a 2251 if(isDdaugh) return 6; //d immediate
2252 return 7; //d chain
a2ad7da1 2253 }
4d978a4a 2254 else return 5; //b, not d
a2ad7da1 2255 }
2256 else return 8; //b hadronization
2257 }
2258 }
4d978a4a 2259 else if(!isDdaugh && !isDchaindaugh && !isBdaugh && !isBchaindaugh) return 0; //no HF decay
2260 //in this case, it's !isQuarkFound, but not in 100% cases it's a non HF particle!
2261 //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...
2262
2263 return -1; //some problem spotted
2264}
2265//________________________________________________________________________
2266Bool_t AliAnalysisTaskSED0Correlations::IsDDaughter(AliAODMCParticle* d, AliAODMCParticle* track, TClonesArray* mcArray) const {
2267
2268 //Daughter removal in MCKine case
2269 Bool_t isDaughter = kFALSE;
2270 Int_t labelD0 = d->GetLabel();
2271
2272 Int_t mother = track->GetMother();
2273
2274 //Loop on the mothers to find the D0 label (it must be the trigger D0, not a generic D0!)
2275 while (mother > 0){
2276 AliAODMCParticle* mcMoth = dynamic_cast<AliAODMCParticle*>(mcArray->At(mother)); //it's the mother of the track!
2277 if (mcMoth){
2278 if (mcMoth->GetLabel() == labelD0) isDaughter = kTRUE;
2279 mother = mcMoth->GetMother(); //goes back by one
2280 } else{
2281 AliError("Failed casting the mother particle!");
2282 break;
2283 }
2284 }
2285
2286 return isDaughter;
a2ad7da1 2287}
2288
a2ad7da1 2289//________________________________________________________________________
2290Int_t AliAnalysisTaskSED0Correlations::PtBinCorr(Double_t pt) const {
2291 //
2292 //give the pt bin where the pt lies.
2293 //
2294 Int_t ptbin=-1;
2295 if(pt<fBinLimsCorr.at(0)) return ptbin; //out of bounds
2296
2297 Int_t i = 0;
2298 while(pt>fBinLimsCorr.at(i)) {ptbin=i; i++;}
2299
2300 return ptbin;
2301}
2302
a2ad7da1 2303//---------------------------------------------------------------------------
2304Bool_t AliAnalysisTaskSED0Correlations::SelectV0(AliAODv0* v0, AliAODVertex *vtx, Int_t opt, Int_t idArrayV0[][2]) const
2305{
2306 //
2307 // Selection for K0 hypotheses
2308 // options: 1 = selects mass invariant about 3 sigma inside the peak + threshold of 0.3 GeV
2309 // 2 = no previous selections
2310
2311 if(!fCutsTracks->IsKZeroSelected(v0,vtx)) return kFALSE;
2312
2313 AliAODTrack *v0Daug1 = (AliAODTrack*)v0->GetDaughter(0);
2314 AliAODTrack *v0Daug2 = (AliAODTrack*)v0->GetDaughter(1);
2315
2316 if(opt==1) { //additional cuts for correlations (V0 has to be closer than 3 sigma from K0 mass)
bce70c96 2317 if(TMath::Abs(v0->MassK0Short()-0.4976) > 3*0.004) return kFALSE;
a2ad7da1 2318 }
2319
2320 //This part removes double counting for swapped tracks!
2321 Int_t i = 0; //while loop (until the last-written entry pair of ID!
2322 while(idArrayV0[i][0]!=-2 && idArrayV0[i][1]!=-2) {
2323 if((v0Daug1->GetID()==idArrayV0[i][0] && v0Daug2->GetID()==idArrayV0[i][1])||
2324 (v0Daug1->GetID()==idArrayV0[i][1] && v0Daug2->GetID()==idArrayV0[i][0])) return kFALSE;
2325 i++;
2326 }
2327 idArrayV0[i][0]=v0Daug1->GetID();
2328 idArrayV0[i][1]=v0Daug2->GetID();
2329
2330 return kTRUE;
2331}
2332
241994ab 2333//---------------------------------------------------------------------------
2334Bool_t AliAnalysisTaskSED0Correlations::IsSoftPion_MCKine(AliAODMCParticle* d, AliAODMCParticle* track, TClonesArray* arrayMC) const
2335{
2336 //
2337 // Removes soft pions in Kine
2338
2339 //Daughter removal in MCKine case
2340 Bool_t isSoftPi = kFALSE;
2341 Int_t labelD0 = d->GetLabel();
2342
2343 Int_t mother = track->GetMother();
38c832ca 2344 if(mother<0) return isSoftPi; //safety check
241994ab 2345
2346 AliAODMCParticle* mcMoth = dynamic_cast<AliAODMCParticle*>(arrayMC->At(mother)); //it's the mother of the track!
38c832ca 2347 if(!mcMoth){
2348 return isSoftPi;
2349 }
241994ab 2350 if(TMath::Abs(mcMoth->GetPdgCode())==413 && mcMoth->GetNDaughters()==2) { //mother is D* with 2 daughs
2351 Int_t labdau1 = mcMoth->GetDaughter(0);
2352 Int_t labdau2 = mcMoth->GetDaughter(1);
2353 AliAODMCParticle* dau1 = dynamic_cast<AliAODMCParticle*>(arrayMC->At(labdau1));
2354 AliAODMCParticle* dau2 = dynamic_cast<AliAODMCParticle*>(arrayMC->At(labdau2));
2355 if(!dau1 || !dau2) return isSoftPi; //safety check
2356 if(dau1->GetLabel()==labelD0 || dau2->GetLabel()==labelD0) { //one of the daughs is the D0 trigger
2357 if((TMath::Abs(dau1->GetPdgCode())==421 && TMath::Abs(dau2->GetPdgCode())==211)||(TMath::Abs(dau1->GetPdgCode())==211 && TMath::Abs(dau2->GetPdgCode())==421)) {
2358 isSoftPi = kTRUE; //ok, soft pion was found
2359 return isSoftPi;
2360 }
2361 }
2362 }
2363
2364 return isSoftPi;
2365}
2366
a2ad7da1 2367//________________________________________________________________________
2368void AliAnalysisTaskSED0Correlations::PrintBinsAndLimits() {
2369
2370 cout << "--------------------------\n";
2371 cout << "PtBins = " << fNPtBinsCorr << "\n";
2372 cout << "PtBin limits--------------\n";
2373 for (int i=0; i<fNPtBinsCorr; i++) {
2374 cout << "Bin "<<i+1<<" = "<<fBinLimsCorr.at(i)<<" to "<<fBinLimsCorr.at(i+1)<<"\n";
2375 }
2376 cout << "\n--------------------------\n";
2377 cout << "PtBin tresh. tracks low---\n";
2378 for (int i=0; i<fNPtBinsCorr; i++) {
2379 cout << fPtThreshLow.at(i) << "; ";
2380 }
2381 cout << "PtBin tresh. tracks up----\n";
2382 for (int i=0; i<fNPtBinsCorr; i++) {
2383 cout << fPtThreshUp.at(i) << "; ";
2384 }
2385 cout << "\n--------------------------\n";
8c2d7467 2386 cout << "D0 Eta cut for Correl = "<<fEtaForCorrel<<"\n";
2387 cout << "--------------------------\n";
4d978a4a 2388 cout << "MC Truth = "<<fReadMC<<" - MC Reco: Trk = "<<fRecoTr<<", D0 = "<<fRecoD0<<"\n";
2389 cout << "--------------------------\n";
2390 cout << "Sel of Event tpye (PP,GS,FE,...)= "<<fSelEvType<<"\n";
a2ad7da1 2391 cout << "--------------------------\n";
bce70c96 2392 cout << "Ev Mixing = "<<fMixing<<"\n";
2393 cout << "--------------------------\n";
241994ab 2394 cout << "ME thresh axis = "<<fMEAxisThresh<<"\n";
2395 cout << "--------------------------\n";
2396 cout << "Soft Pi Cut = "<<fSoftPiCut<<"\n";
2397 cout << "--------------------------\n";
a2ad7da1 2398}
2399