Updated treatment of TOF PID in QA task (Francesco+Pietro)
[u/mrichter/AliRoot.git] / PWGJE / AliAnalysisTaskJetChem.cxx
CommitLineData
96c271c1 1/*************************************************************************
2 * *
3 * *
4 * Task for Jet Chemistry Analysis in PWG4 Jet Task Force Train *
5 * *
6 * *
7 *************************************************************************/
8
9b2de807 9/**************************************************************************
10 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
11 * *
12 * Author: The ALICE Off-line Project. *
13 * Contributors are mentioned in the code where appropriate. *
14 * *
15 * Permission to use, copy, modify and distribute this software and its *
16 * documentation strictly for non-commercial purposes is hereby granted *
17 * without fee, provided that the above copyright notice appears in all *
18 * copies and that both the copyright notice and this permission notice *
19 * appear in the supporting documentation. The authors make no claims *
20 * about the suitability of this software for any purpose. It is *
21 * provided "as is" without express or implied warranty. *
22 **************************************************************************/
23
96c271c1 24/* $Id: */
9b2de807 25
96c271c1 26#include "TH2F.h"
27#include "TH3F.h"
28#include "TProfile.h"
29#include "THnSparse.h"
9b2de807 30
96c271c1 31#include "AliAnalysisHelperJetTasks.h"
9b2de807 32#include "AliAnalysisManager.h"
96c271c1 33#include "AliAODHandler.h"
34#include "AliAODInputHandler.h"
35#include "AliESDEvent.h"
9b2de807 36#include "AliGenPythiaEventHeader.h"
96c271c1 37#include "AliGenHijingEventHeader.h"
9280dfa6 38
96c271c1 39#include "AliAODEvent.h"
40#include "AliAODJet.h"
41#include "AliAODv0.h"
42#include "AliAODTrack.h"
9b2de807 43
44
96c271c1 45#include "AliPID.h"
46#include "AliExternalTrackParam.h"
9b2de807 47
96c271c1 48#include "AliAnalysisTaskJetChem.h"
9b2de807 49
50ClassImp(AliAnalysisTaskJetChem)
51
96c271c1 52//____________________________________________________________________________
53AliAnalysisTaskJetChem::AliAnalysisTaskJetChem()
54 : AliAnalysisTaskFragmentationFunction()
55 ,fK0Type(0)
56 ,fFilterMaskK0(0)
57 ,fListK0s(0)
58 ,fV0QAK0(0)
59 ,fFFHistosRecCutsK0Evt(0)
60 ,fFFHistosIMK0AllEvt(0)
61 ,fFFHistosIMK0Jet(0)
62 ,fFFHistosIMK0Cone(0)
63 ,fFFHistosPhiCorrIMK0(0)
64 ,fFFIMNBinsJetPt(0)
65 ,fFFIMJetPtMin(0)
66 ,fFFIMJetPtMax(0)
67 ,fFFIMNBinsInvM(0)
68 ,fFFIMInvMMin(0)
69 ,fFFIMInvMMax(0)
70 ,fFFIMNBinsPt(0)
71 ,fFFIMPtMin(0)
72 ,fFFIMPtMax(0)
73 ,fFFIMNBinsXi(0)
74 ,fFFIMXiMin(0)
75 ,fFFIMXiMax(0)
76 ,fFFIMNBinsZ(0)
77 ,fFFIMZMin(0)
78 ,fFFIMZMax(0)
79 ,fPhiCorrIMNBinsPt(0)
80 ,fPhiCorrIMPtMin(0)
81 ,fPhiCorrIMPtMax(0)
82 ,fPhiCorrIMNBinsPhi(0)
83 ,fPhiCorrIMPhiMin(0)
84 ,fPhiCorrIMPhiMax(0)
85 ,fPhiCorrIMNBinsInvM(0)
86 ,fPhiCorrIMInvMMin(0)
87 ,fPhiCorrIMInvMMax(0)
88 ,fh1K0Mult(0)
89 ,fh1dPhiJetK0(0)
9b2de807 90{
96c271c1 91 // default constructor
9b2de807 92}
93
96c271c1 94//__________________________________________________________________________________________
95AliAnalysisTaskJetChem::AliAnalysisTaskJetChem(const char *name)
96 : AliAnalysisTaskFragmentationFunction(name)
97 ,fK0Type(0)
98 ,fFilterMaskK0(0)
99 ,fListK0s(0)
100 ,fV0QAK0(0)
101 ,fFFHistosRecCutsK0Evt(0)
102 ,fFFHistosIMK0AllEvt(0)
103 ,fFFHistosIMK0Jet(0)
104 ,fFFHistosIMK0Cone(0)
105 ,fFFHistosPhiCorrIMK0(0)
106 ,fFFIMNBinsJetPt(0)
107 ,fFFIMJetPtMin(0)
108 ,fFFIMJetPtMax(0)
109 ,fFFIMNBinsInvM(0)
110 ,fFFIMInvMMin(0)
111 ,fFFIMInvMMax(0)
112 ,fFFIMNBinsPt(0)
113 ,fFFIMPtMin(0)
114 ,fFFIMPtMax(0)
115 ,fFFIMNBinsXi(0)
116 ,fFFIMXiMin(0)
117 ,fFFIMXiMax(0)
118 ,fFFIMNBinsZ(0)
119 ,fFFIMZMin(0)
120 ,fFFIMZMax(0)
121 ,fPhiCorrIMNBinsPt(0)
122 ,fPhiCorrIMPtMin(0)
123 ,fPhiCorrIMPtMax(0)
124 ,fPhiCorrIMNBinsPhi(0)
125 ,fPhiCorrIMPhiMin(0)
126 ,fPhiCorrIMPhiMax(0)
127 ,fPhiCorrIMNBinsInvM(0)
128 ,fPhiCorrIMInvMMin(0)
129 ,fPhiCorrIMInvMMax(0)
130 ,fh1K0Mult(0)
131 ,fh1dPhiJetK0(0)
9b2de807 132{
96c271c1 133 // constructor
134
135 DefineOutput(1,TList::Class());
9b2de807 136}
137
96c271c1 138//__________________________________________________________________________________________________________________________
139AliAnalysisTaskJetChem::AliAnalysisTaskJetChem(const AliAnalysisTaskJetChem &copy)
140 : AliAnalysisTaskFragmentationFunction()
141 ,fK0Type(copy.fK0Type)
142 ,fFilterMaskK0(copy.fFilterMaskK0)
143 ,fListK0s(copy.fListK0s)
144 ,fV0QAK0(copy.fV0QAK0)
145 ,fFFHistosRecCutsK0Evt(copy.fFFHistosRecCutsK0Evt)
146 ,fFFHistosIMK0AllEvt(copy.fFFHistosIMK0AllEvt)
147 ,fFFHistosIMK0Jet(copy.fFFHistosIMK0Jet)
148 ,fFFHistosIMK0Cone(copy.fFFHistosIMK0Cone)
149 ,fFFHistosPhiCorrIMK0(copy.fFFHistosPhiCorrIMK0)
150 ,fFFIMNBinsJetPt(copy.fFFIMNBinsJetPt)
151 ,fFFIMJetPtMin(copy.fFFIMJetPtMin)
152 ,fFFIMJetPtMax(copy.fFFIMJetPtMax)
153 ,fFFIMNBinsInvM(copy.fFFIMNBinsInvM)
154 ,fFFIMInvMMin(copy.fFFIMInvMMin)
155 ,fFFIMInvMMax(copy.fFFIMInvMMax)
156 ,fFFIMNBinsPt(copy.fFFIMNBinsPt)
157 ,fFFIMPtMin(copy.fFFIMPtMin)
158 ,fFFIMPtMax(copy.fFFIMPtMax)
159 ,fFFIMNBinsXi(copy.fFFIMNBinsXi)
160 ,fFFIMXiMin(copy.fFFIMXiMin)
161 ,fFFIMXiMax(copy.fFFIMXiMax)
162 ,fFFIMNBinsZ(copy.fFFIMNBinsZ)
163 ,fFFIMZMin(copy.fFFIMZMin)
164 ,fFFIMZMax(copy.fFFIMZMax)
165 ,fPhiCorrIMNBinsPt(copy.fPhiCorrIMNBinsPt)
166 ,fPhiCorrIMPtMin(copy.fPhiCorrIMPtMin)
167 ,fPhiCorrIMPtMax(copy.fPhiCorrIMPtMax)
168 ,fPhiCorrIMNBinsPhi(copy.fPhiCorrIMNBinsPhi)
169 ,fPhiCorrIMPhiMin(copy.fPhiCorrIMPhiMin)
170 ,fPhiCorrIMPhiMax(copy.fPhiCorrIMPhiMax)
171 ,fPhiCorrIMNBinsInvM(copy.fPhiCorrIMNBinsInvM)
172 ,fPhiCorrIMInvMMin(copy.fPhiCorrIMInvMMin)
173 ,fPhiCorrIMInvMMax(copy.fPhiCorrIMInvMMax)
174 ,fh1K0Mult(copy.fh1K0Mult)
175 ,fh1dPhiJetK0(copy.fh1dPhiJetK0)
9b2de807 176{
96c271c1 177 // copy constructor
9b2de807 178
96c271c1 179}
9b2de807 180
96c271c1 181// _________________________________________________________________________________________________________________________________
182AliAnalysisTaskJetChem& AliAnalysisTaskJetChem::operator=(const AliAnalysisTaskJetChem& o)
183{
184 // assignment
185
186 if(this!=&o){
187 AliAnalysisTaskFragmentationFunction::operator=(o);
188
189 fK0Type = o.fK0Type;
190 fFilterMaskK0 = o.fFilterMaskK0;
191 fListK0s = o.fListK0s;
192 fV0QAK0 = o.fV0QAK0;
193 fFFHistosRecCutsK0Evt = o.fFFHistosRecCutsK0Evt;
194 fFFHistosIMK0AllEvt = o.fFFHistosIMK0AllEvt;
195 fFFHistosIMK0Jet = o.fFFHistosIMK0Jet;
196 fFFHistosIMK0Cone = o.fFFHistosIMK0Cone;
197 fFFHistosPhiCorrIMK0 = o.fFFHistosPhiCorrIMK0;
198 fFFIMNBinsJetPt = o.fFFIMNBinsJetPt;
199 fFFIMJetPtMin = o.fFFIMJetPtMin;
200 fFFIMJetPtMax = o.fFFIMJetPtMax;
201 fFFIMNBinsPt = o.fFFIMNBinsPt;
202 fFFIMPtMin = o.fFFIMPtMin;
203 fFFIMPtMax = o.fFFIMPtMax;
204 fFFIMNBinsXi = o.fFFIMNBinsXi;
205 fFFIMXiMin = o.fFFIMXiMin;
206 fFFIMXiMax = o.fFFIMXiMax;
207 fFFIMNBinsZ = o.fFFIMNBinsZ;
208 fFFIMZMin = o.fFFIMZMin;
209 fFFIMZMax = o.fFFIMZMax;
210 fPhiCorrIMNBinsPt = o.fPhiCorrIMNBinsPt;
211 fPhiCorrIMPtMin = o.fPhiCorrIMPtMin;
212 fPhiCorrIMPtMax = o.fPhiCorrIMPtMax;
213 fPhiCorrIMNBinsPhi = o.fPhiCorrIMNBinsPhi;
214 fPhiCorrIMPhiMin = o.fPhiCorrIMPhiMin;
215 fPhiCorrIMPhiMax = o.fPhiCorrIMPhiMax;
216 fPhiCorrIMNBinsInvM = o.fPhiCorrIMNBinsInvM;
217 fPhiCorrIMInvMMin = o.fPhiCorrIMInvMMin;
218 fPhiCorrIMInvMMax = o.fPhiCorrIMInvMMax;
219 fh1K0Mult = o.fh1K0Mult;
220 fh1dPhiJetK0 = o.fh1dPhiJetK0;
9b2de807 221 }
9b2de807 222
96c271c1 223 return *this;
9b2de807 224}
225
96c271c1 226//_______________________________________________
227AliAnalysisTaskJetChem::~AliAnalysisTaskJetChem()
9b2de807 228{
96c271c1 229 // destructor
9b2de807 230
96c271c1 231 if(fListK0s) delete fListK0s;
96c271c1 232}
9b2de807 233
96c271c1 234//________________________________________________________________________________________________________________________________
235AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::AliFragFuncHistosInvMass(const char* name,
236 Int_t nJetPt, Float_t jetPtMin, Float_t jetPtMax,
237 Int_t nInvMass, Float_t invMassMin, Float_t invMassMax,
238 Int_t nPt, Float_t ptMin, Float_t ptMax,
239 Int_t nXi, Float_t xiMin, Float_t xiMax,
240 Int_t nZ , Float_t zMin , Float_t zMax )
241 : TObject()
242 ,fNBinsJetPt(nJetPt)
243 ,fJetPtMin(jetPtMin)
244 ,fJetPtMax(jetPtMax)
245 ,fNBinsInvMass(nInvMass)
246 ,fInvMassMin(invMassMin)
247 ,fInvMassMax(invMassMax)
248 ,fNBinsPt(nPt)
249 ,fPtMin(ptMin)
250 ,fPtMax(ptMax)
251 ,fNBinsXi(nXi)
252 ,fXiMin(xiMin)
253 ,fXiMax(xiMax)
254 ,fNBinsZ(nZ)
255 ,fZMin(zMin)
256 ,fZMax(zMax)
257 ,fh3TrackPt(0)
258 ,fh3Xi(0)
259 ,fh3Z(0)
260 ,fh1JetPt(0)
261 ,fNameFF(name)
262{
263 // default constructor
9b2de807 264
96c271c1 265}
9b2de807 266
96c271c1 267//______________________________________________________________________________________________________________
268AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::AliFragFuncHistosInvMass(const AliFragFuncHistosInvMass& copy)
269 : TObject()
270 ,fNBinsJetPt(copy.fNBinsJetPt)
271 ,fJetPtMin(copy.fJetPtMin)
272 ,fJetPtMax(copy.fJetPtMax)
273 ,fNBinsInvMass(copy.fNBinsInvMass)
274 ,fInvMassMin(copy.fInvMassMin)
275 ,fInvMassMax(copy.fInvMassMax)
276 ,fNBinsPt(copy.fNBinsPt)
277 ,fPtMin(copy.fPtMin)
278 ,fPtMax(copy.fPtMax)
279 ,fNBinsXi(copy.fNBinsXi)
280 ,fXiMin(copy.fXiMin)
281 ,fXiMax(copy.fXiMax)
282 ,fNBinsZ(copy.fNBinsZ)
283 ,fZMin(copy.fZMin)
284 ,fZMax(copy.fZMax)
285 ,fh3TrackPt(copy.fh3TrackPt)
286 ,fh3Xi(copy.fh3Xi)
287 ,fh3Z(copy.fh3Z)
288 ,fh1JetPt(copy.fh1JetPt)
289 ,fNameFF(copy.fNameFF)
290{
291 // copy constructor
292}
9280dfa6 293
96c271c1 294//______________________________________________________________________________________________________________________________________________________________________
295AliAnalysisTaskJetChem::AliFragFuncHistosInvMass& AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::operator=(const AliAnalysisTaskJetChem::AliFragFuncHistosInvMass& o)
296{
297 // assignment
298
299 if(this!=&o){
300 TObject::operator=(o);
301 fNBinsJetPt = o.fNBinsJetPt;
302 fJetPtMin = o.fJetPtMin;
303 fJetPtMax = o.fJetPtMax;
304 fNBinsInvMass = o.fNBinsInvMass;
305 fInvMassMin = o.fInvMassMin;
306 fInvMassMax = o.fInvMassMax;
307 fNBinsPt = o.fNBinsPt;
308 fPtMin = o.fPtMin;
309 fPtMax = o.fPtMax;
310 fNBinsXi = o.fNBinsXi;
311 fXiMin = o.fXiMin;
312 fXiMax = o.fXiMax;
313 fNBinsZ = o.fNBinsZ;
314 fZMin = o.fZMin;
315 fZMax = o.fZMax;
316 fh3TrackPt = o.fh3TrackPt;
317 fh3Xi = o.fh3Xi;
318 fh3Z = o.fh3Z;
319 fh1JetPt = o.fh1JetPt;
320 fNameFF = o.fNameFF;
9280dfa6 321 }
9280dfa6 322
96c271c1 323 return *this;
324}
9b2de807 325
96c271c1 326//___________________________________________________________________________
327AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::~AliFragFuncHistosInvMass()
328{
329 // destructor
9b2de807 330
96c271c1 331 if(fh1JetPt) delete fh1JetPt;
332 if(fh3TrackPt) delete fh3TrackPt;
333 if(fh3Xi) delete fh3Xi;
334 if(fh3Z) delete fh3Z;
335}
9b2de807 336
96c271c1 337//_________________________________________________________________
338void AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::DefineHistos()
339{
340 // book FF histos
9b2de807 341
96c271c1 342 fh1JetPt = new TH1F(Form("fh1FFJetPtIM%s", fNameFF.Data()),"",fNBinsJetPt,fJetPtMin,fJetPtMax);
343 fh3TrackPt = new TH3F(Form("fh3FFTrackPtIM%s",fNameFF.Data()),"",fNBinsJetPt, fJetPtMin, fJetPtMax, fNBinsInvMass, fInvMassMin, fInvMassMax, fNBinsPt, fPtMin, fPtMax);
344 fh3Xi = new TH3F(Form("fh3FFXiIM%s", fNameFF.Data()),"", fNBinsJetPt, fJetPtMin, fJetPtMax, fNBinsInvMass, fInvMassMin, fInvMassMax, fNBinsXi, fXiMin, fXiMax);
345 fh3Z = new TH3F(Form("fh3FFZIM%s", fNameFF.Data()),"", fNBinsJetPt, fJetPtMin, fJetPtMax, fNBinsInvMass, fInvMassMin, fInvMassMax, fNBinsZ, fZMin, fZMax);
9b2de807 346
96c271c1 347 AliAnalysisTaskFragmentationFunction::SetProperties(fh1JetPt, "p_{t} (GeV/c)", "entries");
348 AliAnalysisTaskJetChem::SetProperties(fh3TrackPt,"jet p_{t} (GeV/c)","inv Mass (GeV/c^2)","p_{t} (GeV/c)");
349 AliAnalysisTaskJetChem::SetProperties(fh3Xi,"jet p_{t} (GeV/c)","inv Mass (GeV/c^2)","#xi");
350 AliAnalysisTaskJetChem::SetProperties(fh3Z,"jet p_{t} (GeV/c)","inv Mass (GeV/c^2)","z");
351}
9b2de807 352
96c271c1 353//________________________________________________________________________________________________________________________________
354void AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::FillFF(Float_t trackPt, Float_t invM, Float_t jetPt, Bool_t incrementJetPt)
355{
356 // fill FF
9280dfa6 357
96c271c1 358 if(incrementJetPt) fh1JetPt->Fill(jetPt);
359 fh3TrackPt->Fill(jetPt,invM,trackPt);
9b2de807 360
96c271c1 361 Double_t z = 0.;
362 if(jetPt>0) z = trackPt / jetPt;
363 Double_t xi = 0;
364 if(z>0) xi = TMath::Log(1/z);
9280dfa6 365
96c271c1 366 fh3Xi->Fill(jetPt,invM,xi);
367 fh3Z->Fill(jetPt,invM,z);
9b2de807 368}
369
96c271c1 370//___________________________________________________________________________________
371void AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::AddToOutput(TList* list) const
372{
373 // add histos to list
9f77f0ce 374
96c271c1 375 list->Add(fh1JetPt);
9b2de807 376
96c271c1 377 list->Add(fh3TrackPt);
378 list->Add(fh3Xi);
379 list->Add(fh3Z);
380}
9f77f0ce 381
96c271c1 382// ---
383
384
385//_______________________________________________________________________________________________________
386AliAnalysisTaskJetChem::AliFragFuncHistosPhiCorrInvMass::AliFragFuncHistosPhiCorrInvMass(const char* name,
387 Int_t nPt, Float_t ptMin, Float_t ptMax,
388 Int_t nPhi, Float_t phiMin, Float_t phiMax,
389 Int_t nInvMass, Float_t invMassMin, Float_t invMassMax)
390 : TObject()
391 ,fNBinsPt(nPt)
392 ,fPtMin(ptMin)
393 ,fPtMax(ptMax)
394 ,fNBinsPhi(nPhi)
395 ,fPhiMin(phiMin)
396 ,fPhiMax(phiMax)
397 ,fNBinsInvMass(nInvMass)
398 ,fInvMassMin(invMassMin)
399 ,fInvMassMax(invMassMax)
400 ,fh3PhiCorr(0)
401 ,fNamePhiCorr(name)
402{
403 // default constructor
404}
9b2de807 405
96c271c1 406//____________________________________________________________________________________________________________________________________
407AliAnalysisTaskJetChem::AliFragFuncHistosPhiCorrInvMass::AliFragFuncHistosPhiCorrInvMass(const AliFragFuncHistosPhiCorrInvMass& copy)
408 : TObject()
409 ,fNBinsPt(copy.fNBinsPt)
410 ,fPtMin(copy.fPtMin)
411 ,fPtMax(copy.fPtMax)
412 ,fNBinsPhi(copy.fNBinsPhi)
413 ,fPhiMin(copy.fPhiMin)
414 ,fPhiMax(copy.fPhiMax)
415 ,fNBinsInvMass(copy.fNBinsInvMass)
416 ,fInvMassMin(copy.fInvMassMin)
417 ,fInvMassMax(copy.fInvMassMax)
418 ,fh3PhiCorr(copy.fh3PhiCorr)
419 ,fNamePhiCorr(copy.fNamePhiCorr)
420{
421 // copy constructor
9b2de807 422}
423
96c271c1 424//________________________________________________________________________________________________________________________________________________________________________
425AliAnalysisTaskJetChem::AliFragFuncHistosPhiCorrInvMass& AliAnalysisTaskJetChem::AliFragFuncHistosPhiCorrInvMass::operator=(const AliAnalysisTaskJetChem::AliFragFuncHistosPhiCorrInvMass& o)
9b2de807 426{
96c271c1 427 // assignment
428
429 if(this!=&o){
430 TObject::operator=(o);
431 fNBinsPt = o.fNBinsPt;
432 fPtMin = o.fPtMin;
433 fPtMax = o.fPtMax;
434 fNBinsPhi = o.fNBinsPhi;
435 fPhiMin = o.fPhiMin;
436 fPhiMax = o.fPhiMax;
437 fNBinsInvMass = o.fNBinsInvMass;
438 fInvMassMin = o.fInvMassMin;
439 fInvMassMax = o.fInvMassMax;
440
441 fh3PhiCorr = o.fh3PhiCorr;
442 fNamePhiCorr = o.fNamePhiCorr;
443 }
9b2de807 444
96c271c1 445 return *this;
9b2de807 446}
447
96c271c1 448//_________________________________________________________________________________________
449AliAnalysisTaskJetChem::AliFragFuncHistosPhiCorrInvMass::~AliFragFuncHistosPhiCorrInvMass()
9b2de807 450{
96c271c1 451 // destructor
9b2de807 452
96c271c1 453 if(fh3PhiCorr) delete fh3PhiCorr;
9b2de807 454}
455
96c271c1 456//__________________________________________________________________________
457void AliAnalysisTaskJetChem::AliFragFuncHistosPhiCorrInvMass::DefineHistos()
458{
459 // book jet QA histos
9b2de807 460
96c271c1 461 fh3PhiCorr = new TH3F(Form("fh3PhiCorrIM%s", fNamePhiCorr.Data()),
462 Form("%s: p_{t} - #phi - m_{inv} distribution",fNamePhiCorr.Data()),
463 fNBinsPt, fPtMin, fPtMax,
464 fNBinsPhi, fPhiMin, fPhiMax,
465 fNBinsInvMass, fInvMassMin, fInvMassMax);
9b2de807 466
96c271c1 467 AliAnalysisTaskJetChem::SetProperties(fh3PhiCorr, "p_{t} (GeV/c)", "#phi", "m_{inv} (GeV/c^2)");
9b2de807 468}
469
96c271c1 470//___________________________________________________________________________________________________________
471void AliAnalysisTaskJetChem::AliFragFuncHistosPhiCorrInvMass::FillPhiCorr(Float_t pt, Float_t phi, Float_t invM)
9b2de807 472{
96c271c1 473 // fill jet QA histos
9b2de807 474
96c271c1 475 fh3PhiCorr->Fill(pt, phi, invM);
476}
9b2de807 477
96c271c1 478//______________________________________________________________________________________________
479void AliAnalysisTaskJetChem::AliFragFuncHistosPhiCorrInvMass::AddToOutput(TList* list) const
480{
481 // add histos to list
9b2de807 482
96c271c1 483 list->Add(fh3PhiCorr);
484}
9b2de807 485
96c271c1 486//____________________________________________________
487void AliAnalysisTaskJetChem::UserCreateOutputObjects()
488{
489 // create output objects
9b2de807 490
96c271c1 491 if(fDebug > 1) Printf("AliAnalysisTaskFragmentationFunction::UserCreateOutputObjects()");
492
493 // create list of tracks and jets
494
495 fTracksRecCuts = new TList();
496 fTracksRecCuts->SetOwner(kFALSE);
9b2de807 497
96c271c1 498 fJetsRecCuts = new TList();
499 fJetsRecCuts->SetOwner(kFALSE);
9b2de807 500
96c271c1 501 fListK0s = new TList();
502 fListK0s->SetOwner(kFALSE);
9b2de807 503
96c271c1 504 //
505 // Create histograms / output container
506 //
9b2de807 507
96c271c1 508 OpenFile(1);
509 fCommonHistList = new TList();
9b2de807 510
96c271c1 511 Bool_t oldStatus = TH1::AddDirectoryStatus();
512 TH1::AddDirectory(kFALSE);
513
9b2de807 514
96c271c1 515 // histograms inherited from AliAnalysisTaskFragmentationFunction
9b2de807 516
96c271c1 517 fh1EvtSelection = new TH1F("fh1EvtSelection", "Event Selection", 6, -0.5, 5.5);
518 fh1EvtSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
519 fh1EvtSelection->GetXaxis()->SetBinLabel(2,"event selection: rejected");
520 fh1EvtSelection->GetXaxis()->SetBinLabel(3,"event class: rejected");
521 fh1EvtSelection->GetXaxis()->SetBinLabel(4,"vertex Ncontr: rejected");
522 fh1EvtSelection->GetXaxis()->SetBinLabel(5,"vertex z: rejected");
523 fh1EvtSelection->GetXaxis()->SetBinLabel(6,"vertex type: rejected");
9b2de807 524
96c271c1 525 fh1EvtCent = new TH1F("fh1EvtCent","centrality",100,0.,100.);
526 fh1VertexNContributors = new TH1F("fh1VertexNContributors", "Vertex N contributors", 11,-.5, 10.5);
527 fh1VertexZ = new TH1F("fh1VertexZ", "Vertex z distribution", 30, -15., 15.);
528 fh1Xsec = new TProfile("fh1Xsec","xsec from pyxsec.root",1,0,1);
529 fh1Xsec->GetXaxis()->SetBinLabel(1,"<#sigma>");
530 fh1Trials = new TH1F("fh1Trials","trials from pyxsec.root",1,0,1);
531 fh1Trials->GetXaxis()->SetBinLabel(1,"#sum{ntrials}");
532 fh1PtHard = new TH1F("fh1PtHard","PYTHIA Pt hard;p_{T,hard}",350,-.5,349.5);
533 fh1PtHardTrials = new TH1F("fh1PtHardTrials","PYTHIA Pt hard weight with trials;p_{T,hard}",350,-.5,349.5);
9b2de807 534
96c271c1 535 fh1nRecJetsCuts = new TH1F("fh1nRecJetsCuts","reconstructed jets per event",10,-0.5,9.5);
9b2de807 536
9b2de807 537
96c271c1 538 // histograms jetChem proper
9b2de807 539
96c271c1 540 fh1EvtMult = new TH1F("fh1EvtMult","multiplicity",1200,0.,12000.);
541 fh1K0Mult = new TH1F("fh1K0Mult","K0 multiplicity",500,0.,500.);
542 fh1dPhiJetK0 = new TH1F("fh1dPhiJetK0","",640,-1,5.4);
9b2de807 543
9b2de807 544
9b2de807 545
96c271c1 546 fFFHistosRecCuts = new AliFragFuncHistos("RecCuts", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
547 fFFNBinsPt, fFFPtMin, fFFPtMax,
548 fFFNBinsXi, fFFXiMin, fFFXiMax,
549 fFFNBinsZ , fFFZMin , fFFZMax);
9b2de807 550
96c271c1 551 fV0QAK0 = new AliFragFuncQATrackHistos("V0QAK0",fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
552 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
553 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
554 fQATrackHighPtThreshold);
9b2de807 555
96c271c1 556 fFFHistosRecCutsK0Evt = new AliFragFuncHistos("RecCutsK0Evt", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
557 fFFNBinsPt, fFFPtMin, fFFPtMax,
558 fFFNBinsXi, fFFXiMin, fFFXiMax,
559 fFFNBinsZ , fFFZMin , fFFZMax);
560
561
562 fFFHistosIMK0AllEvt = new AliFragFuncHistosInvMass("K0AllEvt", fFFIMNBinsJetPt, fFFIMJetPtMin, fFFIMJetPtMax,
563 fFFIMNBinsInvM,fFFIMInvMMin,fFFIMInvMMax,
564 fFFIMNBinsPt, fFFIMPtMin, fFFIMPtMax,
565 fFFIMNBinsXi, fFFIMXiMin, fFFIMXiMax,
566 fFFIMNBinsZ , fFFIMZMin , fFFIMZMax);
567
568 fFFHistosIMK0Jet = new AliFragFuncHistosInvMass("K0Jet", fFFIMNBinsJetPt, fFFIMJetPtMin, fFFIMJetPtMax,
569 fFFIMNBinsInvM,fFFIMInvMMin,fFFIMInvMMax,
570 fFFIMNBinsPt, fFFIMPtMin, fFFIMPtMax,
571 fFFIMNBinsXi, fFFIMXiMin, fFFIMXiMax,
572 fFFIMNBinsZ , fFFIMZMin , fFFIMZMax);
573
574
575 fFFHistosIMK0Cone = new AliFragFuncHistosInvMass("K0Cone", fFFIMNBinsJetPt, fFFIMJetPtMin, fFFIMJetPtMax,
576 fFFIMNBinsInvM,fFFIMInvMMin,fFFIMInvMMax,
577 fFFIMNBinsPt, fFFIMPtMin, fFFIMPtMax,
578 fFFIMNBinsXi, fFFIMXiMin, fFFIMXiMax,
579 fFFIMNBinsZ , fFFIMZMin , fFFIMZMax);
580
581 fFFHistosPhiCorrIMK0 = new AliFragFuncHistosPhiCorrInvMass("K0",fPhiCorrIMNBinsPt, fPhiCorrIMPtMin, fPhiCorrIMPtMax,
582 fPhiCorrIMNBinsPhi, fPhiCorrIMPhiMin, fPhiCorrIMPhiMax,
583 fPhiCorrIMNBinsInvM , fPhiCorrIMInvMMin , fPhiCorrIMInvMMax);
584
585
586 fV0QAK0->DefineHistos();
587 fFFHistosRecCuts->DefineHistos();
588 fFFHistosRecCutsK0Evt->DefineHistos();
589 fFFHistosIMK0AllEvt->DefineHistos();
590 fFFHistosIMK0Jet->DefineHistos();
591 fFFHistosIMK0Cone->DefineHistos();
592 fFFHistosPhiCorrIMK0->DefineHistos();
9b2de807 593
96c271c1 594 const Int_t saveLevel = 5;
595 if(saveLevel>0){
9b2de807 596
96c271c1 597 fCommonHistList->Add(fh1EvtSelection);
598 fCommonHistList->Add(fh1EvtCent);
599 fCommonHistList->Add(fh1VertexNContributors);
600 fCommonHistList->Add(fh1VertexZ);
601 fCommonHistList->Add(fh1Xsec);
602 fCommonHistList->Add(fh1Trials);
603 fCommonHistList->Add(fh1PtHard);
604 fCommonHistList->Add(fh1PtHardTrials);
605 fCommonHistList->Add(fh1nRecJetsCuts);
606 fCommonHistList->Add(fh1EvtMult);
607 fCommonHistList->Add(fh1K0Mult);
608 fCommonHistList->Add(fh1dPhiJetK0);
609
610 fV0QAK0->AddToOutput(fCommonHistList);
611 fFFHistosRecCuts->AddToOutput(fCommonHistList);
612 fFFHistosRecCutsK0Evt->AddToOutput(fCommonHistList);
613
614 fFFHistosIMK0AllEvt->AddToOutput(fCommonHistList);
615 fFFHistosIMK0Jet->AddToOutput(fCommonHistList);
616 fFFHistosIMK0Cone->AddToOutput(fCommonHistList);
617 fFFHistosPhiCorrIMK0->AddToOutput(fCommonHistList);
618 }
619
620 // =========== Switch on Sumw2 for all histos ===========
621 for (Int_t i=0; i<fCommonHistList->GetEntries(); ++i){
622 TH1 *h1 = dynamic_cast<TH1*>(fCommonHistList->At(i));
623 if (h1) h1->Sumw2();
624 else{
625 THnSparse *hnSparse = dynamic_cast<THnSparse*>(fCommonHistList->At(i));
626 if(hnSparse) hnSparse->Sumw2();
9b2de807 627 }
628 }
96c271c1 629
630 TH1::AddDirectory(oldStatus);
9b2de807 631}
632
96c271c1 633//_______________________________________________
634void AliAnalysisTaskJetChem::UserExec(Option_t *)
9b2de807 635{
96c271c1 636 // Main loop
637 // Called for each event
9b2de807 638
96c271c1 639 if(fDebug > 1) Printf("Analysis event #%5d", (Int_t) fEntry);
9b2de807 640
96c271c1 641 AliInputEventHandler* inputHandler = (AliInputEventHandler*)
642 ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
643 if(inputHandler->IsEventSelected() & AliVEvent::kMB){
644 if(fDebug > 1) Printf(" Trigger Selection: event ACCEPTED ... ");
645 fh1EvtSelection->Fill(1.);
646 } else {
647 fh1EvtSelection->Fill(0.);
648 if(inputHandler->InheritsFrom("AliESDInputHandler") && fUsePhysicsSelection){ // PhysicsSelection only with ESD input
649 if (fDebug > 1 ) Printf(" Trigger Selection: event REJECTED ... ");
650 PostData(1, fCommonHistList);
651 return;
652 }
9b2de807 653 }
654
96c271c1 655 fESD = dynamic_cast<AliESDEvent*>(InputEvent());
656 if(!fESD){
657 if(fDebug>3) Printf("%s:%d ESDEvent not found in the input", (char*)__FILE__,__LINE__);
658 }
9b2de807 659
96c271c1 660 fMCEvent = MCEvent();
661 if(!fMCEvent){
662 if(fDebug>3) Printf("%s:%d MCEvent not found in the input", (char*)__FILE__,__LINE__);
663 }
9b2de807 664
96c271c1 665 // get AOD event from input/ouput
666 TObject* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
667 if( handler && handler->InheritsFrom("AliAODInputHandler") ) {
668 fAOD = ((AliAODInputHandler*)handler)->GetEvent();
669 if(fUseAODInputJets) fAODJets = fAOD;
670 if (fDebug > 1) Printf("%s:%d AOD event from input", (char*)__FILE__,__LINE__);
671 }
672 else {
673 handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
674 if( handler && handler->InheritsFrom("AliAODHandler") ) {
675 fAOD = ((AliAODHandler*)handler)->GetAOD();
676 fAODJets = fAOD;
677 if (fDebug > 1) Printf("%s:%d AOD event from output", (char*)__FILE__,__LINE__);
9b2de807 678 }
9b2de807 679 }
9b2de807 680
96c271c1 681 if(!fAODJets && !fUseAODInputJets){ // case we have AOD in input & output and want jets from output
682 TObject* outHandler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
683 if( outHandler && outHandler->InheritsFrom("AliAODHandler") ) {
684 fAODJets = ((AliAODHandler*)outHandler)->GetAOD();
685 if (fDebug > 1) Printf("%s:%d jets from output AOD", (char*)__FILE__,__LINE__);
686 }
9b2de807 687 }
688
96c271c1 689 if(fNonStdFile.Length()!=0){
690 // case we have an AOD extension - fetch the jets from the extended output
691
692 AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
693 fAODExtension = (aodH?aodH->GetExtension(fNonStdFile.Data()):0);
694 if(!fAODExtension){
695 if(fDebug>1)Printf("AODExtension not found for %s",fNonStdFile.Data());
696 }
9b2de807 697 }
698
96c271c1 699 if(!fAOD){
700 Printf("%s:%d AODEvent not found", (char*)__FILE__,__LINE__);
701 return;
702 }
703 if(!fAODJets){
704 Printf("%s:%d AODEvent with jet branch not found", (char*)__FILE__,__LINE__);
705 return;
9b2de807 706 }
707
9b2de807 708
96c271c1 709 // event selection *****************************************
9b2de807 710
96c271c1 711 Double_t centPercent = -1;
712 if(fEventClass>0){
713 Int_t cl = 0;
c500e0a0 714 if(handler && handler->InheritsFrom("AliAODInputHandler")){
96c271c1 715 // since it is not supported by the helper task define own classes
716 centPercent = fAOD->GetHeader()->GetCentrality();
717 cl = 1;
718 if(centPercent>10) cl = 2;
719 if(centPercent>30) cl = 3;
720 if(centPercent>50) cl = 4;
9b2de807 721 }
96c271c1 722 else {
723 cl = AliAnalysisHelperJetTasks::EventClass();
c500e0a0 724 if(fESD) centPercent = fESD->GetCentrality()->GetCentralityPercentile("V0M"); // OB added
9b2de807 725 }
96c271c1 726
727 if(cl!=fEventClass){
728 // event not in selected event class, reject event
729 if (fDebug > 1) Printf("%s:%d event not in selected event class: event REJECTED ...",(char*)__FILE__,__LINE__);
730 fh1EvtSelection->Fill(2.);
731 PostData(1, fCommonHistList);
732 return;
9b2de807 733 }
734 }
9b2de807 735
96c271c1 736 // *** vertex cut ***
737 AliAODVertex* primVtx = fAOD->GetPrimaryVertex();
738 Int_t nTracksPrim = primVtx->GetNContributors();
739 fh1VertexNContributors->Fill(nTracksPrim);
9b2de807 740
96c271c1 741 if (fDebug > 1) Printf("%s:%d primary vertex selection: %d", (char*)__FILE__,__LINE__,nTracksPrim);
742 if(!nTracksPrim){
743 if (fDebug > 1) Printf("%s:%d primary vertex selection: event REJECTED...",(char*)__FILE__,__LINE__);
744 fh1EvtSelection->Fill(3.);
745 PostData(1, fCommonHistList);
746 return;
747 }
9b2de807 748
96c271c1 749 fh1VertexZ->Fill(primVtx->GetZ());
9b2de807 750
96c271c1 751 if(TMath::Abs(primVtx->GetZ())>fMaxVertexZ){
752 if (fDebug > 1) Printf("%s:%d primary vertex z = %f: event REJECTED...",(char*)__FILE__,__LINE__,primVtx->GetZ());
753 fh1EvtSelection->Fill(4.);
754 PostData(1, fCommonHistList);
755 return;
756 }
9b2de807 757
96c271c1 758 TString primVtxName(primVtx->GetName());
9b2de807 759
96c271c1 760 if(primVtxName.CompareTo("TPCVertex",TString::kIgnoreCase) == 1){
761 if (fDebug > 1) Printf("%s:%d primary vertex selection: TPC vertex, event REJECTED...",(char*)__FILE__,__LINE__);
762 fh1EvtSelection->Fill(5.);
763 PostData(1, fCommonHistList);
9b2de807 764 return;
765 }
766
96c271c1 767 if (fDebug > 1) Printf("%s:%d event ACCEPTED ...",(char*)__FILE__,__LINE__);
768 fh1EvtSelection->Fill(0.);
769 fh1EvtCent->Fill(centPercent);
770
9b2de807 771
96c271c1 772 //___ get MC information __________________________________________________________________
9b2de807 773
96c271c1 774 Double_t ptHard = 0.;
775 Double_t nTrials = 1; // trials for MC trigger weight for real data
9b2de807 776
96c271c1 777 if(fMCEvent){
778 AliGenEventHeader* genHeader = fMCEvent->GenEventHeader();
779 AliGenPythiaEventHeader* pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
780 AliGenHijingEventHeader* hijingGenHeader = 0x0;
9b2de807 781
96c271c1 782 if(pythiaGenHeader){
783 if(fDebug>3) Printf("%s:%d pythiaGenHeader found", (char*)__FILE__,__LINE__);
784 nTrials = pythiaGenHeader->Trials();
785 ptHard = pythiaGenHeader->GetPtHard();
9b2de807 786
96c271c1 787 fh1PtHard->Fill(ptHard);
788 fh1PtHardTrials->Fill(ptHard,nTrials);
9280dfa6 789
9280dfa6 790
96c271c1 791 } else { // no pythia, hijing?
9280dfa6 792
96c271c1 793 if(fDebug>3) Printf("%s:%d no pythiaGenHeader found", (char*)__FILE__,__LINE__);
9280dfa6 794
96c271c1 795 hijingGenHeader = dynamic_cast<AliGenHijingEventHeader*>(genHeader);
796 if(!hijingGenHeader){
797 Printf("%s:%d no pythiaGenHeader or hjingGenHeader found", (char*)__FILE__,__LINE__);
798 } else {
799 if(fDebug>3) Printf("%s:%d hijingGenHeader found", (char*)__FILE__,__LINE__);
800 }
801 }
9280dfa6 802
96c271c1 803 fh1Trials->Fill("#sum{ntrials}",fAvgTrials);
9b2de807 804 }
805
9b2de807 806
96c271c1 807 //____ fetch jets ______________________________________________________________
9b2de807 808
96c271c1 809 Int_t nJCuts = GetListOfJets(fJetsRecCuts, kJetsRecAcceptance);
810 Int_t nRecJetsCuts = 0;
811 if(nJCuts>=0) nRecJetsCuts = fJetsRecCuts->GetEntries();
812 if(fDebug>2)Printf("%s:%d Selected Rec jets after cuts: %d %d",(char*)__FILE__,__LINE__,nJCuts,nRecJetsCuts);
813 if(nRecJetsCuts != nJCuts) Printf("%s:%d Mismatch selected Rec jets after cuts: %d %d",(char*)__FILE__,__LINE__,nJCuts,nRecJetsCuts);
814 fh1nRecJetsCuts->Fill(nRecJetsCuts);
9280dfa6 815
816
96c271c1 817 //____ fetch particles __________________________________________________________
9b2de807 818
96c271c1 819 Int_t nTCuts = GetListOfTracks(fTracksRecCuts, kTrackAODCuts);
820 if(fDebug>2)Printf("%s:%d Selected Rec tracks after cuts: %d %d",(char*)__FILE__,__LINE__,nTCuts,fTracksRecCuts->GetEntries());
821 if(fTracksRecCuts->GetEntries() != nTCuts)
822 Printf("%s:%d Mismatch selected Rec tracks after cuts: %d %d",(char*)__FILE__,__LINE__,nTCuts,fTracksRecCuts->GetEntries());
823 fh1EvtMult->Fill(fTracksRecCuts->GetEntries());
9b2de807 824
96c271c1 825 Int_t nK0s = GetListOfK0s(fListK0s,fK0Type);
826 if(fDebug>2)Printf("%s:%d Selected Rec tracks after cuts: %d %d",(char*)__FILE__,__LINE__,nK0s,fListK0s->GetEntries());
827 if(nK0s != fListK0s->GetEntries()) Printf("%s:%d Mismatch selected K0s: %d %d",(char*)__FILE__,__LINE__,nK0s,fListK0s->GetEntries());
828 fh1K0Mult->Fill(fListK0s->GetEntries());
9b2de807 829
96c271c1 830 // ___ V0 QA + K0 pt spectra all events _______________________________________________
9b2de807 831
96c271c1 832 if(fListK0s->GetEntries()>0){
9b2de807 833
96c271c1 834 for(Int_t it=0; it<fListK0s->GetSize(); ++it){ // loop all K0s
835
836 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListK0s->At(it));
f98082ec 837 if(!v0) continue;
838
96c271c1 839 Float_t trackPt = v0->Pt();
840 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
841 Double_t invM = v0->MassK0Short();
842 Double_t jetPt = fFFIMJetPtMin; // assign pro forma jet energy
9b2de807 843
96c271c1 844 fV0QAK0->FillTrackQA(v0->Eta(), TVector2::Phi_0_2pi(v0->Phi()), v0->Pt());
845 fFFHistosIMK0AllEvt->FillFF(trackPt, invM, jetPt, incrementJetPt);
846 }
9b2de807 847 }
848
9b2de807 849
96c271c1 850 //____ fill FF histos __________________________________________________________
9b2de807 851
96c271c1 852 for(Int_t ij=0; ij<nRecJetsCuts; ++ij){
9b2de807 853
c500e0a0 854 AliAODJet* jet = (AliAODJet*) (fJetsRecCuts->At(ij));
6daac008 855 Double_t jetPt = jet->Pt();
856
96c271c1 857 if(ij==0){ // leading jet
858
859 TList* jettracklist = new TList();
860 Double_t sumPt = 0.;
6daac008 861 Bool_t isBadJet = kFALSE;
96c271c1 862
863 if(GetFFRadius()<=0){
6daac008 864 GetJetTracksTrackrefs(jettracklist, jet, GetFFMaxTrackPt(), isBadJet);
96c271c1 865 } else {
6daac008 866 GetJetTracksPointing(fTracksRecCuts, jettracklist, jet, GetFFRadius(), sumPt, GetFFMaxTrackPt(), isBadJet);
96c271c1 867 }
6daac008 868
869 if(GetFFMinNTracks()>0 && jettracklist->GetSize() <= GetFFMinNTracks()) isBadJet = kTRUE;
870 if(isBadJet) continue;
871
96c271c1 872
873 for(Int_t it=0; it<jettracklist->GetSize(); ++it){
9b2de807 874
c500e0a0 875 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(jettracklist->At(it));
876 if(!trackVP)continue;
877
96c271c1 878 Float_t trackPt = trackVP->Pt();
879
880 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
881
882 fFFHistosRecCuts->FillFF(trackPt, jetPt, incrementJetPt);
883 if(nK0s>0) fFFHistosRecCutsK0Evt->FillFF(trackPt, jetPt, incrementJetPt);
884 }
885
886 delete jettracklist;
9b2de807 887
6daac008 888 // ---- K0s ----
96c271c1 889
6daac008 890 // fQAJetHistosRecCutsLeading->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt() );
9b2de807 891
96c271c1 892 for(Int_t it=0; it<fListK0s->GetSize(); ++it){ // loop all K0s
9b2de807 893
96c271c1 894 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListK0s->At(it));
f98082ec 895 if(!v0) continue;
896
96c271c1 897 Double_t v0Mom[3];
898 v0->PxPyPz(v0Mom);
899 TVector3 v0MomVect(v0Mom);
9b2de807 900
96c271c1 901 Double_t dPhiJetK0 = (jet->MomentumVector()->Vect()).DeltaPhi(v0MomVect);
902
903 Float_t trackPt = v0->Pt();
904 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
905 Double_t invM = v0->MassK0Short();
906
907 fFFHistosIMK0Jet->FillFF(trackPt, invM, jetPt, incrementJetPt);
908 fFFHistosPhiCorrIMK0->FillPhiCorr(trackPt,TVector2::Phi_0_2pi(dPhiJetK0),invM);
9b2de807 909
96c271c1 910 if(dPhiJetK0<fh1dPhiJetK0->GetXaxis()->GetXmin()) dPhiJetK0 += 2*TMath::Pi();
911 fh1dPhiJetK0->Fill(dPhiJetK0);
9b2de807 912
96c271c1 913 }
6daac008 914
96c271c1 915 if(fListK0s->GetSize() == 0){ // no K0: increment jet pt spectrum
9b2de807 916
96c271c1 917 Bool_t incrementJetPt = kTRUE;
918 fFFHistosIMK0Jet->FillFF(-1, -1, jetPt, incrementJetPt);
919 }
920
9b2de807 921
96c271c1 922 TList* jetConeK0list = new TList();
6daac008 923 Double_t sumPtK0 = 0.;
924 Bool_t isBadJetK0 = kFALSE; // dummy, do not use
925
926 GetJetTracksPointing(fListK0s, jetConeK0list, jet, GetFFRadius(), sumPtK0, GetFFMaxTrackPt(), isBadJetK0);
927
928
96c271c1 929 if(fDebug>2)Printf("%s:%d nK0s total: %d, in jet cone: %d,FFRadius %f ",(char*)__FILE__,__LINE__,nK0s,jetConeK0list->GetEntries(),GetFFRadius());
930
931 for(Int_t it=0; it<jetConeK0list->GetSize(); ++it){ // loop K0s in jet cone
c500e0a0 932
96c271c1 933 AliAODv0* v0 = dynamic_cast<AliAODv0*>(jetConeK0list->At(it));
c500e0a0 934 if(!v0) continue;
9b2de807 935
96c271c1 936 Double_t invM = v0->MassK0Short();
937 Float_t trackPt = v0->Pt();
938 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
939
940 //std::cout<<" trackPt "<<trackPt<<" invM "<<invM<<std::endl;
941 fFFHistosIMK0Cone->FillFF(trackPt, invM, jetPt, incrementJetPt);
942 }
943 if(jetConeK0list->GetSize() == 0){ // no K0: increment jet pt spectrum
9b2de807 944
96c271c1 945 Bool_t incrementJetPt = kTRUE;
946 fFFHistosIMK0Cone->FillFF(-1, -1, jetPt, incrementJetPt);
947 }
9b2de807 948
96c271c1 949 delete jetConeK0list;
950 }
951 }
952
953 fTracksRecCuts->Clear();
954 fJetsRecCuts->Clear();
955 fListK0s->Clear();
9b2de807 956
96c271c1 957 //Post output data.
958 PostData(1, fCommonHistList);
9b2de807 959}
960
96c271c1 961// ____________________________________________________________________________________________
962void AliAnalysisTaskJetChem::SetProperties(TH3F* h,const char* x, const char* y, const char* z)
963{
964 //Set properties of histos (x,y and z title)
965
966 h->SetXTitle(x);
967 h->SetYTitle(y);
968 h->SetZTitle(z);
969 h->GetXaxis()->SetTitleColor(1);
970 h->GetYaxis()->SetTitleColor(1);
971 h->GetZaxis()->SetTitleColor(1);
9b2de807 972}
973
96c271c1 974// ____________________________________________________________________________________________________________________________________________
9b2de807 975Bool_t AliAnalysisTaskJetChem::IsAccepteddEdx(const Double_t mom,const Double_t signal, AliPID::EParticleType n, const Double_t cutnSig) const{
976
977 // apply TPC dE/dx cut similar as in AliTPCpidESD
978 // note: AliTPCidESD uses inner track param for momentum - not avaiable on AOD,
979 // so we use global track momentum
980 // should use separate parametrisation for MC and data, but probably ALEPH param & 7% resolution used here anyway not the last word
981
982
983 const Double_t kBBMIP(50.);
984 const Double_t kBBRes(0.07);
985 //const Double_t kBBRange(5.);
986 const Double_t kBBp1(0.76176e-1);
987 const Double_t kBBp2(10.632);
988 const Double_t kBBp3(0.13279e-4);
989 const Double_t kBBp4(1.8631);
990 const Double_t kBBp5(1.9479);
991
992 Double_t mass=AliPID::ParticleMass(n);
993 Double_t betaGamma = mom/mass;
994
995 const Float_t kmeanCorrection =0.1;
996 Double_t bb = AliExternalTrackParam::BetheBlochAleph(betaGamma,kBBp1,kBBp2,kBBp3,kBBp4,kBBp5);
997 Double_t meanCorrection =(1+(bb-1)*kmeanCorrection);
998 Double_t bethe = bb * meanCorrection; // expected
999 Double_t sigma = bethe * kBBRes;
1000
1001
1002 Double_t dedx = signal/kBBMIP; // measured
1003
1004 Double_t nSig = (TMath::Abs(dedx - bethe))/sigma;
1005
1006 if(nSig > cutnSig) return kFALSE;
1007
1008 return kTRUE;
1009}
1010
96c271c1 1011//___________________________________________________________________
1012Bool_t AliAnalysisTaskJetChem::IsK0InvMass(const Double_t mass) const
1013{
1014 // K0 mass ? Use FF histo limits
1015
1016 if(fFFIMInvMMin <= mass && mass <= fFFIMInvMMax) return kTRUE;
9b2de807 1017
96c271c1 1018 return kFALSE;
1019}
9280dfa6 1020
96c271c1 1021//_____________________________________________________________________________________
1022Int_t AliAnalysisTaskJetChem::GetListOfK0s(TList *list, const Int_t type)
1023{
1024 // fill list of V0s selected according to type
9b2de807 1025
96c271c1 1026 if(!list){
1027 if(fDebug>1) Printf("%s:%d no input list", (char*)__FILE__,__LINE__);
1028 return -1;
1029 }
1030
1031 if(type==kTrackUndef) return 0;
9b2de807 1032
1033 for(int i=0; i<fAOD->GetNumberOfV0s(); i++){ // loop over V0s
1034
1035 AliAODv0* v0 = fAOD->GetV0(i);
1036
1037 Bool_t isOnFly = v0->GetOnFlyStatus();
9b2de807 1038
96c271c1 1039 if(!isOnFly && (type == kOnFly || type == kOnFlyPID || type == kOnFlydEdx || type == kOnFlyPrim)) continue;
1040 if( isOnFly && (type == kOffl || type == kOfflPID || type == kOffldEdx || type == kOfflPrim)) continue;
656457c9 1041
96c271c1 1042 Double_t massK0 = v0->MassK0Short();
656457c9 1043
96c271c1 1044 if(!(IsK0InvMass(massK0))) continue; // moved invMass cut for HI - otherwise too slow
656457c9 1045
96c271c1 1046 if(type == kOnFlyPID || type == kOfflPID){
1047
1048 AliAODTrack *trackPos = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(0)); // slow
1049 AliAODTrack *trackNeg = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(1)); // slow
1050
1051 // AOD pid - cuts strongly into signal
9280dfa6 1052
96c271c1 1053 AliAODTrack::AODTrkPID_t mpPIDNeg = trackNeg->GetMostProbablePID();
1054 AliAODTrack::AODTrkPID_t mpPIDPos = trackPos->GetMostProbablePID();
1055
1056 if(!( (mpPIDNeg == AliAODTrack::kPion) && (mpPIDPos == AliAODTrack::kPion) ) ) continue;
1057 }
1058
1059 if(type == kOnFlydEdx || type == kOffldEdx){
9b2de807 1060
656457c9 1061 AliAODTrack *trackPos = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(0)); // slow
1062 AliAODTrack *trackNeg = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(1)); // slow
9b2de807 1063
9b2de807 1064 AliAODPid* aodPidPos = trackPos->GetDetPid();
1065 AliAODPid* aodPidNeg = trackNeg->GetDetPid();
1066
1067 Double_t dEdxPos = aodPidPos->GetTPCsignal();
1068 Double_t dEdxNeg = aodPidNeg->GetTPCsignal();
1069
1070 Double_t momPos = trackPos->P();
1071 Double_t momNeg = trackNeg->P();
9b2de807 1072
96c271c1 1073 Int_t cutnSigdEdx = 2;
1074 if(! (IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,cutnSigdEdx)) && (IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,cutnSigdEdx)) ) continue;
9b2de807 1075
96c271c1 1076 }
9b2de807 1077
96c271c1 1078 if(type == kOnFlyPrim || type == kOfflPrim){
9b2de807 1079
96c271c1 1080 AliAODTrack *trackPos = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(0)); // slow
1081 AliAODTrack *trackNeg = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(1)); // slow
9b2de807 1082
96c271c1 1083 //std::cout<<" filer map trackPos "<<trackPos->GetFilterMap()<<" trackNeg "<<trackNeg->GetFilterMap()<<std::endl;
9b2de807 1084
96c271c1 1085 if((fFilterMaskK0>0) && !(trackPos->TestFilterBit(fFilterMaskK0))) continue;
1086 if((fFilterMaskK0>0) && !(trackNeg->TestFilterBit(fFilterMaskK0))) continue;
9b2de807 1087 }
9b2de807 1088
96c271c1 1089 list->Add(v0);
9b2de807 1090 }
9280dfa6 1091
96c271c1 1092 Int_t nK0s = list->GetSize();
9280dfa6 1093
96c271c1 1094 return nK0s;
9280dfa6 1095}
1096
9280dfa6 1097