]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGJE/AliAnalysisTaskFragmentationFunction.cxx
nicer printout
[u/mrichter/AliRoot.git] / PWGJE / AliAnalysisTaskFragmentationFunction.cxx
CommitLineData
ffa175ff 1// *************************************************************************
2// * *
3// * Task for Fragmentation Function Analysis in PWG4 Jet Task Force Train *
4// * *
5// *************************************************************************
6
7
8/**************************************************************************
9 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
10 * *
11 * Author: The ALICE Off-line Project. *
12 * Contributors are mentioned in the code where appropriate. *
13 * *
14 * Permission to use, copy, modify and distribute this software and its *
15 * documentation strictly for non-commercial purposes is hereby granted *
16 * without fee, provided that the above copyright notice appears in all *
17 * copies and that both the copyright notice and this permission notice *
18 * appear in the supporting documentation. The authors make no claims *
19 * about the suitability of this software for any purpose. It is *
20 * provided "as is" without express or implied warranty. *
21 **************************************************************************/
22
23/* $Id: */
24
25#include "TList.h"
26#include "TH1F.h"
27#include "TH2F.h"
28#include "TH3F.h"
29#include "TString.h"
30#include "THnSparse.h"
31#include "TProfile.h"
32#include "TFile.h"
33#include "TKey.h"
34#include "TRandom3.h"
35#include "TAxis.h"
36
37#include "AliAODInputHandler.h"
38#include "AliAODHandler.h"
39#include "AliESDEvent.h"
40#include "AliAODMCParticle.h"
41#include "AliAODJet.h"
42#include "AliAODJetEventBackground.h"
43#include "AliGenPythiaEventHeader.h"
44#include "AliGenHijingEventHeader.h"
45#include "AliInputEventHandler.h"
46
47#include "AliAnalysisHelperJetTasks.h"
48#include "AliAnalysisManager.h"
49#include "AliAnalysisTaskSE.h"
50#include "AliVParticle.h"
51#include "AliVEvent.h"
52
53#include "AliAnalysisTaskFragmentationFunction.h"
c2aad3ae 54using std::cout;
55using std::endl;
56using std::cerr;
ffa175ff 57
58ClassImp(AliAnalysisTaskFragmentationFunction)
59
60//____________________________________________________________________________
61AliAnalysisTaskFragmentationFunction::AliAnalysisTaskFragmentationFunction()
62 : AliAnalysisTaskSE()
63 ,fESD(0)
64 ,fAOD(0)
65 ,fAODJets(0)
66 ,fAODExtension(0)
67 ,fNonStdFile("")
68 ,fBranchRecJets("jets")
69 ,fBranchRecBckgClusters("")
70 ,fBranchGenJets("")
71 ,fBranchEmbeddedJets("")
72 ,fTrackTypeGen(0)
73 ,fJetTypeGen(0)
74 ,fJetTypeRecEff(0)
75 ,fUseAODInputJets(kTRUE)
76 ,fFilterMask(0)
77 ,fUsePhysicsSelection(kTRUE)
78 ,fEvtSelectionMask(0)
79 ,fEventClass(0)
80 ,fMaxVertexZ(10)
81 ,fTrackPtCut(0)
82 ,fTrackEtaMin(0)
83 ,fTrackEtaMax(0)
84 ,fTrackPhiMin(0)
85 ,fTrackPhiMax(0)
86 ,fUseExtraTracks(0)
87 ,fUseExtraTracksBgr(0)
88 ,fCutFractionPtEmbedded(0)
89 ,fUseEmbeddedJetAxis(0)
90 ,fUseEmbeddedJetPt(0)
91 ,fJetPtCut(0)
92 ,fJetEtaMin(0)
93 ,fJetEtaMax(0)
94 ,fJetPhiMin(0)
95 ,fJetPhiMax(0)
96 ,fFFRadius(0)
97 ,fFFMinLTrackPt(-1)
98 ,fFFMaxTrackPt(-1)
99 ,fFFMinnTracks(0)
100 ,fFFBckgRadius(0)
101 ,fBckgMode(0)
102 ,fQAMode(0)
103 ,fFFMode(0)
104 ,fEffMode(0)
105 ,fJSMode(0)
106 ,fAvgTrials(0)
107 ,fTracksRecCuts(0)
108 ,fTracksGen(0)
109 ,fTracksAODMCCharged(0)
110 ,fTracksAODMCChargedSecNS(0)
111 ,fTracksAODMCChargedSecS(0)
112 ,fTracksRecQualityCuts(0)
113 ,fJetsRec(0)
114 ,fJetsRecCuts(0)
115 ,fJetsGen(0)
116 ,fJetsRecEff(0)
117 ,fJetsEmbedded(0)
118 ,fBckgJetsRec(0)
119 ,fBckgJetsRecCuts(0)
120 ,fBckgJetsGen(0)
121 ,fQATrackHistosRecCuts(0)
122 ,fQATrackHistosGen(0)
123 ,fQAJetHistosRec(0)
124 ,fQAJetHistosRecCuts(0)
125 ,fQAJetHistosRecCutsLeading(0)
126 ,fQAJetHistosGen(0)
127 ,fQAJetHistosGenLeading(0)
128 ,fQAJetHistosRecEffLeading(0)
129 ,fFFHistosRecCuts(0)
37f35567 130 ,fFFHistosRecCutsInc(0)
131 ,fFFHistosRecLeadingTrack(0)
ffa175ff 132 ,fFFHistosGen(0)
37f35567 133 ,fFFHistosGenInc(0)
134 ,fFFHistosGenLeadingTrack(0)
ffa175ff 135 ,fQATrackHighPtThreshold(0)
136 ,fFFNBinsJetPt(0)
137 ,fFFJetPtMin(0)
138 ,fFFJetPtMax(0)
139 ,fFFNBinsPt(0)
140 ,fFFPtMin(0)
141 ,fFFPtMax(0)
142 ,fFFNBinsXi(0)
143 ,fFFXiMin(0)
144 ,fFFXiMax(0)
145 ,fFFNBinsZ(0)
146 ,fFFZMin(0)
147 ,fFFZMax(0)
148 ,fQAJetNBinsPt(0)
149 ,fQAJetPtMin(0)
150 ,fQAJetPtMax(0)
151 ,fQAJetNBinsEta(0)
152 ,fQAJetEtaMin(0)
153 ,fQAJetEtaMax(0)
154 ,fQAJetNBinsPhi(0)
155 ,fQAJetPhiMin(0)
156 ,fQAJetPhiMax(0)
157 ,fQATrackNBinsPt(0)
158 ,fQATrackPtMin(0)
159 ,fQATrackPtMax(0)
160 ,fQATrackNBinsEta(0)
161 ,fQATrackEtaMin(0)
162 ,fQATrackEtaMax(0)
163 ,fQATrackNBinsPhi(0)
164 ,fQATrackPhiMin(0)
165 ,fQATrackPhiMax(0)
166 ,fCommonHistList(0)
167 ,fh1EvtSelection(0)
168 ,fh1VertexNContributors(0)
169 ,fh1VertexZ(0)
170 ,fh1EvtMult(0)
171 ,fh1EvtCent(0)
172 ,fh1Xsec(0)
173 ,fh1Trials(0)
174 ,fh1PtHard(0)
175 ,fh1PtHardTrials(0)
176 ,fh1nRecJetsCuts(0)
177 ,fh1nGenJets(0)
178 ,fh1nRecEffJets(0)
179 ,fh1nEmbeddedJets(0)
180 ,fh1nRecBckgJetsCuts(0)
181 ,fh1nGenBckgJets(0)
182 ,fh2PtRecVsGenPrim(0)
183 ,fh2PtRecVsGenSec(0)
184 ,fQATrackHistosRecEffGen(0)
185 ,fQATrackHistosRecEffRec(0)
186 ,fQATrackHistosSecRecNS(0)
187 ,fQATrackHistosSecRecS(0)
188 ,fQATrackHistosSecRecSsc(0)
189 ,fFFHistosRecEffRec(0)
190 ,fFFHistosSecRecNS(0)
191 ,fFFHistosSecRecS(0)
192 ,fFFHistosSecRecSsc(0)
193 // Background
194 ,fh1BckgMult0(0)
195 ,fh1BckgMult1(0)
196 ,fh1BckgMult2(0)
197 ,fh1BckgMult3(0)
198 ,fh1BckgMult4(0)
199 ,fh1FractionPtEmbedded(0)
200 ,fh1IndexEmbedded(0)
201 ,fh2DeltaPtVsJetPtEmbedded(0)
202 ,fh2DeltaPtVsRecJetPtEmbedded(0)
203 ,fh1DeltaREmbedded(0)
204 ,fQABckgHisto0RecCuts(0)
205 ,fQABckgHisto0Gen(0)
206 ,fQABckgHisto1RecCuts(0)
207 ,fQABckgHisto1Gen(0)
208 ,fQABckgHisto2RecCuts(0)
209 ,fQABckgHisto2Gen(0)
210 ,fQABckgHisto3RecCuts(0)
211 ,fQABckgHisto3Gen(0)
212 ,fQABckgHisto4RecCuts(0)
213 ,fQABckgHisto4Gen(0)
214 ,fFFBckgHisto0RecCuts(0)
215 ,fFFBckgHisto0Gen(0)
216 ,fFFBckgHisto1RecCuts(0)
217 ,fFFBckgHisto1Gen(0)
218 ,fFFBckgHisto2RecCuts(0)
219 ,fFFBckgHisto2Gen(0)
220 ,fFFBckgHisto3RecCuts(0)
221 ,fFFBckgHisto3Gen(0)
222 ,fFFBckgHisto4RecCuts(0)
223 ,fFFBckgHisto4Gen(0)
224 ,fFFBckgHisto0RecEffRec(0)
225 ,fFFBckgHisto0SecRecNS(0)
226 ,fFFBckgHisto0SecRecS(0)
227 ,fFFBckgHisto0SecRecSsc(0)
228 // jet shape
229 ,fProNtracksLeadingJet(0)
230 ,fProDelR80pcPt(0)
231 ,fProNtracksLeadingJetGen(0)
232 ,fProDelR80pcPtGen(0)
233 ,fProNtracksLeadingJetBgrPerp2(0)
234 ,fProNtracksLeadingJetRecPrim(0)
235 ,fProDelR80pcPtRecPrim(0)
236 ,fProNtracksLeadingJetRecSecNS(0)
237 ,fProNtracksLeadingJetRecSecS(0)
238 ,fProNtracksLeadingJetRecSecSsc(0)
239
240 ,fRandom(0)
241{
242 // default constructor
243 fBckgType[0] = 0;
244 fBckgType[1] = 0;
245 fBckgType[2] = 0;
246 fBckgType[3] = 0;
247 fBckgType[4] = 0;
248
249 for(Int_t ii=0; ii<5; ii++){
250 fProDelRPtSum[ii] = 0;
251 fProDelRPtSumGen[ii] = 0;
252 fProDelRPtSumBgrPerp2[ii] = 0;
253 fProDelRPtSumRecPrim[ii] = 0;
254 fProDelRPtSumRecSecNS[ii] = 0;
255 fProDelRPtSumRecSecS[ii] = 0;
256 fProDelRPtSumRecSecSsc[ii] = 0;
257 }
258}
259
260//_______________________________________________________________________________________________
261AliAnalysisTaskFragmentationFunction::AliAnalysisTaskFragmentationFunction(const char *name)
262 : AliAnalysisTaskSE(name)
263 ,fESD(0)
264 ,fAOD(0)
265 ,fAODJets(0)
266 ,fAODExtension(0)
267 ,fNonStdFile("")
268 ,fBranchRecJets("jets")
269 ,fBranchRecBckgClusters("")
270 ,fBranchGenJets("")
271 ,fBranchEmbeddedJets("")
272 ,fTrackTypeGen(0)
273 ,fJetTypeGen(0)
274 ,fJetTypeRecEff(0)
275 ,fUseAODInputJets(kTRUE)
276 ,fFilterMask(0)
277 ,fUsePhysicsSelection(kTRUE)
278 ,fEvtSelectionMask(0)
279 ,fEventClass(0)
280 ,fMaxVertexZ(10)
281 ,fTrackPtCut(0)
282 ,fTrackEtaMin(0)
283 ,fTrackEtaMax(0)
284 ,fTrackPhiMin(0)
285 ,fTrackPhiMax(0)
286 ,fUseExtraTracks(0)
287 ,fUseExtraTracksBgr(0)
288 ,fCutFractionPtEmbedded(0)
289 ,fUseEmbeddedJetAxis(0)
290 ,fUseEmbeddedJetPt(0)
291 ,fJetPtCut(0)
292 ,fJetEtaMin(0)
293 ,fJetEtaMax(0)
294 ,fJetPhiMin(0)
295 ,fJetPhiMax(0)
296 ,fFFRadius(0)
297 ,fFFMinLTrackPt(-1)
298 ,fFFMaxTrackPt(-1)
299 ,fFFMinnTracks(0)
300 ,fFFBckgRadius(0)
301 ,fBckgMode(0)
302 ,fQAMode(0)
303 ,fFFMode(0)
304 ,fEffMode(0)
305 ,fJSMode(0)
306 ,fAvgTrials(0)
307 ,fTracksRecCuts(0)
308 ,fTracksGen(0)
309 ,fTracksAODMCCharged(0)
310 ,fTracksAODMCChargedSecNS(0)
311 ,fTracksAODMCChargedSecS(0)
312 ,fTracksRecQualityCuts(0)
313 ,fJetsRec(0)
314 ,fJetsRecCuts(0)
315 ,fJetsGen(0)
316 ,fJetsRecEff(0)
317 ,fJetsEmbedded(0)
318 ,fBckgJetsRec(0)
319 ,fBckgJetsRecCuts(0)
320 ,fBckgJetsGen(0)
321 ,fQATrackHistosRecCuts(0)
322 ,fQATrackHistosGen(0)
323 ,fQAJetHistosRec(0)
324 ,fQAJetHistosRecCuts(0)
325 ,fQAJetHistosRecCutsLeading(0)
326 ,fQAJetHistosGen(0)
327 ,fQAJetHistosGenLeading(0)
328 ,fQAJetHistosRecEffLeading(0)
329 ,fFFHistosRecCuts(0)
37f35567 330 ,fFFHistosRecCutsInc(0)
331 ,fFFHistosRecLeadingTrack(0)
ffa175ff 332 ,fFFHistosGen(0)
37f35567 333 ,fFFHistosGenInc(0)
334 ,fFFHistosGenLeadingTrack(0)
ffa175ff 335 ,fQATrackHighPtThreshold(0)
336 ,fFFNBinsJetPt(0)
337 ,fFFJetPtMin(0)
338 ,fFFJetPtMax(0)
339 ,fFFNBinsPt(0)
340 ,fFFPtMin(0)
341 ,fFFPtMax(0)
342 ,fFFNBinsXi(0)
343 ,fFFXiMin(0)
344 ,fFFXiMax(0)
345 ,fFFNBinsZ(0)
346 ,fFFZMin(0)
347 ,fFFZMax(0)
348 ,fQAJetNBinsPt(0)
349 ,fQAJetPtMin(0)
350 ,fQAJetPtMax(0)
351 ,fQAJetNBinsEta(0)
352 ,fQAJetEtaMin(0)
353 ,fQAJetEtaMax(0)
354 ,fQAJetNBinsPhi(0)
355 ,fQAJetPhiMin(0)
356 ,fQAJetPhiMax(0)
357 ,fQATrackNBinsPt(0)
358 ,fQATrackPtMin(0)
359 ,fQATrackPtMax(0)
360 ,fQATrackNBinsEta(0)
361 ,fQATrackEtaMin(0)
362 ,fQATrackEtaMax(0)
363 ,fQATrackNBinsPhi(0)
364 ,fQATrackPhiMin(0)
365 ,fQATrackPhiMax(0)
366 ,fCommonHistList(0)
367 ,fh1EvtSelection(0)
368 ,fh1VertexNContributors(0)
369 ,fh1VertexZ(0)
370 ,fh1EvtMult(0)
371 ,fh1EvtCent(0)
372 ,fh1Xsec(0)
373 ,fh1Trials(0)
374 ,fh1PtHard(0)
375 ,fh1PtHardTrials(0)
376 ,fh1nRecJetsCuts(0)
377 ,fh1nGenJets(0)
378 ,fh1nRecEffJets(0)
379 ,fh1nEmbeddedJets(0)
380 ,fh1nRecBckgJetsCuts(0)
381 ,fh1nGenBckgJets(0)
382 ,fh2PtRecVsGenPrim(0)
383 ,fh2PtRecVsGenSec(0)
384 ,fQATrackHistosRecEffGen(0)
385 ,fQATrackHistosRecEffRec(0)
386 ,fQATrackHistosSecRecNS(0)
387 ,fQATrackHistosSecRecS(0)
388 ,fQATrackHistosSecRecSsc(0)
389 ,fFFHistosRecEffRec(0)
390 ,fFFHistosSecRecNS(0)
391 ,fFFHistosSecRecS(0)
392 ,fFFHistosSecRecSsc(0)
393 // Background
394 ,fh1BckgMult0(0)
395 ,fh1BckgMult1(0)
396 ,fh1BckgMult2(0)
397 ,fh1BckgMult3(0)
398 ,fh1BckgMult4(0)
399 ,fh1FractionPtEmbedded(0)
400 ,fh1IndexEmbedded(0)
401 ,fh2DeltaPtVsJetPtEmbedded(0)
402 ,fh2DeltaPtVsRecJetPtEmbedded(0)
403 ,fh1DeltaREmbedded(0)
404 ,fQABckgHisto0RecCuts(0)
405 ,fQABckgHisto0Gen(0)
406 ,fQABckgHisto1RecCuts(0)
407 ,fQABckgHisto1Gen(0)
408 ,fQABckgHisto2RecCuts(0)
409 ,fQABckgHisto2Gen(0)
410 ,fQABckgHisto3RecCuts(0)
411 ,fQABckgHisto3Gen(0)
412 ,fQABckgHisto4RecCuts(0)
413 ,fQABckgHisto4Gen(0)
414 ,fFFBckgHisto0RecCuts(0)
415 ,fFFBckgHisto0Gen(0)
416 ,fFFBckgHisto1RecCuts(0)
417 ,fFFBckgHisto1Gen(0)
418 ,fFFBckgHisto2RecCuts(0)
419 ,fFFBckgHisto2Gen(0)
420 ,fFFBckgHisto3RecCuts(0)
421 ,fFFBckgHisto3Gen(0)
422 ,fFFBckgHisto4RecCuts(0)
423 ,fFFBckgHisto4Gen(0)
424 ,fFFBckgHisto0RecEffRec(0)
425 ,fFFBckgHisto0SecRecNS(0)
426 ,fFFBckgHisto0SecRecS(0)
427 ,fFFBckgHisto0SecRecSsc(0)
428 // jet shape
429 ,fProNtracksLeadingJet(0)
430 ,fProDelR80pcPt(0)
431 ,fProNtracksLeadingJetGen(0)
432 ,fProDelR80pcPtGen(0)
433 ,fProNtracksLeadingJetBgrPerp2(0)
434 ,fProNtracksLeadingJetRecPrim(0)
435 ,fProDelR80pcPtRecPrim(0)
436 ,fProNtracksLeadingJetRecSecNS(0)
437 ,fProNtracksLeadingJetRecSecS(0)
438 ,fProNtracksLeadingJetRecSecSsc(0)
439 ,fRandom(0)
440{
441 // constructor
442 fBckgType[0] = 0;
443 fBckgType[1] = 0;
444 fBckgType[2] = 0;
445 fBckgType[3] = 0;
446 fBckgType[4] = 0;
447
448 for(Int_t ii=0; ii<5; ii++){
449 fProDelRPtSum[ii] = 0;
450 fProDelRPtSumGen[ii] = 0;
451 fProDelRPtSumBgrPerp2[ii] = 0;
452 fProDelRPtSumRecPrim[ii] = 0;
453 fProDelRPtSumRecSecNS[ii] = 0;
454 fProDelRPtSumRecSecS[ii] = 0;
455 fProDelRPtSumRecSecSsc[ii] = 0;
456 }
457
458 DefineOutput(1,TList::Class());
459}
460
461//__________________________________________________________________________________________________________________________
462AliAnalysisTaskFragmentationFunction::AliAnalysisTaskFragmentationFunction(const AliAnalysisTaskFragmentationFunction &copy)
463 : AliAnalysisTaskSE()
464 ,fESD(copy.fESD)
465 ,fAOD(copy.fAOD)
466 ,fAODJets(copy.fAODJets)
467 ,fAODExtension(copy.fAODExtension)
468 ,fNonStdFile(copy.fNonStdFile)
469 ,fBranchRecJets(copy.fBranchRecJets)
470 ,fBranchRecBckgClusters(copy.fBranchRecBckgClusters)
471 ,fBranchGenJets(copy.fBranchGenJets)
472 ,fBranchEmbeddedJets(copy.fBranchEmbeddedJets)
473 ,fTrackTypeGen(copy.fTrackTypeGen)
474 ,fJetTypeGen(copy.fJetTypeGen)
475 ,fJetTypeRecEff(copy.fJetTypeRecEff)
476 ,fUseAODInputJets(copy.fUseAODInputJets)
477 ,fFilterMask(copy.fFilterMask)
478 ,fUsePhysicsSelection(copy.fUsePhysicsSelection)
479 ,fEvtSelectionMask(copy.fEvtSelectionMask)
480 ,fEventClass(copy.fEventClass)
481 ,fMaxVertexZ(copy.fMaxVertexZ)
482 ,fTrackPtCut(copy.fTrackPtCut)
483 ,fTrackEtaMin(copy.fTrackEtaMin)
484 ,fTrackEtaMax(copy.fTrackEtaMax)
485 ,fTrackPhiMin(copy.fTrackPhiMin)
486 ,fTrackPhiMax(copy.fTrackPhiMax)
487 ,fUseExtraTracks(copy.fUseExtraTracks)
488 ,fUseExtraTracksBgr(copy.fUseExtraTracksBgr)
489 ,fCutFractionPtEmbedded(copy.fCutFractionPtEmbedded)
490 ,fUseEmbeddedJetAxis(copy.fUseEmbeddedJetAxis)
491 ,fUseEmbeddedJetPt(copy.fUseEmbeddedJetPt)
492 ,fJetPtCut(copy.fJetPtCut)
493 ,fJetEtaMin(copy.fJetEtaMin)
494 ,fJetEtaMax(copy.fJetEtaMax)
495 ,fJetPhiMin(copy.fJetPhiMin)
496 ,fJetPhiMax(copy.fJetPhiMax)
497 ,fFFRadius(copy.fFFRadius)
498 ,fFFMinLTrackPt(copy.fFFMinLTrackPt)
499 ,fFFMaxTrackPt(copy.fFFMaxTrackPt)
500 ,fFFMinnTracks(copy.fFFMinnTracks)
501 ,fFFBckgRadius(copy.fFFBckgRadius)
502 ,fBckgMode(copy.fBckgMode)
503 ,fQAMode(copy.fQAMode)
504 ,fFFMode(copy.fFFMode)
505 ,fEffMode(copy.fEffMode)
506 ,fJSMode(copy.fJSMode)
507 ,fAvgTrials(copy.fAvgTrials)
508 ,fTracksRecCuts(copy.fTracksRecCuts)
509 ,fTracksGen(copy.fTracksGen)
510 ,fTracksAODMCCharged(copy.fTracksAODMCCharged)
511 ,fTracksAODMCChargedSecNS(copy.fTracksAODMCChargedSecNS)
512 ,fTracksAODMCChargedSecS(copy.fTracksAODMCChargedSecS)
513 ,fTracksRecQualityCuts(copy.fTracksRecQualityCuts)
514 ,fJetsRec(copy.fJetsRec)
515 ,fJetsRecCuts(copy.fJetsRecCuts)
516 ,fJetsGen(copy.fJetsGen)
517 ,fJetsRecEff(copy.fJetsRecEff)
518 ,fJetsEmbedded(copy.fJetsEmbedded)
519 ,fBckgJetsRec(copy.fBckgJetsRec)
520 ,fBckgJetsRecCuts(copy.fBckgJetsRecCuts)
521 ,fBckgJetsGen(copy.fBckgJetsGen)
522 ,fQATrackHistosRecCuts(copy.fQATrackHistosRecCuts)
523 ,fQATrackHistosGen(copy.fQATrackHistosGen)
524 ,fQAJetHistosRec(copy.fQAJetHistosRec)
525 ,fQAJetHistosRecCuts(copy.fQAJetHistosRecCuts)
526 ,fQAJetHistosRecCutsLeading(copy.fQAJetHistosRecCutsLeading)
527 ,fQAJetHistosGen(copy.fQAJetHistosGen)
528 ,fQAJetHistosGenLeading(copy.fQAJetHistosGenLeading)
529 ,fQAJetHistosRecEffLeading(copy.fQAJetHistosRecEffLeading)
530 ,fFFHistosRecCuts(copy.fFFHistosRecCuts)
37f35567 531 ,fFFHistosRecCutsInc(copy.fFFHistosRecCutsInc)
532 ,fFFHistosRecLeadingTrack(copy.fFFHistosRecLeadingTrack)
ffa175ff 533 ,fFFHistosGen(copy.fFFHistosGen)
37f35567 534 ,fFFHistosGenInc(copy.fFFHistosGenInc)
535 ,fFFHistosGenLeadingTrack(copy.fFFHistosGenLeadingTrack)
ffa175ff 536 ,fQATrackHighPtThreshold(copy.fQATrackHighPtThreshold)
537 ,fFFNBinsJetPt(copy.fFFNBinsJetPt)
538 ,fFFJetPtMin(copy.fFFJetPtMin)
539 ,fFFJetPtMax(copy.fFFJetPtMax)
540 ,fFFNBinsPt(copy.fFFNBinsPt)
541 ,fFFPtMin(copy.fFFPtMin)
542 ,fFFPtMax(copy.fFFPtMax)
543 ,fFFNBinsXi(copy.fFFNBinsXi)
544 ,fFFXiMin(copy.fFFXiMin)
545 ,fFFXiMax(copy.fFFXiMax)
546 ,fFFNBinsZ(copy.fFFNBinsZ)
547 ,fFFZMin(copy.fFFZMin)
548 ,fFFZMax(copy.fFFZMax)
549 ,fQAJetNBinsPt(copy.fQAJetNBinsPt)
550 ,fQAJetPtMin(copy.fQAJetPtMin)
551 ,fQAJetPtMax(copy.fQAJetPtMax)
552 ,fQAJetNBinsEta(copy.fQAJetNBinsEta)
553 ,fQAJetEtaMin(copy.fQAJetEtaMin)
554 ,fQAJetEtaMax(copy.fQAJetEtaMax)
555 ,fQAJetNBinsPhi(copy.fQAJetNBinsPhi)
556 ,fQAJetPhiMin(copy.fQAJetPhiMin)
557 ,fQAJetPhiMax(copy.fQAJetPhiMax)
558 ,fQATrackNBinsPt(copy.fQATrackNBinsPt)
559 ,fQATrackPtMin(copy.fQATrackPtMin)
560 ,fQATrackPtMax(copy.fQATrackPtMax)
561 ,fQATrackNBinsEta(copy.fQATrackNBinsEta)
562 ,fQATrackEtaMin(copy.fQATrackEtaMin)
563 ,fQATrackEtaMax(copy.fQATrackEtaMax)
564 ,fQATrackNBinsPhi(copy.fQATrackNBinsPhi)
565 ,fQATrackPhiMin(copy.fQATrackPhiMin)
566 ,fQATrackPhiMax(copy.fQATrackPhiMax)
567 ,fCommonHistList(copy.fCommonHistList)
568 ,fh1EvtSelection(copy.fh1EvtSelection)
569 ,fh1VertexNContributors(copy.fh1VertexNContributors)
570 ,fh1VertexZ(copy.fh1VertexZ)
571 ,fh1EvtMult(copy.fh1EvtMult)
572 ,fh1EvtCent(copy.fh1EvtCent)
573 ,fh1Xsec(copy.fh1Xsec)
574 ,fh1Trials(copy.fh1Trials)
575 ,fh1PtHard(copy.fh1PtHard)
576 ,fh1PtHardTrials(copy.fh1PtHardTrials)
577 ,fh1nRecJetsCuts(copy.fh1nRecJetsCuts)
578 ,fh1nGenJets(copy.fh1nGenJets)
579 ,fh1nRecEffJets(copy.fh1nRecEffJets)
580 ,fh1nEmbeddedJets(copy.fh1nEmbeddedJets)
581 ,fh1nRecBckgJetsCuts(copy.fh1nRecBckgJetsCuts)
582 ,fh1nGenBckgJets(copy.fh1nGenBckgJets)
583 ,fh2PtRecVsGenPrim(copy.fh2PtRecVsGenPrim)
584 ,fh2PtRecVsGenSec(copy.fh2PtRecVsGenSec)
585 ,fQATrackHistosRecEffGen(copy.fQATrackHistosRecEffGen)
586 ,fQATrackHistosRecEffRec(copy.fQATrackHistosRecEffRec)
587 ,fQATrackHistosSecRecNS(copy.fQATrackHistosSecRecNS)
588 ,fQATrackHistosSecRecS(copy.fQATrackHistosSecRecS)
589 ,fQATrackHistosSecRecSsc(copy.fQATrackHistosSecRecSsc)
590 ,fFFHistosRecEffRec(copy.fFFHistosRecEffRec)
591 ,fFFHistosSecRecNS(copy.fFFHistosSecRecNS)
592 ,fFFHistosSecRecS(copy.fFFHistosSecRecS)
593 ,fFFHistosSecRecSsc(copy.fFFHistosSecRecSsc)
594 // Background
595 ,fh1BckgMult0(copy.fh1BckgMult0)
596 ,fh1BckgMult1(copy.fh1BckgMult1)
597 ,fh1BckgMult2(copy.fh1BckgMult2)
598 ,fh1BckgMult3(copy.fh1BckgMult3)
599 ,fh1BckgMult4(copy.fh1BckgMult4)
600 ,fh1FractionPtEmbedded(copy.fh1FractionPtEmbedded)
601 ,fh1IndexEmbedded(copy.fh1IndexEmbedded)
602 ,fh2DeltaPtVsJetPtEmbedded(copy.fh2DeltaPtVsJetPtEmbedded)
603 ,fh2DeltaPtVsRecJetPtEmbedded(copy.fh2DeltaPtVsRecJetPtEmbedded)
604 ,fh1DeltaREmbedded(copy.fh1DeltaREmbedded)
605 ,fQABckgHisto0RecCuts(copy.fQABckgHisto0RecCuts)
606 ,fQABckgHisto0Gen(copy.fQABckgHisto0Gen)
607 ,fQABckgHisto1RecCuts(copy.fQABckgHisto1RecCuts)
608 ,fQABckgHisto1Gen(copy.fQABckgHisto1Gen)
609 ,fQABckgHisto2RecCuts(copy.fQABckgHisto2RecCuts)
610 ,fQABckgHisto2Gen(copy.fQABckgHisto2Gen)
611 ,fQABckgHisto3RecCuts(copy.fQABckgHisto3RecCuts)
612 ,fQABckgHisto3Gen(copy.fQABckgHisto3Gen)
613 ,fQABckgHisto4RecCuts(copy.fQABckgHisto4RecCuts)
614 ,fQABckgHisto4Gen(copy.fQABckgHisto4Gen)
615 ,fFFBckgHisto0RecCuts(copy.fFFBckgHisto0RecCuts)
616 ,fFFBckgHisto0Gen(copy.fFFBckgHisto0Gen)
617 ,fFFBckgHisto1RecCuts(copy.fFFBckgHisto1RecCuts)
618 ,fFFBckgHisto1Gen(copy.fFFBckgHisto1Gen)
619 ,fFFBckgHisto2RecCuts(copy.fFFBckgHisto2RecCuts)
620 ,fFFBckgHisto2Gen(copy.fFFBckgHisto2Gen)
621 ,fFFBckgHisto3RecCuts(copy.fFFBckgHisto3RecCuts)
622 ,fFFBckgHisto3Gen(copy.fFFBckgHisto3Gen)
623 ,fFFBckgHisto4RecCuts(copy.fFFBckgHisto4RecCuts)
624 ,fFFBckgHisto4Gen(copy.fFFBckgHisto4Gen)
625 ,fFFBckgHisto0RecEffRec(copy.fFFBckgHisto0RecEffRec)
626 ,fFFBckgHisto0SecRecNS(copy.fFFBckgHisto0SecRecNS)
627 ,fFFBckgHisto0SecRecS(copy.fFFBckgHisto0SecRecS)
628 ,fFFBckgHisto0SecRecSsc(copy.fFFBckgHisto0SecRecSsc)
629 // jet shape
630 ,fProNtracksLeadingJet(copy.fProNtracksLeadingJet)
631 ,fProDelR80pcPt(copy.fProDelR80pcPt)
632 ,fProNtracksLeadingJetGen(copy.fProNtracksLeadingJetGen)
633 ,fProDelR80pcPtGen(copy.fProDelR80pcPtGen)
634 ,fProNtracksLeadingJetBgrPerp2(copy.fProNtracksLeadingJetBgrPerp2)
635 ,fProNtracksLeadingJetRecPrim(copy.fProNtracksLeadingJetRecPrim)
636 ,fProDelR80pcPtRecPrim(copy.fProDelR80pcPtRecPrim)
637 ,fProNtracksLeadingJetRecSecNS(copy.fProNtracksLeadingJetRecSecNS)
638 ,fProNtracksLeadingJetRecSecS(copy.fProNtracksLeadingJetRecSecS)
639 ,fProNtracksLeadingJetRecSecSsc(copy.fProNtracksLeadingJetRecSecSsc)
640 ,fRandom(copy.fRandom)
641{
642 // copy constructor
643 fBckgType[0] = copy.fBckgType[0];
644 fBckgType[1] = copy.fBckgType[1];
645 fBckgType[2] = copy.fBckgType[2];
646 fBckgType[3] = copy.fBckgType[3];
647 fBckgType[4] = copy.fBckgType[4];
648
649
650 for(Int_t ii=0; ii<5; ii++){
651 fProDelRPtSum[ii] = copy.fProDelRPtSum[ii];
652 fProDelRPtSumGen[ii] = copy.fProDelRPtSumGen[ii];
653 fProDelRPtSumBgrPerp2[ii] = copy.fProDelRPtSumBgrPerp2[ii];
654 fProDelRPtSumRecPrim[ii] = copy.fProDelRPtSumRecPrim[ii];
655 fProDelRPtSumRecSecNS[ii] = copy.fProDelRPtSumRecSecNS[ii];
656 fProDelRPtSumRecSecS[ii] = copy.fProDelRPtSumRecSecS[ii];
657 fProDelRPtSumRecSecSsc[ii] = copy.fProDelRPtSumRecSecSsc[ii];
658 }
659}
660
661// _________________________________________________________________________________________________________________________________
662AliAnalysisTaskFragmentationFunction& AliAnalysisTaskFragmentationFunction::operator=(const AliAnalysisTaskFragmentationFunction& o)
663{
664 // assignment
665
666 if(this!=&o){
667
668 AliAnalysisTaskSE::operator=(o);
669 fESD = o.fESD;
670 fAOD = o.fAOD;
671 fAODJets = o.fAODJets;
672 fAODExtension = o.fAODExtension;
673 fNonStdFile = o.fNonStdFile;
674 fBranchRecJets = o.fBranchRecJets;
675 fBranchRecBckgClusters = o.fBranchRecBckgClusters;
676 fBranchGenJets = o.fBranchGenJets;
677 fBranchEmbeddedJets = o.fBranchEmbeddedJets;
678 fTrackTypeGen = o.fTrackTypeGen;
679 fJetTypeGen = o.fJetTypeGen;
680 fJetTypeRecEff = o.fJetTypeRecEff;
681 fUseAODInputJets = o.fUseAODInputJets;
682 fFilterMask = o.fFilterMask;
683 fUsePhysicsSelection = o.fUsePhysicsSelection;
684 fEvtSelectionMask = o.fEvtSelectionMask;
685 fEventClass = o.fEventClass;
686 fMaxVertexZ = o.fMaxVertexZ;
687 fTrackPtCut = o.fTrackPtCut;
688 fTrackEtaMin = o.fTrackEtaMin;
689 fTrackEtaMax = o.fTrackEtaMax;
690 fTrackPhiMin = o.fTrackPhiMin;
691 fTrackPhiMax = o.fTrackPhiMax;
692 fUseExtraTracks = o.fUseExtraTracks;
693 fUseExtraTracksBgr = o.fUseExtraTracksBgr;
694 fCutFractionPtEmbedded = o.fCutFractionPtEmbedded;
695 fUseEmbeddedJetAxis = o.fUseEmbeddedJetAxis;
696 fUseEmbeddedJetPt = o.fUseEmbeddedJetPt;
697 fJetPtCut = o.fJetPtCut;
698 fJetEtaMin = o.fJetEtaMin;
699 fJetEtaMax = o.fJetEtaMax;
700 fJetPhiMin = o.fJetPhiMin;
701 fJetPhiMax = o.fJetPhiMin;
702 fFFRadius = o.fFFRadius;
703 fFFMinLTrackPt = o.fFFMinLTrackPt;
704 fFFMaxTrackPt = o.fFFMaxTrackPt;
705 fFFMinnTracks = o.fFFMinnTracks;
706 fFFBckgRadius = o.fFFBckgRadius;
707 fBckgMode = o.fBckgMode;
708 fQAMode = o.fQAMode;
709 fFFMode = o.fFFMode;
710 fEffMode = o.fEffMode;
711 fJSMode = o.fJSMode;
712 fBckgType[0] = o.fBckgType[0];
713 fBckgType[1] = o.fBckgType[1];
714 fBckgType[2] = o.fBckgType[2];
715 fBckgType[3] = o.fBckgType[3];
716 fBckgType[4] = o.fBckgType[4];
717 fAvgTrials = o.fAvgTrials;
718 fTracksRecCuts = o.fTracksRecCuts;
719 fTracksGen = o.fTracksGen;
720 fTracksAODMCCharged = o.fTracksAODMCCharged;
721 fTracksAODMCChargedSecNS = o.fTracksAODMCChargedSecNS;
722 fTracksAODMCChargedSecS = o.fTracksAODMCChargedSecS;
723 fTracksRecQualityCuts = o.fTracksRecQualityCuts;
724 fJetsRec = o.fJetsRec;
725 fJetsRecCuts = o.fJetsRecCuts;
726 fJetsGen = o.fJetsGen;
727 fJetsRecEff = o.fJetsRecEff;
728 fJetsEmbedded = o.fJetsEmbedded;
729 fBckgJetsRec = o.fBckgJetsRec;
730 fBckgJetsRecCuts = o.fBckgJetsRecCuts;
731 fBckgJetsGen = o.fBckgJetsGen;
732 fQATrackHistosRecCuts = o.fQATrackHistosRecCuts;
733 fQATrackHistosGen = o.fQATrackHistosGen;
734 fQAJetHistosRec = o.fQAJetHistosRec;
735 fQAJetHistosRecCuts = o.fQAJetHistosRecCuts;
736 fQAJetHistosRecCutsLeading = o.fQAJetHistosRecCutsLeading;
737 fQAJetHistosGen = o.fQAJetHistosGen;
738 fQAJetHistosGenLeading = o.fQAJetHistosGenLeading;
739 fQAJetHistosRecEffLeading = o.fQAJetHistosRecEffLeading;
740 fFFHistosRecCuts = o.fFFHistosRecCuts;
37f35567 741 fFFHistosRecCutsInc = o.fFFHistosRecCutsInc;
742 fFFHistosRecLeadingTrack = o.fFFHistosRecLeadingTrack;
ffa175ff 743 fFFHistosGen = o.fFFHistosGen;
37f35567 744 fFFHistosGenInc = o.fFFHistosGenInc;
745 fFFHistosGenLeadingTrack = o.fFFHistosGenLeadingTrack;
ffa175ff 746 fQATrackHighPtThreshold = o.fQATrackHighPtThreshold;
747 fFFNBinsJetPt = o.fFFNBinsJetPt;
748 fFFJetPtMin = o.fFFJetPtMin;
749 fFFJetPtMax = o.fFFJetPtMax;
750 fFFNBinsPt = o.fFFNBinsPt;
751 fFFPtMin = o.fFFPtMin;
752 fFFPtMax = o.fFFPtMax;
753 fFFNBinsXi = o.fFFNBinsXi;
754 fFFXiMin = o.fFFXiMin;
755 fFFXiMax = o.fFFXiMax;
756 fFFNBinsZ = o.fFFNBinsZ;
757 fFFZMin = o.fFFZMin;
758 fFFZMax = o.fFFZMax;
759 fQAJetNBinsPt = o.fQAJetNBinsPt;
760 fQAJetPtMin = o.fQAJetPtMin;
761 fQAJetPtMax = o.fQAJetPtMax;
762 fQAJetNBinsEta = o.fQAJetNBinsEta;
763 fQAJetEtaMin = o.fQAJetEtaMin;
764 fQAJetEtaMax = o.fQAJetEtaMax;
765 fQAJetNBinsPhi = o.fQAJetNBinsPhi;
766 fQAJetPhiMin = o.fQAJetPhiMin;
767 fQAJetPhiMax = o.fQAJetPhiMax;
768 fQATrackNBinsPt = o.fQATrackNBinsPt;
769 fQATrackPtMin = o.fQATrackPtMin;
770 fQATrackPtMax = o.fQATrackPtMax;
771 fQATrackNBinsEta = o.fQATrackNBinsEta;
772 fQATrackEtaMin = o.fQATrackEtaMin;
773 fQATrackEtaMax = o.fQATrackEtaMax;
774 fQATrackNBinsPhi = o.fQATrackNBinsPhi;
775 fQATrackPhiMin = o.fQATrackPhiMin;
776 fQATrackPhiMax = o.fQATrackPhiMax;
777 fCommonHistList = o.fCommonHistList;
778 fh1EvtSelection = o.fh1EvtSelection;
779 fh1VertexNContributors = o.fh1VertexNContributors;
780 fh1VertexZ = o.fh1VertexZ;
781 fh1EvtMult = o.fh1EvtMult;
782 fh1EvtCent = o.fh1EvtCent;
783 fh1Xsec = o.fh1Xsec;
784 fh1Trials = o.fh1Trials;
785 fh1PtHard = o.fh1PtHard;
786 fh1PtHardTrials = o.fh1PtHardTrials;
787 fh1nRecJetsCuts = o.fh1nRecJetsCuts;
788 fh1nGenJets = o.fh1nGenJets;
789 fh1nRecEffJets = o.fh1nRecEffJets;
790 fh1nEmbeddedJets = o.fh1nEmbeddedJets;
791 fh2PtRecVsGenPrim = o.fh2PtRecVsGenPrim;
792 fh2PtRecVsGenSec = o.fh2PtRecVsGenSec;
793 fQATrackHistosRecEffGen = o.fQATrackHistosRecEffGen;
794 fQATrackHistosRecEffRec = o.fQATrackHistosRecEffRec;
795 fQATrackHistosSecRecNS = o.fQATrackHistosSecRecNS;
796 fQATrackHistosSecRecS = o.fQATrackHistosSecRecS;
797 fQATrackHistosSecRecSsc = o.fQATrackHistosSecRecSsc;
798 fFFHistosRecEffRec = o.fFFHistosRecEffRec;
799 fFFHistosSecRecNS = o.fFFHistosSecRecNS;
800 fFFHistosSecRecS = o.fFFHistosSecRecS;
801 fFFHistosSecRecSsc = o.fFFHistosSecRecSsc;
802 // Background
803 fh1BckgMult0 = o.fh1BckgMult0;
804 fh1BckgMult1 = o.fh1BckgMult1;
805 fh1BckgMult2 = o.fh1BckgMult2;
806 fh1BckgMult3 = o.fh1BckgMult3;
807 fh1BckgMult4 = o.fh1BckgMult4;
808 fh1FractionPtEmbedded = o.fh1FractionPtEmbedded;
809 fh1IndexEmbedded = o.fh1IndexEmbedded;
810 fh2DeltaPtVsJetPtEmbedded = o.fh2DeltaPtVsJetPtEmbedded;
811 fh2DeltaPtVsRecJetPtEmbedded = o.fh2DeltaPtVsRecJetPtEmbedded;
812 fh1DeltaREmbedded = o.fh1DeltaREmbedded;
813 fQABckgHisto0RecCuts = o.fQABckgHisto0RecCuts;
814 fQABckgHisto0Gen = o.fQABckgHisto0Gen;
815 fQABckgHisto1RecCuts = o.fQABckgHisto1RecCuts;
816 fQABckgHisto1Gen = o.fQABckgHisto1Gen;
817 fQABckgHisto2RecCuts = o.fQABckgHisto2RecCuts;
818 fQABckgHisto2Gen = o.fQABckgHisto2Gen;
819 fQABckgHisto3RecCuts = o.fQABckgHisto3RecCuts;
820 fQABckgHisto3Gen = o.fQABckgHisto3Gen;
821 fQABckgHisto4RecCuts = o.fQABckgHisto4RecCuts;
822 fQABckgHisto4Gen = o.fQABckgHisto4Gen;
823 fFFBckgHisto0RecCuts = o.fFFBckgHisto0RecCuts;
824 fFFBckgHisto0Gen = o.fFFBckgHisto0Gen;
825 fFFBckgHisto1RecCuts = o.fFFBckgHisto1RecCuts;
826 fFFBckgHisto1Gen = o.fFFBckgHisto1Gen;
827 fFFBckgHisto2RecCuts = o.fFFBckgHisto2RecCuts;
828 fFFBckgHisto2Gen = o.fFFBckgHisto2Gen;
829 fFFBckgHisto3RecCuts = o.fFFBckgHisto3RecCuts;
830 fFFBckgHisto3Gen = o.fFFBckgHisto3Gen;
831 fFFBckgHisto4RecCuts = o.fFFBckgHisto4RecCuts;
832 fFFBckgHisto4Gen = o.fFFBckgHisto4Gen;
833 fFFBckgHisto0RecEffRec = o.fFFBckgHisto0RecEffRec;
834 fFFBckgHisto0SecRecNS = o.fFFBckgHisto0SecRecNS;
835 fFFBckgHisto0SecRecS = o.fFFBckgHisto0SecRecS;
836 fFFBckgHisto0SecRecSsc = o.fFFBckgHisto0SecRecSsc;
837 fProNtracksLeadingJet = o.fProNtracksLeadingJet;
838 fProDelR80pcPt = o.fProDelR80pcPt;
839 fProNtracksLeadingJetGen = o.fProNtracksLeadingJetGen;
840 fProDelR80pcPtGen = o.fProDelR80pcPtGen;
841 fProNtracksLeadingJetBgrPerp2 = o.fProNtracksLeadingJetBgrPerp2;
842 fProNtracksLeadingJetRecPrim = o.fProNtracksLeadingJetRecPrim;
843 fProDelR80pcPtRecPrim = o.fProDelR80pcPtRecPrim;
844 fProNtracksLeadingJetRecSecNS = o.fProNtracksLeadingJetRecSecNS;
845 fProNtracksLeadingJetRecSecS = o.fProNtracksLeadingJetRecSecS;
846 fProNtracksLeadingJetRecSecSsc = o.fProNtracksLeadingJetRecSecSsc;
847 fRandom = o.fRandom;
848
849 for(Int_t ii=0; ii<5; ii++){
850 fProDelRPtSum[ii] = o.fProDelRPtSum[ii];
851 fProDelRPtSumGen[ii] = o.fProDelRPtSumGen[ii];
852 fProDelRPtSumBgrPerp2[ii] = o.fProDelRPtSumBgrPerp2[ii];
853 fProDelRPtSumRecPrim[ii] = o.fProDelRPtSumRecPrim[ii];
854 fProDelRPtSumRecSecNS[ii] = o.fProDelRPtSumRecSecNS[ii];
855 fProDelRPtSumRecSecS[ii] = o.fProDelRPtSumRecSecS[ii];
856 fProDelRPtSumRecSecSsc[ii] = o.fProDelRPtSumRecSecSsc[ii];
857 }
858 }
859
860 return *this;
861}
862
863//___________________________________________________________________________
864AliAnalysisTaskFragmentationFunction::~AliAnalysisTaskFragmentationFunction()
865{
866 // destructor
867
868 if(fTracksRecCuts) delete fTracksRecCuts;
869 if(fTracksGen) delete fTracksGen;
870 if(fTracksAODMCCharged) delete fTracksAODMCCharged;
871 if(fTracksAODMCChargedSecNS) delete fTracksAODMCChargedSecNS;
872 if(fTracksAODMCChargedSecS) delete fTracksAODMCChargedSecS;
873 if(fTracksRecQualityCuts) delete fTracksRecQualityCuts;
874 if(fJetsRec) delete fJetsRec;
875 if(fJetsRecCuts) delete fJetsRecCuts;
876 if(fJetsGen) delete fJetsGen;
877 if(fJetsRecEff) delete fJetsRecEff;
878 if(fJetsEmbedded) delete fJetsEmbedded;
879
880 if(fBckgMode &&
881 (fBckgType[0]==kBckgClusters || fBckgType[1]==kBckgClusters || fBckgType[2]==kBckgClusters || fBckgType[3]==kBckgClusters || fBckgType[4]==kBckgClusters ||
882 fBckgType[0]==kBckgClustersOutLeading || fBckgType[1]==kBckgClustersOutLeading || fBckgType[2]==kBckgClustersOutLeading ||
883 fBckgType[3]==kBckgClustersOutLeading || fBckgType[4]==kBckgClustersOutLeading)){
884
885 if(fBckgJetsRec) delete fBckgJetsRec;
886 if(fBckgJetsRecCuts) delete fBckgJetsRecCuts;
887 if(fBckgJetsGen) delete fBckgJetsGen;
888 }
889 if(fRandom) delete fRandom;
890}
891
892//______________________________________________________________________________________________________
893AliAnalysisTaskFragmentationFunction::AliFragFuncHistos::AliFragFuncHistos(const char* name,
894 Int_t nJetPt, Float_t jetPtMin, Float_t jetPtMax,
895 Int_t nPt, Float_t ptMin, Float_t ptMax,
896 Int_t nXi, Float_t xiMin, Float_t xiMax,
897 Int_t nZ , Float_t zMin , Float_t zMax)
898 : TObject()
899 ,fNBinsJetPt(nJetPt)
900 ,fJetPtMin(jetPtMin)
901 ,fJetPtMax(jetPtMax)
902 ,fNBinsPt(nPt)
903 ,fPtMin(ptMin)
904 ,fPtMax(ptMax)
905 ,fNBinsXi(nXi)
906 ,fXiMin(xiMin)
907 ,fXiMax(xiMax)
908 ,fNBinsZ(nZ)
909 ,fZMin(zMin)
910 ,fZMax(zMax)
911 ,fh2TrackPt(0)
912 ,fh2Xi(0)
913 ,fh2Z(0)
914 ,fh1JetPt(0)
915 ,fNameFF(name)
916{
917 // default constructor
918
919}
920
921//___________________________________________________________________________
922AliAnalysisTaskFragmentationFunction::AliFragFuncHistos::AliFragFuncHistos(const AliFragFuncHistos& copy)
923 : TObject()
924 ,fNBinsJetPt(copy.fNBinsJetPt)
925 ,fJetPtMin(copy.fJetPtMin)
926 ,fJetPtMax(copy.fJetPtMax)
927 ,fNBinsPt(copy.fNBinsPt)
928 ,fPtMin(copy.fPtMin)
929 ,fPtMax(copy.fPtMax)
930 ,fNBinsXi(copy.fNBinsXi)
931 ,fXiMin(copy.fXiMin)
932 ,fXiMax(copy.fXiMax)
933 ,fNBinsZ(copy.fNBinsZ)
934 ,fZMin(copy.fZMin)
935 ,fZMax(copy.fZMax)
936 ,fh2TrackPt(copy.fh2TrackPt)
937 ,fh2Xi(copy.fh2Xi)
938 ,fh2Z(copy.fh2Z)
939 ,fh1JetPt(copy.fh1JetPt)
940 ,fNameFF(copy.fNameFF)
941{
942 // copy constructor
943}
944
945//_______________________________________________________________________________________________________________________________________________________________
946AliAnalysisTaskFragmentationFunction::AliFragFuncHistos& AliAnalysisTaskFragmentationFunction::AliFragFuncHistos::operator=(const AliAnalysisTaskFragmentationFunction::AliFragFuncHistos& o)
947{
948 // assignment
949
950 if(this!=&o){
951 TObject::operator=(o);
952 fNBinsJetPt = o.fNBinsJetPt;
953 fJetPtMin = o.fJetPtMin;
954 fJetPtMax = o.fJetPtMax;
955 fNBinsPt = o.fNBinsPt;
956 fPtMin = o.fPtMin;
957 fPtMax = o.fPtMax;
958 fNBinsXi = o.fNBinsXi;
959 fXiMin = o.fXiMin;
960 fXiMax = o.fXiMax;
961 fNBinsZ = o.fNBinsZ;
962 fZMin = o.fZMin;
963 fZMax = o.fZMax;
964 fh2TrackPt = o.fh2TrackPt;
965 fh2Xi = o.fh2Xi;
966 fh2Z = o.fh2Z;
967 fh1JetPt = o.fh1JetPt;
968 fNameFF = o.fNameFF;
969 }
970
971 return *this;
972}
973
974//_________________________________________________________
975AliAnalysisTaskFragmentationFunction::AliFragFuncHistos::~AliFragFuncHistos()
976{
977 // destructor
978
979 if(fh1JetPt) delete fh1JetPt;
980 if(fh2TrackPt) delete fh2TrackPt;
981 if(fh2Xi) delete fh2Xi;
982 if(fh2Z) delete fh2Z;
983}
984
985//_________________________________________________________________
986void AliAnalysisTaskFragmentationFunction::AliFragFuncHistos::DefineHistos()
987{
988 // book FF histos
989
990 fh1JetPt = new TH1F(Form("fh1FFJetPt%s", fNameFF.Data()),"",fNBinsJetPt,fJetPtMin,fJetPtMax);
991 fh2TrackPt = new TH2F(Form("fh2FFTrackPt%s",fNameFF.Data()),"",fNBinsJetPt, fJetPtMin, fJetPtMax,fNBinsPt, fPtMin, fPtMax);
992 fh2Z = new TH2F(Form("fh2FFZ%s",fNameFF.Data()),"",fNBinsJetPt, fJetPtMin, fJetPtMax, fNBinsZ, fZMin, fZMax);
993 fh2Xi = new TH2F(Form("fh2FFXi%s",fNameFF.Data()),"",fNBinsJetPt, fJetPtMin, fJetPtMax, fNBinsXi, fXiMin, fXiMax);
994
995 AliAnalysisTaskFragmentationFunction::SetProperties(fh1JetPt, "p_{T} [GeV/c]", "entries");
996 AliAnalysisTaskFragmentationFunction::SetProperties(fh2TrackPt,"jet p_{T} [GeV/c]","p_{T} [GeV/c]","entries");
997 AliAnalysisTaskFragmentationFunction::SetProperties(fh2Xi,"jet p_{T} [GeV/c]","#xi", "entries");
998 AliAnalysisTaskFragmentationFunction::SetProperties(fh2Z,"jet p_{T} [GeV/c]","z","entries");
999}
1000
1001//_______________________________________________________________________________________________________________
1002void AliAnalysisTaskFragmentationFunction::AliFragFuncHistos::FillFF(Float_t trackPt, Float_t jetPt, Bool_t incrementJetPt, Float_t norm,
1003 Bool_t scaleStrangeness, Float_t scaleFacStrangeness)
1004{
1005 // fill FF
1006
1007 if(incrementJetPt && norm) fh1JetPt->Fill(jetPt,1/norm);
1008 else if(incrementJetPt) fh1JetPt->Fill(jetPt);
1009
1010 // Added for proper normalization of FF background estimation
1011 // when zero track are found in the background region
1012 if((int)trackPt==-1) return;
1013
1014 if(norm)fh2TrackPt->Fill(jetPt,trackPt,1/norm);
1015 else if(scaleStrangeness) fh2TrackPt->Fill(jetPt,trackPt,scaleFacStrangeness);
1016 else fh2TrackPt->Fill(jetPt,trackPt);
1017
1018 Double_t z = 0.;
1019 if(jetPt>0) z = trackPt / jetPt;
1020 Double_t xi = 0;
1021 if(z>0) xi = TMath::Log(1/z);
1022
1023 if(trackPt>(1-1e-06)*jetPt && trackPt<(1+1e-06)*jetPt){ // case z=1 : move entry to last histo bin <1
1024 z = 1-1e-06;
1025 xi = 1e-06;
1026 }
1027
1028
1029 if(norm){
1030 fh2Xi->Fill(jetPt,xi,1/norm);
1031 fh2Z->Fill(jetPt,z,1/norm);
1032 }
1033 else if(scaleStrangeness){
1034 fh2Xi->Fill(jetPt,xi,scaleFacStrangeness);
1035 fh2Z->Fill(jetPt,z,scaleFacStrangeness);
1036 }
1037 else {
1038 fh2Xi->Fill(jetPt,xi);
1039 fh2Z->Fill(jetPt,z);
1040 }
1041}
1042
1043//_________________________________________________________________________________
1044void AliAnalysisTaskFragmentationFunction::AliFragFuncHistos::AddToOutput(TList* list) const
1045{
1046 // add histos to list
1047
1048 list->Add(fh1JetPt);
1049
1050 list->Add(fh2TrackPt);
1051 list->Add(fh2Xi);
1052 list->Add(fh2Z);
1053}
1054
1055//_________________________________________________________________________________________________________
1056AliAnalysisTaskFragmentationFunction::AliFragFuncQAJetHistos::AliFragFuncQAJetHistos(const char* name,
1057 Int_t nPt, Float_t ptMin, Float_t ptMax,
1058 Int_t nEta, Float_t etaMin, Float_t etaMax,
1059 Int_t nPhi, Float_t phiMin, Float_t phiMax)
1060 : TObject()
1061 ,fNBinsPt(nPt)
1062 ,fPtMin(ptMin)
1063 ,fPtMax(ptMax)
1064 ,fNBinsEta(nEta)
1065 ,fEtaMin(etaMin)
1066 ,fEtaMax(etaMax)
1067 ,fNBinsPhi(nPhi)
1068 ,fPhiMin(phiMin)
1069 ,fPhiMax(phiMax)
1070 ,fh2EtaPhi(0)
1071 ,fh1Pt(0)
1072 ,fNameQAJ(name)
1073{
1074 // default constructor
1075}
1076
1077//____________________________________________________________________________________
1078AliAnalysisTaskFragmentationFunction::AliFragFuncQAJetHistos::AliFragFuncQAJetHistos(const AliFragFuncQAJetHistos& copy)
1079 : TObject()
1080 ,fNBinsPt(copy.fNBinsPt)
1081 ,fPtMin(copy.fPtMin)
1082 ,fPtMax(copy.fPtMax)
1083 ,fNBinsEta(copy.fNBinsEta)
1084 ,fEtaMin(copy.fEtaMin)
1085 ,fEtaMax(copy.fEtaMax)
1086 ,fNBinsPhi(copy.fNBinsPhi)
1087 ,fPhiMin(copy.fPhiMin)
1088 ,fPhiMax(copy.fPhiMax)
1089 ,fh2EtaPhi(copy.fh2EtaPhi)
1090 ,fh1Pt(copy.fh1Pt)
1091 ,fNameQAJ(copy.fNameQAJ)
1092{
1093 // copy constructor
1094}
1095
1096//________________________________________________________________________________________________________________________________________________________________________
1097AliAnalysisTaskFragmentationFunction::AliFragFuncQAJetHistos& AliAnalysisTaskFragmentationFunction::AliFragFuncQAJetHistos::operator=(const AliAnalysisTaskFragmentationFunction::AliFragFuncQAJetHistos& o)
1098{
1099 // assignment
1100
1101 if(this!=&o){
1102 TObject::operator=(o);
1103 fNBinsPt = o.fNBinsPt;
1104 fPtMin = o.fPtMin;
1105 fPtMax = o.fPtMax;
1106 fNBinsEta = o.fNBinsEta;
1107 fEtaMin = o.fEtaMin;
1108 fEtaMax = o.fEtaMax;
1109 fNBinsPhi = o.fNBinsPhi;
1110 fPhiMin = o.fPhiMin;
1111 fPhiMax = o.fPhiMax;
1112 fh2EtaPhi = o.fh2EtaPhi;
1113 fh1Pt = o.fh1Pt;
1114 fNameQAJ = o.fNameQAJ;
1115 }
1116
1117 return *this;
1118}
1119
1120//______________________________________________________________
1121AliAnalysisTaskFragmentationFunction::AliFragFuncQAJetHistos::~AliFragFuncQAJetHistos()
1122{
1123 // destructor
1124
1125 if(fh2EtaPhi) delete fh2EtaPhi;
1126 if(fh1Pt) delete fh1Pt;
1127}
1128
1129//____________________________________________________________________
1130void AliAnalysisTaskFragmentationFunction::AliFragFuncQAJetHistos::DefineHistos()
1131{
1132 // book jet QA histos
1133
1134 fh2EtaPhi = new TH2F(Form("fh2JetQAEtaPhi%s", fNameQAJ.Data()), Form("%s: #eta - #phi distribution", fNameQAJ.Data()), fNBinsEta, fEtaMin, fEtaMax, fNBinsPhi, fPhiMin, fPhiMax);
1135 fh1Pt = new TH1F(Form("fh1JetQAPt%s", fNameQAJ.Data()), Form("%s: p_{T} distribution", fNameQAJ.Data()), fNBinsPt, fPtMin, fPtMax);
1136
1137 AliAnalysisTaskFragmentationFunction::SetProperties(fh2EtaPhi, "#eta", "#phi");
1138 AliAnalysisTaskFragmentationFunction::SetProperties(fh1Pt, "p_{T} [GeV/c]", "entries");
1139}
1140
1141//____________________________________________________________________________________________________
1142void AliAnalysisTaskFragmentationFunction::AliFragFuncQAJetHistos::FillJetQA(Float_t eta, Float_t phi, Float_t pt)
1143{
1144 // fill jet QA histos
1145
1146 fh2EtaPhi->Fill( eta, phi);
1147 fh1Pt->Fill( pt );
1148}
1149
1150//____________________________________________________________________________________
1151void AliAnalysisTaskFragmentationFunction::AliFragFuncQAJetHistos::AddToOutput(TList* list) const
1152{
1153 // add histos to list
1154
1155 list->Add(fh2EtaPhi);
1156 list->Add(fh1Pt);
1157}
1158
1159//___________________________________________________________________________________________________________
1160AliAnalysisTaskFragmentationFunction::AliFragFuncQATrackHistos::AliFragFuncQATrackHistos(const char* name,
1161 Int_t nPt, Float_t ptMin, Float_t ptMax,
1162 Int_t nEta, Float_t etaMin, Float_t etaMax,
1163 Int_t nPhi, Float_t phiMin, Float_t phiMax,
1164 Float_t ptThresh)
1165 : TObject()
1166 ,fNBinsPt(nPt)
1167 ,fPtMin(ptMin)
1168 ,fPtMax(ptMax)
1169 ,fNBinsEta(nEta)
1170 ,fEtaMin(etaMin)
1171 ,fEtaMax(etaMax)
1172 ,fNBinsPhi(nPhi)
1173 ,fPhiMin(phiMin)
1174 ,fPhiMax(phiMax)
1175 ,fHighPtThreshold(ptThresh)
1176 ,fh2EtaPhi(0)
1177 ,fh1Pt(0)
1178 ,fh2HighPtEtaPhi(0)
1179 ,fh2PhiPt(0)
1180 ,fNameQAT(name)
1181{
1182 // default constructor
1183}
1184
1185//__________________________________________________________________________________________
1186AliAnalysisTaskFragmentationFunction::AliFragFuncQATrackHistos::AliFragFuncQATrackHistos(const AliFragFuncQATrackHistos& copy)
1187 : TObject()
1188 ,fNBinsPt(copy.fNBinsPt)
1189 ,fPtMin(copy.fPtMin)
1190 ,fPtMax(copy.fPtMax)
1191 ,fNBinsEta(copy.fNBinsEta)
1192 ,fEtaMin(copy.fEtaMin)
1193 ,fEtaMax(copy.fEtaMax)
1194 ,fNBinsPhi(copy.fNBinsPhi)
1195 ,fPhiMin(copy.fPhiMin)
1196 ,fPhiMax(copy.fPhiMax)
1197 ,fHighPtThreshold(copy.fHighPtThreshold)
1198 ,fh2EtaPhi(copy.fh2EtaPhi)
1199 ,fh1Pt(copy.fh1Pt)
1200 ,fh2HighPtEtaPhi(copy.fh2HighPtEtaPhi)
1201 ,fh2PhiPt(copy.fh2PhiPt)
1202 ,fNameQAT(copy.fNameQAT)
1203{
1204 // copy constructor
1205}
1206
1207// _____________________________________________________________________________________________________________________________________________________________________________
1208AliAnalysisTaskFragmentationFunction::AliFragFuncQATrackHistos& AliAnalysisTaskFragmentationFunction::AliFragFuncQATrackHistos::operator=(const AliAnalysisTaskFragmentationFunction::AliFragFuncQATrackHistos& o)
1209{
1210 // assignment
1211
1212 if(this!=&o){
1213 TObject::operator=(o);
1214 fNBinsPt = o.fNBinsPt;
1215 fPtMin = o.fPtMin;
1216 fPtMax = o.fPtMax;
1217 fNBinsEta = o.fNBinsEta;
1218 fEtaMin = o.fEtaMin;
1219 fEtaMax = o.fEtaMax;
1220 fNBinsPhi = o.fNBinsPhi;
1221 fPhiMin = o.fPhiMin;
1222 fPhiMax = o.fPhiMax;
1223 fHighPtThreshold = o.fHighPtThreshold;
1224 fh2EtaPhi = o.fh2EtaPhi;
1225 fh1Pt = o.fh1Pt;
1226 fh2HighPtEtaPhi = o.fh2HighPtEtaPhi;
1227 fh2PhiPt = o.fh2PhiPt;
1228 fNameQAT = o.fNameQAT;
1229 }
1230
1231 return *this;
1232}
1233
1234//___________________________________________________________________
1235AliAnalysisTaskFragmentationFunction::AliFragFuncQATrackHistos::~AliFragFuncQATrackHistos()
1236{
1237 // destructor
1238
1239 if(fh2EtaPhi) delete fh2EtaPhi;
1240 if(fh2HighPtEtaPhi) delete fh2HighPtEtaPhi;
1241 if(fh1Pt) delete fh1Pt;
1242 if(fh2PhiPt) delete fh2PhiPt;
1243}
1244
1245//______________________________________________________________________
1246void AliAnalysisTaskFragmentationFunction::AliFragFuncQATrackHistos::DefineHistos()
1247{
1248 // book track QA histos
1249
1250 fh2EtaPhi = new TH2F(Form("fh2TrackQAEtaPhi%s", fNameQAT.Data()), Form("%s: #eta - #phi distribution", fNameQAT.Data()), fNBinsEta, fEtaMin, fEtaMax, fNBinsPhi, fPhiMin, fPhiMax);
1251 fh2HighPtEtaPhi = new TH2F(Form("fh2TrackQAHighPtEtaPhi%s", fNameQAT.Data()), Form("%s: #eta - #phi distribution for high-p_{T}", fNameQAT.Data()), fNBinsEta, fEtaMin, fEtaMax, fNBinsPhi, fPhiMin, fPhiMax);
1252 fh1Pt = new TH1F(Form("fh1TrackQAPt%s", fNameQAT.Data()), Form("%s: p_{T} distribution", fNameQAT.Data()), fNBinsPt, fPtMin, fPtMax);
1253 fh2PhiPt = new TH2F(Form("fh2TrackQAPhiPt%s", fNameQAT.Data()), Form("%s: #phi - p_{T} distribution", fNameQAT.Data()), fNBinsPhi, fPhiMin, fPhiMax, fNBinsPt, fPtMin, fPtMax);
1254
1255 AliAnalysisTaskFragmentationFunction::SetProperties(fh2EtaPhi, "#eta", "#phi");
1256 AliAnalysisTaskFragmentationFunction::SetProperties(fh2HighPtEtaPhi, "#eta", "#phi");
1257 AliAnalysisTaskFragmentationFunction::SetProperties(fh1Pt, "p_{T} [GeV/c]", "entries");
1258 AliAnalysisTaskFragmentationFunction::SetProperties(fh2PhiPt, "#phi", "p_{T} [GeV/c]");
1259}
1260
1261//________________________________________________________________________________________________________
1262void AliAnalysisTaskFragmentationFunction::AliFragFuncQATrackHistos::FillTrackQA(Float_t eta, Float_t phi, Float_t pt, Bool_t weightPt, Float_t norm,
1263 Bool_t scaleStrangeness, Float_t scaleFacStrangeness)
1264{
1265 // fill track QA histos
1266 Float_t weight = 1.;
1267 if(weightPt) weight = pt;
1268 fh2EtaPhi->Fill( eta, phi, weight);
1269 if(scaleStrangeness) fh2EtaPhi->Fill( eta, phi, scaleFacStrangeness);
1270 if(pt > fHighPtThreshold) fh2HighPtEtaPhi->Fill( eta, phi, weight);
1271 if(pt > fHighPtThreshold && scaleStrangeness) fh2HighPtEtaPhi->Fill( eta, phi, weight);
1272 if(norm) fh1Pt->Fill( pt, 1/norm );
1273 else if(scaleStrangeness) fh1Pt->Fill(pt,scaleFacStrangeness);
1274 else fh1Pt->Fill( pt );
1275
1276 if(scaleFacStrangeness) fh2PhiPt->Fill(phi, pt, scaleFacStrangeness);
1277 else fh2PhiPt->Fill(phi, pt);
1278}
1279
1280//______________________________________________________________________________________
1281void AliAnalysisTaskFragmentationFunction::AliFragFuncQATrackHistos::AddToOutput(TList* list) const
1282{
1283 // add histos to list
1284
1285 list->Add(fh2EtaPhi);
1286 list->Add(fh2HighPtEtaPhi);
1287 list->Add(fh1Pt);
1288 list->Add(fh2PhiPt);
1289}
1290
1291//_________________________________________________________________________________
1292Bool_t AliAnalysisTaskFragmentationFunction::Notify()
1293{
1294 //
1295 // Implemented Notify() to read the cross sections
1296 // and number of trials from pyxsec.root
1297 // (taken from AliAnalysisTaskJetSpectrum2)
1298 //
1299 TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
1300 Float_t xsection = 0;
1301 Float_t ftrials = 1;
1302
1303 fAvgTrials = 1;
1304 if(tree){
1305 TFile *curfile = tree->GetCurrentFile();
1306 if (!curfile) {
1307 Error("Notify","No current file");
1308 return kFALSE;
1309 }
1310 if(!fh1Xsec||!fh1Trials){
1311 Printf("%s%d No Histogram fh1Xsec",(char*)__FILE__,__LINE__);
1312 return kFALSE;
1313 }
1314 AliAnalysisHelperJetTasks::PythiaInfoFromFile(curfile->GetName(),xsection,ftrials);
1315 fh1Xsec->Fill("<#sigma>",xsection);
1316 // construct a poor man average trials
1317 Float_t nEntries = (Float_t)tree->GetTree()->GetEntries();
1318 if(ftrials>=nEntries && nEntries>0.)fAvgTrials = ftrials/nEntries;
1319 }
1320
1321 // Set seed for backg study
1322 fRandom = new TRandom3();
1323 fRandom->SetSeed(0);
1324
1325 return kTRUE;
1326}
1327
1328//__________________________________________________________________
1329void AliAnalysisTaskFragmentationFunction::UserCreateOutputObjects()
1330{
1331 // create output objects
1332
1333 if(fDebug > 1) Printf("AliAnalysisTaskFragmentationFunction::UserCreateOutputObjects()");
1334
1335 // create list of tracks and jets
1336
1337 fTracksRecCuts = new TList();
1338 fTracksRecCuts->SetOwner(kFALSE);
1339
1340 fTracksGen = new TList();
1341 fTracksGen->SetOwner(kFALSE);
1342
1343 fTracksAODMCCharged = new TList();
1344 fTracksAODMCCharged->SetOwner(kFALSE);
1345
1346 fTracksAODMCChargedSecNS = new TList();
1347 fTracksAODMCChargedSecNS->SetOwner(kFALSE);
1348
1349 fTracksAODMCChargedSecS = new TList();
1350 fTracksAODMCChargedSecS->SetOwner(kFALSE);
1351
1352 fTracksRecQualityCuts = new TList();
1353 fTracksRecQualityCuts->SetOwner(kFALSE);
1354
1355 fJetsRec = new TList();
1356 fJetsRec->SetOwner(kFALSE);
1357
1358 fJetsRecCuts = new TList();
1359 fJetsRecCuts->SetOwner(kFALSE);
1360
1361 fJetsGen = new TList();
1362 fJetsGen->SetOwner(kFALSE);
1363
1364 fJetsRecEff = new TList();
1365 fJetsRecEff->SetOwner(kFALSE);
1366
1367 fJetsEmbedded = new TList();
1368 fJetsEmbedded->SetOwner(kFALSE);
1369
1370
1371 if(fBckgMode &&
1372 (fBckgType[0]==kBckgClusters || fBckgType[1]==kBckgClusters || fBckgType[2]==kBckgClusters || fBckgType[3]==kBckgClusters || fBckgType[4]==kBckgClusters ||
1373 fBckgType[0]==kBckgClustersOutLeading || fBckgType[1]==kBckgClustersOutLeading || fBckgType[2]==kBckgClustersOutLeading ||
1374 fBckgType[3]==kBckgClustersOutLeading || fBckgType[4]==kBckgClustersOutLeading)){
1375
1376 fBckgJetsRec = new TList();
1377 fBckgJetsRec->SetOwner(kFALSE);
1378
1379 fBckgJetsRecCuts = new TList();
1380 fBckgJetsRecCuts->SetOwner(kFALSE);
1381
1382 fBckgJetsGen = new TList();
1383 fBckgJetsGen->SetOwner(kFALSE);
1384 }
1385
1386 //
1387 // Create histograms / output container
1388 //
1389
1390 OpenFile(1);
1391 fCommonHistList = new TList();
1392 fCommonHistList->SetOwner(kTRUE);
1393
1394 Bool_t oldStatus = TH1::AddDirectoryStatus();
1395 TH1::AddDirectory(kFALSE);
1396
1397
1398 // Histograms
1399 fh1EvtSelection = new TH1F("fh1EvtSelection", "Event Selection", 6, -0.5, 5.5);
1400 fh1EvtSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
1401 fh1EvtSelection->GetXaxis()->SetBinLabel(2,"event selection: rejected");
1402 fh1EvtSelection->GetXaxis()->SetBinLabel(3,"event class: rejected");
1403 fh1EvtSelection->GetXaxis()->SetBinLabel(4,"vertex Ncontr: rejected");
1404 fh1EvtSelection->GetXaxis()->SetBinLabel(5,"vertex z: rejected");
1405 fh1EvtSelection->GetXaxis()->SetBinLabel(6,"vertex type: rejected");
1406
1407 fh1VertexNContributors = new TH1F("fh1VertexNContributors", "Vertex N contributors", 2500,-.5, 2499.5);
1408 fh1VertexZ = new TH1F("fh1VertexZ", "Vertex z distribution", 30, -15., 15.);
1409 fh1EvtMult = new TH1F("fh1EvtMult","Event multiplicity, track pT cut > 150 MeV/c, |#eta| < 0.9",120,0.,12000.);
1410 fh1EvtCent = new TH1F("fh1EvtCent","centrality",100,0.,100.);
1411
1412 fh1Xsec = new TProfile("fh1Xsec","xsec from pyxsec.root",1,0,1);
1413 fh1Xsec->GetXaxis()->SetBinLabel(1,"<#sigma>");
1414 fh1Trials = new TH1F("fh1Trials","trials from pyxsec.root",1,0,1);
1415 fh1Trials->GetXaxis()->SetBinLabel(1,"#sum{ntrials}");
1416 fh1PtHard = new TH1F("fh1PtHard","PYTHIA Pt hard;p_{T,hard}",350,-.5,349.5);
1417 fh1PtHardTrials = new TH1F("fh1PtHardTrials","PYTHIA Pt hard weight with trials;p_{T,hard}",350,-.5,349.5);
1418
1419 fh1nRecJetsCuts = new TH1F("fh1nRecJetsCuts","reconstructed jets per event",10,-0.5,9.5);
1420 fh1nGenJets = new TH1F("fh1nGenJets","generated jets per event",10,-0.5,9.5);
1421 fh1nRecEffJets = new TH1F("fh1nRecEffJets","reconstruction effiency: jets per event",10,-0.5,9.5);
1422 fh1nEmbeddedJets = new TH1F("fh1nEmbeddedJets","embedded jets per event",10,-0.5,9.5);
1423
1424 fh2PtRecVsGenPrim = new TH2F("fh2PtRecVsGenPrim","rec vs gen pt",fQATrackNBinsPt,fQATrackPtMin,fQATrackPtMax,fQATrackNBinsPt,fQATrackPtMin,fQATrackPtMax);
1425 fh2PtRecVsGenSec = new TH2F("fh2PtRecVsGenSec","rec vs gen pt",fQATrackNBinsPt,fQATrackPtMin,fQATrackPtMax,fQATrackNBinsPt,fQATrackPtMin,fQATrackPtMax);
1426
1427 // embedding
1428 if(fBranchEmbeddedJets.Length()){
1429 fh1FractionPtEmbedded = new TH1F("fh1FractionPtEmbedded","",200,0,2);
1430 fh1IndexEmbedded = new TH1F("fh1IndexEmbedded","",11,-1,10);
1431 fh2DeltaPtVsJetPtEmbedded = new TH2F("fh2DeltaPtVsJetPtEmbedded","",250,0,250,200,-100,100);
1432 fh2DeltaPtVsRecJetPtEmbedded = new TH2F("fh2DeltaPtVsRecJetPtEmbedded","",250,0,250,200,-100,100);
1433 fh1DeltaREmbedded = new TH1F("fh1DeltaREmbedded","",50,0,0.5);
1434 fh1nEmbeddedJets = new TH1F("fh1nEmbeddedJets","embedded jets per event",10,-0.5,9.5);
1435 }
1436
1437
1438 if(fQAMode){
1439 if(fQAMode&1){ // track QA
1440 fQATrackHistosRecCuts = new AliFragFuncQATrackHistos("RecCuts", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1441 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1442 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1443 fQATrackHighPtThreshold);
1444 fQATrackHistosGen = new AliFragFuncQATrackHistos("Gen", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1445 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1446 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1447 fQATrackHighPtThreshold);
1448 }
1449
1450 if(fQAMode&2){ // jet QA
1451 fQAJetHistosRec = new AliFragFuncQAJetHistos("Rec", fQAJetNBinsPt, fQAJetPtMin, fQAJetPtMax,
1452 fQAJetNBinsEta, fQAJetEtaMin, fQAJetEtaMax,
1453 fQAJetNBinsPhi, fQAJetPhiMin, fQAJetPhiMax);
1454 fQAJetHistosRecCuts = new AliFragFuncQAJetHistos("RecCuts", fQAJetNBinsPt, fQAJetPtMin, fQAJetPtMax,
1455 fQAJetNBinsEta, fQAJetEtaMin, fQAJetEtaMax,
1456 fQAJetNBinsPhi, fQAJetPhiMin, fQAJetPhiMax);
1457 fQAJetHistosRecCutsLeading = new AliFragFuncQAJetHistos("RecCutsLeading", fQAJetNBinsPt, fQAJetPtMin, fQAJetPtMax,
1458 fQAJetNBinsEta, fQAJetEtaMin, fQAJetEtaMax,
1459 fQAJetNBinsPhi, fQAJetPhiMin, fQAJetPhiMax);
1460 fQAJetHistosGen = new AliFragFuncQAJetHistos("Gen", fQAJetNBinsPt, fQAJetPtMin, fQAJetPtMax,
1461 fQAJetNBinsEta, fQAJetEtaMin, fQAJetEtaMax,
1462 fQAJetNBinsPhi, fQAJetPhiMin, fQAJetPhiMax);
1463 fQAJetHistosGenLeading = new AliFragFuncQAJetHistos("GenLeading", fQAJetNBinsPt, fQAJetPtMin, fQAJetPtMax,
1464 fQAJetNBinsEta, fQAJetEtaMin, fQAJetEtaMax,
1465 fQAJetNBinsPhi, fQAJetPhiMin, fQAJetPhiMax);
1466 if(fEffMode) fQAJetHistosRecEffLeading = new AliFragFuncQAJetHistos("RecEffLeading", fQAJetNBinsPt, fQAJetPtMin, fQAJetPtMax,
1467 fQAJetNBinsEta, fQAJetEtaMin, fQAJetEtaMax,fQAJetNBinsPhi, fQAJetPhiMin, fQAJetPhiMax);
1468 }
1469 } // end: QA
1470
1471 if(fFFMode){
1472
1473 fFFHistosRecCuts = new AliFragFuncHistos("RecCuts", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1474 fFFNBinsPt, fFFPtMin, fFFPtMax,
1475 fFFNBinsXi, fFFXiMin, fFFXiMax,
1476 fFFNBinsZ , fFFZMin , fFFZMax );
1477
37f35567 1478
1479 fFFHistosRecCutsInc = new AliFragFuncHistos("RecCutsInc", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
ffa175ff 1480 fFFNBinsPt, fFFPtMin, fFFPtMax,
1481 fFFNBinsXi, fFFXiMin, fFFXiMax,
37f35567 1482 fFFNBinsZ , fFFZMin , fFFZMax );
1483
1484
1485 fFFHistosRecLeadingTrack = new AliFragFuncHistos("RecLeadingTrack", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1486 fFFNBinsPt, fFFPtMin, fFFPtMax,
1487 fFFNBinsXi, fFFXiMin, fFFXiMax,
1488 fFFNBinsZ , fFFZMin , fFFZMax );
1489
1490 fFFHistosGen = new AliFragFuncHistos("Gen", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1491 fFFNBinsPt, fFFPtMin, fFFPtMax,
1492 fFFNBinsXi, fFFXiMin, fFFXiMax,
1493 fFFNBinsZ , fFFZMin , fFFZMax);
1494
1495 fFFHistosGenInc = new AliFragFuncHistos("GenInc", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1496 fFFNBinsPt, fFFPtMin, fFFPtMax,
1497 fFFNBinsXi, fFFXiMin, fFFXiMax,
1498 fFFNBinsZ , fFFZMin , fFFZMax);
1499
1500 fFFHistosGenLeadingTrack = new AliFragFuncHistos("GenLeadingTrack", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1501 fFFNBinsPt, fFFPtMin, fFFPtMax,
1502 fFFNBinsXi, fFFXiMin, fFFXiMax,
1503 fFFNBinsZ , fFFZMin , fFFZMax);
1504
ffa175ff 1505 } // end: FF
1506
1507 // efficiency
1508
1509 if(fEffMode){
1510 if(fQAMode&1){
1511 fQATrackHistosRecEffGen = new AliFragFuncQATrackHistos("RecEffGen", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1512 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1513 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1514 fQATrackHighPtThreshold);
1515
1516 fQATrackHistosRecEffRec = new AliFragFuncQATrackHistos("RecEffRec", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1517 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1518 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1519 fQATrackHighPtThreshold);
1520
1521 fQATrackHistosSecRecNS = new AliFragFuncQATrackHistos("SecRecNS", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1522 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1523 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1524 fQATrackHighPtThreshold);
1525
1526 fQATrackHistosSecRecS = new AliFragFuncQATrackHistos("SecRecS", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1527 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1528 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1529 fQATrackHighPtThreshold);
1530
1531 fQATrackHistosSecRecSsc = new AliFragFuncQATrackHistos("SecRecSsc", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1532 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1533 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1534 fQATrackHighPtThreshold);
1535
1536 }
1537 if(fFFMode){
1538 fFFHistosRecEffRec = new AliFragFuncHistos("RecEffRec", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1539 fFFNBinsPt, fFFPtMin, fFFPtMax,
1540 fFFNBinsXi, fFFXiMin, fFFXiMax,
1541 fFFNBinsZ , fFFZMin , fFFZMax);
1542
1543 fFFHistosSecRecNS = new AliFragFuncHistos("SecRecNS", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1544 fFFNBinsPt, fFFPtMin, fFFPtMax,
1545 fFFNBinsXi, fFFXiMin, fFFXiMax,
1546 fFFNBinsZ , fFFZMin , fFFZMax);
1547
1548 fFFHistosSecRecS = new AliFragFuncHistos("SecRecS", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1549 fFFNBinsPt, fFFPtMin, fFFPtMax,
1550 fFFNBinsXi, fFFXiMin, fFFXiMax,
1551 fFFNBinsZ , fFFZMin , fFFZMax);
1552
1553 fFFHistosSecRecSsc = new AliFragFuncHistos("SecRecSsc", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1554 fFFNBinsPt, fFFPtMin, fFFPtMax,
1555 fFFNBinsXi, fFFXiMin, fFFXiMax,
1556 fFFNBinsZ , fFFZMin , fFFZMax);
1557
1558 }
1559 } // end: efficiency
1560
1561 // Background
1562 if(fBckgMode){
1563 if(fBckgType[0]==kBckgNone){
1564 AliError("no bgr method selected !");
1565 }
1566
1567 TString title[5];
1568 for(Int_t i=0; i<5; i++){
1569 if(fBckgType[i]==kBckgPerp) title[i]="Perp";
1570 else if(fBckgType[i]==kBckgPerp2) title[i]="Perp2";
1571 else if(fBckgType[i]==kBckgPerp2Area) title[i]="Perp2Area";
1572 else if(fBckgType[i]==kBckgPerpWindow) title[i]="PerpW";
1573 else if(fBckgType[i]==kBckgASide) title[i]="ASide";
1574 else if(fBckgType[i]==kBckgASideWindow) title[i]="ASideW";
1575 else if(fBckgType[i]==kBckgOutLJ) title[i]="OutLeadingJet";
1576 else if(fBckgType[i]==kBckgOut2J) title[i]="Out2Jets";
1577 else if(fBckgType[i]==kBckgOut3J) title[i]="Out3Jets";
1578 else if(fBckgType[i]==kBckgOutAJ) title[i]="AllJets";
1579 else if(fBckgType[i]==kBckgOutLJStat) title[i]="OutLeadingJetStat";
1580 else if(fBckgType[i]==kBckgOut2JStat) title[i]="Out2JetsStat";
1581 else if(fBckgType[i]==kBckgOut3JStat) title[i]="Out3JetsStat";
1582 else if(fBckgType[i]==kBckgOutAJStat) title[i]="AllJetsStat";
1583 else if(fBckgType[i]==kBckgClustersOutLeading) title[i]="OutClusters";
1584 else if(fBckgType[i]==kBckgClusters) title[i]="MedianClusters";
1585 else if(fBckgType[i]==kBckgNone) title[i]="";
1586 else printf("Please chose background method number %d!",i);
1587 }
1588
1589
1590 if(fBckgType[0]==kBckgClusters || fBckgType[1]==kBckgClusters || fBckgType[2]==kBckgClusters || fBckgType[3]==kBckgClusters || fBckgType[4]==kBckgClusters ||
1591 fBckgType[0]==kBckgClustersOutLeading || fBckgType[1]==kBckgClustersOutLeading || fBckgType[2]==kBckgClustersOutLeading ||
1592 fBckgType[3]==kBckgClustersOutLeading || fBckgType[4]==kBckgClustersOutLeading){
1593
1594 fh1nRecBckgJetsCuts = new TH1F("fh1nRecBckgJetsCuts","reconstructed background jets per event",10,-0.5,9.5);
1595 fh1nGenBckgJets = new TH1F("fh1nGenBckgJets","generated background jets per event",10,-0.5,9.5);
1596 }
1597
1598
1599 fh1BckgMult0 = new TH1F("fh1BckgMult0","bckg mult "+title[0],500,0,500);
1600 if(fBckgType[1] != kBckgNone) fh1BckgMult1 = new TH1F("fh1BckgMult1","bckg mult "+title[1],500,0,500);
1601 if(fBckgType[2] != kBckgNone) fh1BckgMult2 = new TH1F("fh1BckgMult2","bckg mult "+title[2],500,0,500);
1602 if(fBckgType[3] != kBckgNone) fh1BckgMult3 = new TH1F("fh1BckgMult3","bckg mult "+title[3],500,0,500);
1603 if(fBckgType[4] != kBckgNone) fh1BckgMult4 = new TH1F("fh1BckgMult4","bckg mult "+title[4],500,0,500);
1604
1605
1606 if(fQAMode&1){
1607 fQABckgHisto0RecCuts = new AliFragFuncQATrackHistos("Bckg"+title[0]+"RecCuts", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1608 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1609 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1610 fQATrackHighPtThreshold);
1611 fQABckgHisto0Gen = new AliFragFuncQATrackHistos("Bckg"+title[0]+"Gen", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1612 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1613 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1614 fQATrackHighPtThreshold);
1615
1616 if(fBckgType[1] != kBckgNone){
1617 fQABckgHisto1RecCuts = new AliFragFuncQATrackHistos("Bckg"+title[1]+"RecCuts", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1618 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1619 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1620 fQATrackHighPtThreshold);
1621 fQABckgHisto1Gen = new AliFragFuncQATrackHistos("Bckg"+title[1]+"Gen", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1622 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1623 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1624 fQATrackHighPtThreshold);
1625 }
1626 if(fBckgType[2] != kBckgNone){
1627 fQABckgHisto2RecCuts = new AliFragFuncQATrackHistos("Bckg"+title[2]+"RecCuts", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1628 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1629 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1630 fQATrackHighPtThreshold);
1631 fQABckgHisto2Gen = new AliFragFuncQATrackHistos("Bckg"+title[2]+"Gen", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1632 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1633 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1634 fQATrackHighPtThreshold);
1635 }
1636 if(fBckgType[3] != kBckgNone){
1637 fQABckgHisto3RecCuts = new AliFragFuncQATrackHistos("Bckg"+title[3]+"RecCuts", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1638 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1639 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1640 fQATrackHighPtThreshold);
1641 fQABckgHisto3Gen = new AliFragFuncQATrackHistos("Bckg"+title[3]+"Gen", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1642 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1643 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1644 fQATrackHighPtThreshold);
1645 }
1646 if(fBckgType[4] != kBckgNone){
1647 fQABckgHisto4RecCuts = new AliFragFuncQATrackHistos("Bckg"+title[4]+"RecCuts", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1648 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1649 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1650 fQATrackHighPtThreshold);
1651 fQABckgHisto4Gen = new AliFragFuncQATrackHistos("Bckg"+title[4]+"Gen", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1652 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1653 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1654 fQATrackHighPtThreshold);
1655 }
1656 } // end: background QA
1657
1658 if(fFFMode){
1659 fFFBckgHisto0RecCuts = new AliFragFuncHistos("Bckg"+title[0]+"RecCuts", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1660 fFFNBinsPt, fFFPtMin, fFFPtMax,
1661 fFFNBinsXi, fFFXiMin, fFFXiMax,
1662 fFFNBinsZ , fFFZMin , fFFZMax);
1663
1664 fFFBckgHisto0Gen = new AliFragFuncHistos("Bckg"+title[0]+"Gen", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1665 fFFNBinsPt, fFFPtMin, fFFPtMax,
1666 fFFNBinsXi, fFFXiMin, fFFXiMax,
1667 fFFNBinsZ , fFFZMin , fFFZMax);
1668
1669 if(fBckgType[1] != kBckgNone){
1670 fFFBckgHisto1RecCuts = new AliFragFuncHistos("Bckg"+title[1]+"RecCuts", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1671 fFFNBinsPt, fFFPtMin, fFFPtMax,
1672 fFFNBinsXi, fFFXiMin, fFFXiMax,
1673 fFFNBinsZ , fFFZMin , fFFZMax);
1674 fFFBckgHisto1Gen = new AliFragFuncHistos("Bckg"+title[1]+"Gen", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1675 fFFNBinsPt, fFFPtMin, fFFPtMax,
1676 fFFNBinsXi, fFFXiMin, fFFXiMax,
1677 fFFNBinsZ , fFFZMin , fFFZMax);
1678 }
1679 if(fBckgType[2] != kBckgNone){
1680 fFFBckgHisto2RecCuts = new AliFragFuncHistos("Bckg"+title[2]+"RecCuts", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1681 fFFNBinsPt, fFFPtMin, fFFPtMax,
1682 fFFNBinsXi, fFFXiMin, fFFXiMax,
1683 fFFNBinsZ , fFFZMin , fFFZMax);
1684
1685 fFFBckgHisto2Gen = new AliFragFuncHistos("Bckg"+title[2]+"Gen", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1686 fFFNBinsPt, fFFPtMin, fFFPtMax,
1687 fFFNBinsXi, fFFXiMin, fFFXiMax,
1688 fFFNBinsZ , fFFZMin , fFFZMax);
1689 }
1690 if(fBckgType[3] != kBckgNone){
1691 fFFBckgHisto3RecCuts = new AliFragFuncHistos("Bckg"+title[3]+"RecCuts", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1692 fFFNBinsPt, fFFPtMin, fFFPtMax,
1693 fFFNBinsXi, fFFXiMin, fFFXiMax,
1694 fFFNBinsZ , fFFZMin , fFFZMax);
1695
1696 fFFBckgHisto3Gen = new AliFragFuncHistos("Bckg"+title[3]+"Gen", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1697 fFFNBinsPt, fFFPtMin, fFFPtMax,
1698 fFFNBinsXi, fFFXiMin, fFFXiMax,
1699 fFFNBinsZ , fFFZMin , fFFZMax);
1700 }
1701 if(fBckgType[4] != kBckgNone){
1702 fFFBckgHisto4RecCuts = new AliFragFuncHistos("Bckg"+title[4]+"RecCuts", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1703 fFFNBinsPt, fFFPtMin, fFFPtMax,
1704 fFFNBinsXi, fFFXiMin, fFFXiMax,
1705 fFFNBinsZ , fFFZMin , fFFZMax);
1706
1707 fFFBckgHisto4Gen = new AliFragFuncHistos("Bckg"+title[4]+"Gen", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1708 fFFNBinsPt, fFFPtMin, fFFPtMax,
1709 fFFNBinsXi, fFFXiMin, fFFXiMax,
1710 fFFNBinsZ , fFFZMin , fFFZMax);
1711 }
1712 if(fEffMode){
1713 fFFBckgHisto0RecEffRec = new AliFragFuncHistos("Bckg"+title[0]+"RecEffRec", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1714 fFFNBinsPt, fFFPtMin, fFFPtMax,
1715 fFFNBinsXi, fFFXiMin, fFFXiMax,
1716 fFFNBinsZ , fFFZMin , fFFZMax);
1717
1718 fFFBckgHisto0SecRecNS = new AliFragFuncHistos("Bckg"+title[0]+"SecRecNS", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1719 fFFNBinsPt, fFFPtMin, fFFPtMax,
1720 fFFNBinsXi, fFFXiMin, fFFXiMax,
1721 fFFNBinsZ , fFFZMin , fFFZMax);
1722
1723 fFFBckgHisto0SecRecS = new AliFragFuncHistos("Bckg"+title[0]+"SecRecS", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1724 fFFNBinsPt, fFFPtMin, fFFPtMax,
1725 fFFNBinsXi, fFFXiMin, fFFXiMax,
1726 fFFNBinsZ , fFFZMin , fFFZMax);
1727
1728 fFFBckgHisto0SecRecSsc = new AliFragFuncHistos("Bckg"+title[0]+"SecRecSsc", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1729 fFFNBinsPt, fFFPtMin, fFFPtMax,
1730 fFFNBinsXi, fFFXiMin, fFFXiMax,
1731 fFFNBinsZ , fFFZMin , fFFZMax);
1732
1733 }
1734 } // end: background FF
1735
1736
1737 } // end: background
1738
1739
1740 // ____________ define histograms ____________________
1741
1742 if(fQAMode){
1743 if(fQAMode&1){ // track QA
1744 fQATrackHistosRecCuts->DefineHistos();
1745 fQATrackHistosGen->DefineHistos();
1746 }
1747
1748 if(fQAMode&2){ // jet QA
1749 fQAJetHistosRec->DefineHistos();
1750 fQAJetHistosRecCuts->DefineHistos();
1751 fQAJetHistosRecCutsLeading->DefineHistos();
1752 fQAJetHistosGen->DefineHistos();
1753 fQAJetHistosGenLeading->DefineHistos();
1754 if(fEffMode) fQAJetHistosRecEffLeading->DefineHistos();
1755 }
1756 }
1757
1758 if(fFFMode){
1759 fFFHistosRecCuts->DefineHistos();
37f35567 1760 fFFHistosRecCutsInc->DefineHistos();
1761 fFFHistosRecLeadingTrack->DefineHistos();
ffa175ff 1762 fFFHistosGen->DefineHistos();
37f35567 1763 fFFHistosGenInc->DefineHistos();
1764 fFFHistosGenLeadingTrack->DefineHistos();
ffa175ff 1765 }
1766
1767 if(fEffMode){
1768 if(fQAMode&1){
1769 fQATrackHistosRecEffGen->DefineHistos();
1770 fQATrackHistosRecEffRec->DefineHistos();
1771 fQATrackHistosSecRecNS->DefineHistos();
1772 fQATrackHistosSecRecS->DefineHistos();
1773 fQATrackHistosSecRecSsc->DefineHistos();
1774 }
1775 if(fFFMode){
1776 fFFHistosRecEffRec->DefineHistos();
1777 fFFHistosSecRecNS->DefineHistos();
1778 fFFHistosSecRecS->DefineHistos();
1779 fFFHistosSecRecSsc->DefineHistos();
1780 }
1781 } // end: efficiency
1782
1783 // Background
1784 if(fBckgMode){
1785 if(fFFMode){
1786 fFFBckgHisto0RecCuts->DefineHistos();
1787 fFFBckgHisto0Gen->DefineHistos();
1788 if(fBckgType[1] != kBckgNone) fFFBckgHisto1RecCuts->DefineHistos();
1789 if(fBckgType[1] != kBckgNone) fFFBckgHisto1Gen->DefineHistos();
1790 if(fBckgType[2] != kBckgNone) fFFBckgHisto2RecCuts->DefineHistos();
1791 if(fBckgType[2] != kBckgNone) fFFBckgHisto2Gen->DefineHistos();
1792 if(fBckgType[3] != kBckgNone) fFFBckgHisto3RecCuts->DefineHistos();
1793 if(fBckgType[3] != kBckgNone) fFFBckgHisto3Gen->DefineHistos();
1794 if(fBckgType[4] != kBckgNone) fFFBckgHisto4RecCuts->DefineHistos();
1795 if(fBckgType[4] != kBckgNone) fFFBckgHisto4Gen->DefineHistos();
1796
1797 if(fEffMode){
1798 fFFBckgHisto0RecEffRec->DefineHistos();
1799 fFFBckgHisto0SecRecNS->DefineHistos();
1800 fFFBckgHisto0SecRecS->DefineHistos();
1801 fFFBckgHisto0SecRecSsc->DefineHistos();
1802 }
1803 }
1804
1805 if(fQAMode&1){
1806 fQABckgHisto0RecCuts->DefineHistos();
1807 fQABckgHisto0Gen->DefineHistos();
1808 if(fBckgType[1] != kBckgNone) fQABckgHisto1RecCuts->DefineHistos();
1809 if(fBckgType[1] != kBckgNone) fQABckgHisto1Gen->DefineHistos();
1810 if(fBckgType[2] != kBckgNone) fQABckgHisto2RecCuts->DefineHistos();
1811 if(fBckgType[2] != kBckgNone) fQABckgHisto2Gen->DefineHistos();
1812 if(fBckgType[3] != kBckgNone) fQABckgHisto3RecCuts->DefineHistos();
1813 if(fBckgType[3] != kBckgNone) fQABckgHisto3Gen->DefineHistos();
1814 if(fBckgType[4] != kBckgNone) fQABckgHisto4RecCuts->DefineHistos();
1815 if(fBckgType[4] != kBckgNone) fQABckgHisto4Gen->DefineHistos();
1816 }
1817 } // end: background
1818
1819
1820 Bool_t genJets = (fJetTypeGen != kJetsUndef) ? kTRUE : kFALSE;
1821 Bool_t genTracks = (fTrackTypeGen != kTrackUndef) ? kTRUE : kFALSE;
1822 Bool_t recJetsEff = (fJetTypeRecEff != kJetsUndef) ? kTRUE : kFALSE;
1823
1824 fCommonHistList->Add(fh1EvtSelection);
1825 fCommonHistList->Add(fh1EvtMult);
1826 fCommonHistList->Add(fh1EvtCent);
1827 fCommonHistList->Add(fh1VertexNContributors);
1828 fCommonHistList->Add(fh1VertexZ);
1829 fCommonHistList->Add(fh1nRecJetsCuts);
1830 fCommonHistList->Add(fh1Xsec);
1831 fCommonHistList->Add(fh1Trials);
1832 fCommonHistList->Add(fh1PtHard);
1833 fCommonHistList->Add(fh1PtHardTrials);
1834
1835 if(genJets) fCommonHistList->Add(fh1nGenJets);
1836
1837 // FF histograms
1838 if(fFFMode){
1839 fFFHistosRecCuts->AddToOutput(fCommonHistList);
37f35567 1840 fFFHistosRecCutsInc->AddToOutput(fCommonHistList);
1841 fFFHistosRecLeadingTrack->AddToOutput(fCommonHistList);
1842
ffa175ff 1843 if(genJets && genTracks){
1844 fFFHistosGen->AddToOutput(fCommonHistList);
37f35567 1845 fFFHistosGenInc->AddToOutput(fCommonHistList);
1846 fFFHistosGenLeadingTrack->AddToOutput(fCommonHistList);
ffa175ff 1847 }
1848 }
1849
1850 // Background
1851 if(fBckgMode){
1852 if(fFFMode){
1853 fFFBckgHisto0RecCuts->AddToOutput(fCommonHistList);
1854 if(fBckgType[1] != kBckgNone) fFFBckgHisto1RecCuts->AddToOutput(fCommonHistList);
1855 if(fBckgType[2] != kBckgNone) fFFBckgHisto2RecCuts->AddToOutput(fCommonHistList);
1856 if(fBckgType[3] != kBckgNone) fFFBckgHisto3RecCuts->AddToOutput(fCommonHistList);
1857 if(fBckgType[4] != kBckgNone) fFFBckgHisto4RecCuts->AddToOutput(fCommonHistList);
1858
1859 if(genJets && genTracks){
1860 fFFBckgHisto0Gen->AddToOutput(fCommonHistList);
1861 if(fBckgType[1] != kBckgNone) fFFBckgHisto1Gen->AddToOutput(fCommonHistList);
1862 if(fBckgType[2] != kBckgNone) fFFBckgHisto2Gen->AddToOutput(fCommonHistList);
1863 if(fBckgType[3] != kBckgNone) fFFBckgHisto3Gen->AddToOutput(fCommonHistList);
1864 if(fBckgType[4] != kBckgNone) fFFBckgHisto4Gen->AddToOutput(fCommonHistList);
1865 }
1866
1867 if(fEffMode){
1868 fFFBckgHisto0RecEffRec->AddToOutput(fCommonHistList);
1869 fFFBckgHisto0SecRecNS->AddToOutput(fCommonHistList);
1870 fFFBckgHisto0SecRecS->AddToOutput(fCommonHistList);
1871 fFFBckgHisto0SecRecSsc->AddToOutput(fCommonHistList);
1872 }
1873 }
1874
1875 if(fQAMode&1){
1876 fQABckgHisto0RecCuts->AddToOutput(fCommonHistList);
1877 if(fBckgType[1] != kBckgNone) fQABckgHisto1RecCuts->AddToOutput(fCommonHistList);
1878 if(fBckgType[2] != kBckgNone) fQABckgHisto2RecCuts->AddToOutput(fCommonHistList);
1879 if(fBckgType[3] != kBckgNone) fQABckgHisto3RecCuts->AddToOutput(fCommonHistList);
1880 if(fBckgType[4] != kBckgNone) fQABckgHisto4RecCuts->AddToOutput(fCommonHistList);
1881 if(genJets && genTracks){
1882 fQABckgHisto0Gen->AddToOutput(fCommonHistList);
1883 if(fBckgType[1] != kBckgNone) fQABckgHisto1Gen->AddToOutput(fCommonHistList);
1884 if(fBckgType[2] != kBckgNone) fQABckgHisto2Gen->AddToOutput(fCommonHistList);
1885 if(fBckgType[3] != kBckgNone) fQABckgHisto3Gen->AddToOutput(fCommonHistList);
1886 if(fBckgType[4] != kBckgNone) fQABckgHisto4Gen->AddToOutput(fCommonHistList);
1887 }
1888 }
1889
1890 if(fh1BckgMult0) fCommonHistList->Add(fh1BckgMult0);
1891 if(fBckgType[1] != kBckgNone) fCommonHistList->Add(fh1BckgMult1);
1892 if(fBckgType[2] != kBckgNone) fCommonHistList->Add(fh1BckgMult2);
1893 if(fBckgType[3] != kBckgNone) fCommonHistList->Add(fh1BckgMult3);
1894 if(fBckgType[4] != kBckgNone) fCommonHistList->Add(fh1BckgMult4);
1895 }
1896
1897
1898 if(fBranchEmbeddedJets.Length()){
1899 fCommonHistList->Add(fh1FractionPtEmbedded);
1900 fCommonHistList->Add(fh1IndexEmbedded);
1901 fCommonHistList->Add(fh2DeltaPtVsJetPtEmbedded);
1902 fCommonHistList->Add(fh2DeltaPtVsRecJetPtEmbedded);
1903 fCommonHistList->Add(fh1DeltaREmbedded);
1904 fCommonHistList->Add(fh1nEmbeddedJets);
1905 }
1906
1907
1908 // QA
1909 if(fQAMode){
1910 if(fQAMode&1){ // track QA
1911 fQATrackHistosRecCuts->AddToOutput(fCommonHistList);
1912 if(genTracks) fQATrackHistosGen->AddToOutput(fCommonHistList);
1913 }
1914
1915 if(fQAMode&2){ // jet QA
1916 fQAJetHistosRec->AddToOutput(fCommonHistList);
1917 fQAJetHistosRecCuts->AddToOutput(fCommonHistList);
1918 fQAJetHistosRecCutsLeading->AddToOutput(fCommonHistList);
1919 if(recJetsEff && fEffMode) fQAJetHistosRecEffLeading->AddToOutput(fCommonHistList);
1920 if(genJets){
1921 fQAJetHistosGen->AddToOutput(fCommonHistList);
1922 fQAJetHistosGenLeading->AddToOutput(fCommonHistList);
1923 }
1924 }
1925 }
1926
1927 if(fBckgMode &&
1928 (fBckgType[0]==kBckgClusters || fBckgType[1]==kBckgClusters || fBckgType[2]==kBckgClusters || fBckgType[3]==kBckgClusters || fBckgType[4]==kBckgClusters ||
1929 fBckgType[0]==kBckgClustersOutLeading || fBckgType[1]==kBckgClustersOutLeading || fBckgType[2]==kBckgClustersOutLeading ||
1930 fBckgType[3]==kBckgClustersOutLeading || fBckgType[4]==kBckgClustersOutLeading)) {
1931 fCommonHistList->Add(fh1nRecBckgJetsCuts);
1932 if(genJets) fCommonHistList->Add(fh1nGenBckgJets);
1933 }
1934
1935
1936 if(fEffMode && recJetsEff && genTracks){
1937 if(fQAMode&1){
1938 fQATrackHistosRecEffGen->AddToOutput(fCommonHistList);
1939 fQATrackHistosRecEffRec->AddToOutput(fCommonHistList);
1940 fQATrackHistosSecRecNS->AddToOutput(fCommonHistList);
1941 fQATrackHistosSecRecS->AddToOutput(fCommonHistList);
1942 fQATrackHistosSecRecSsc->AddToOutput(fCommonHistList);
1943 }
1944 if(fFFMode){
1945 fFFHistosRecEffRec->AddToOutput(fCommonHistList);
1946 fFFHistosSecRecNS->AddToOutput(fCommonHistList);
1947 fFFHistosSecRecS->AddToOutput(fCommonHistList);
1948 fFFHistosSecRecSsc->AddToOutput(fCommonHistList);
1949 }
1950 fCommonHistList->Add(fh1nRecEffJets);
1951 fCommonHistList->Add(fh2PtRecVsGenPrim);
1952 fCommonHistList->Add(fh2PtRecVsGenSec);
1953 }
1954
1955 // jet shape
1956 if(fJSMode){
1957
1958 fProNtracksLeadingJet = new TProfile("AvgNoOfTracksLeadingJet","AvgNoOfTracksLeadingJet",100,0,250,0,50);
1959 fProDelR80pcPt = new TProfile("AvgdelR80pcPt","AvgdelR80pcPt",100,0,250,0,1);
1960
1961 if(genJets && genTracks){
1962 fProNtracksLeadingJetGen = new TProfile("AvgNoOfTracksLeadingJetGen","AvgNoOfTracksLeadingJetGen",100,0,250,0,50);
1963 fProDelR80pcPtGen = new TProfile("AvgdelR80pcPtGen","AvgdelR80pcPtGen",100,0,250,0,1);
1964 }
1965
1966 if(fBckgMode)
1967 fProNtracksLeadingJetBgrPerp2 = new TProfile("AvgNoOfTracksLeadingJetBgrPerp2","AvgNoOfTracksLeadingJetBgrPerp2",100,0,250,0,50);
1968
1969 if(fEffMode){
1970 fProNtracksLeadingJetRecPrim = new TProfile("AvgNoOfTracksLeadingJetRecPrim","AvgNoOfTracksLeadingJetRecPrim",100,0,250,0,50);
1971 fProDelR80pcPtRecPrim = new TProfile("AvgdelR80pcPtRecPrim","AvgdelR80pcPtRecPrim",100,0,250,0,1);
1972 fProNtracksLeadingJetRecSecNS = new TProfile("AvgNoOfTracksLeadingJetRecSecNS","AvgNoOfTracksLeadingJetRecSecNS",100,0,250,0,50);
1973 fProNtracksLeadingJetRecSecS = new TProfile("AvgNoOfTracksLeadingJetRecSecS","AvgNoOfTracksLeadingJetRecSecS",100,0,250,0,50);
1974 fProNtracksLeadingJetRecSecSsc = new TProfile("AvgNoOfTracksLeadingJetRecSecSsc","AvgNoOfTracksLeadingJetRecSecSsc",100,0,250,0,50);
1975 }
1976
1977 TString strTitJS;
1978 for(Int_t ii=0; ii<5; ii++){
1979 if(ii==0)strTitJS = "_JetPt20to30";
1980 if(ii==1)strTitJS = "_JetPt30to40";
1981 if(ii==2)strTitJS = "_JetPt40to60";
1982 if(ii==3)strTitJS = "_JetPt60to80";
1983 if(ii==4)strTitJS = "_JetPt80to100";
1984
fd5c68ba 1985 fProDelRPtSum[ii] = new TProfile(Form("AvgPtSumDelR%s",strTitJS.Data()),Form("AvgPtSumDelR%s",strTitJS.Data()),50,0,1,0,250);
ffa175ff 1986 if(genJets && genTracks)
fd5c68ba 1987 fProDelRPtSumGen[ii] = new TProfile(Form("AvgPtSumDelRGen%s",strTitJS.Data()),Form("AvgPtSumDelRGen%s",strTitJS.Data()),50,0,1,0,250);
ffa175ff 1988 if(fBckgMode)
fd5c68ba 1989 fProDelRPtSumBgrPerp2[ii] = new TProfile(Form("AvgPtSumDelRBgrPerp2%s",strTitJS.Data()),Form("AvgPtSumDelRBgrPerp2%s",strTitJS.Data()),50,0,1,0,250);
ffa175ff 1990 if(fEffMode){
fd5c68ba 1991 fProDelRPtSumRecPrim[ii] = new TProfile(Form("AvgPtSumDelRRecPrim%s",strTitJS.Data()),Form("AvgPtSumDelRRecPrim%s",strTitJS.Data()),50,0,1,0,250);
1992 fProDelRPtSumRecSecNS[ii] = new TProfile(Form("AvgPtSumDelRRecSecNS%s",strTitJS.Data()),Form("AvgPtSumDelRRecSecNS%s",strTitJS.Data()),50,0,1,0,250);
1993 fProDelRPtSumRecSecS[ii] = new TProfile(Form("AvgPtSumDelRRecSecS%s",strTitJS.Data()),Form("AvgPtSumDelRRecSecS%s",strTitJS.Data()),50,0,1,0,250);
1994 fProDelRPtSumRecSecSsc[ii] = new TProfile(Form("AvgPtSumDelRRecSecSsc%s",strTitJS.Data()),Form("AvgPtSumDelRRecSecSsc%s",strTitJS.Data()),50,0,1,0,250);
ffa175ff 1995 }
1996 }
1997
1998 fCommonHistList->Add(fProNtracksLeadingJet);
1999 fCommonHistList->Add(fProDelR80pcPt);
2000 for(int ii=0; ii<5; ii++) fCommonHistList->Add(fProDelRPtSum[ii]);
2001
2002 if(genJets && genTracks){
2003 fCommonHistList->Add(fProNtracksLeadingJetGen);
2004 fCommonHistList->Add(fProDelR80pcPtGen);
2005 for(Int_t ii=0; ii<5; ii++) fCommonHistList->Add(fProDelRPtSumGen[ii]);
2006 }
2007
2008 if(fBckgMode){
2009 fCommonHistList->Add(fProNtracksLeadingJetBgrPerp2);
2010 for(Int_t ii=0; ii<5; ii++) fCommonHistList->Add(fProDelRPtSumBgrPerp2[ii]);
2011 }
2012
2013 if(fEffMode){
2014 fCommonHistList->Add(fProNtracksLeadingJetRecPrim);
2015 fCommonHistList->Add(fProDelR80pcPtRecPrim);
2016 for(Int_t ii=0; ii<5; ii++) fCommonHistList->Add(fProDelRPtSumRecPrim[ii]);
2017
2018 fCommonHistList->Add(fProNtracksLeadingJetRecSecNS);
2019 for(Int_t ii=0; ii<5; ii++) fCommonHistList->Add(fProDelRPtSumRecSecNS[ii]);
2020
2021 fCommonHistList->Add(fProNtracksLeadingJetRecSecS);
2022 for(Int_t ii=0; ii<5; ii++) fCommonHistList->Add(fProDelRPtSumRecSecS[ii]);
2023
2024 fCommonHistList->Add(fProNtracksLeadingJetRecSecSsc);
2025 for(Int_t ii=0; ii<5; ii++) fCommonHistList->Add(fProDelRPtSumRecSecSsc[ii]);
2026 }
2027 }
2028
2029 // =========== Switch on Sumw2 for all histos ===========
2030 for (Int_t i=0; i<fCommonHistList->GetEntries(); ++i){
2031 TH1 *h1 = dynamic_cast<TH1*>(fCommonHistList->At(i));
2032 if (h1) h1->Sumw2();
2033 else{
2034 THnSparse *hnSparse = dynamic_cast<THnSparse*>(fCommonHistList->At(i));
2035 if(hnSparse) hnSparse->Sumw2();
2036 }
2037 }
2038
2039 TH1::AddDirectory(oldStatus);
2040
2041 PostData(1, fCommonHistList);
2042}
2043
2044//_______________________________________________
2045void AliAnalysisTaskFragmentationFunction::Init()
2046{
2047 // Initialization
2048 if(fDebug > 1) Printf("AliAnalysisTaskFragmentationFunction::Init()");
2049
2050}
2051
2052//_____________________________________________________________
2053void AliAnalysisTaskFragmentationFunction::UserExec(Option_t *)
2054{
2055 // Main loop
2056 // Called for each event
2057 if(fDebug > 1) Printf("AliAnalysisTaskFragmentationFunction::UserExec()");
2058
2059
2060 if(fDebug > 1) Printf("Analysis event #%5d", (Int_t) fEntry);
2061
2062 // Trigger selection
2063 AliInputEventHandler* inputHandler = (AliInputEventHandler*)
2064 ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
2065
2066 if(!(inputHandler->IsEventSelected() & fEvtSelectionMask)){
2067 fh1EvtSelection->Fill(1.);
2068 if (fDebug > 1 ) Printf(" Trigger Selection: event REJECTED ... ");
2069 PostData(1, fCommonHistList);
2070 return;
2071 }
2072
2073 fESD = dynamic_cast<AliESDEvent*>(InputEvent());
2074 if(!fESD){
2075 if(fDebug>3) Printf("%s:%d ESDEvent not found in the input", (char*)__FILE__,__LINE__);
2076 }
2077
2078 fMCEvent = MCEvent();
2079 if(!fMCEvent){
2080 if(fDebug>3) Printf("%s:%d MCEvent not found in the input", (char*)__FILE__,__LINE__);
2081 }
2082
2083 // get AOD event from input/ouput
2084 TObject* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
2085 if( handler && handler->InheritsFrom("AliAODInputHandler") ) {
2086 fAOD = ((AliAODInputHandler*)handler)->GetEvent();
2087 if(fUseAODInputJets) fAODJets = fAOD;
2088 if (fDebug > 1) Printf("%s:%d AOD event from input", (char*)__FILE__,__LINE__);
2089 }
2090 else {
2091 handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
2092 if( handler && handler->InheritsFrom("AliAODHandler") ) {
2093 fAOD = ((AliAODHandler*)handler)->GetAOD();
2094 fAODJets = fAOD;
2095 if (fDebug > 1) Printf("%s:%d AOD event from output", (char*)__FILE__,__LINE__);
2096 }
2097 }
2098
2099 if(!fAODJets && !fUseAODInputJets){ // case we have AOD in input & output and want jets from output
2100 TObject* outHandler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
2101 if( outHandler && outHandler->InheritsFrom("AliAODHandler") ) {
2102 fAODJets = ((AliAODHandler*)outHandler)->GetAOD();
2103 if (fDebug > 1) Printf("%s:%d jets from output AOD", (char*)__FILE__,__LINE__);
2104 }
2105 }
2106
2107 if(fNonStdFile.Length()!=0){
2108 // case we have an AOD extension - fetch the jets from the extended output
2109
2110 AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
2111 fAODExtension = (aodH?aodH->GetExtension(fNonStdFile.Data()):0);
2112 if(!fAODExtension){
2113 if(fDebug>1)Printf("AODExtension not found for %s",fNonStdFile.Data());
2114 }
2115 }
2116
2117 if(!fAOD){
2118 Printf("%s:%d AODEvent not found", (char*)__FILE__,__LINE__);
2119 return;
2120 }
2121 if(!fAODJets){
2122 Printf("%s:%d AODEvent with jet branch not found", (char*)__FILE__,__LINE__);
2123 return;
2124 }
2125
2126
2127 // event selection **************************************************
2128 // *** event class ***
2129 Double_t centPercent = -1;
2130 if(fEventClass>0){
2131 Int_t cl = 0;
2132 if(handler->InheritsFrom("AliAODInputHandler")){
2133 // since it is not supported by the helper task define own classes
2134 centPercent = fAOD->GetHeader()->GetCentrality();
2135 cl = 1;
2136 if(centPercent>10) cl = 2;
2137 if(centPercent>30) cl = 3;
2138 if(centPercent>50) cl = 4;
2139 }
2140 else {
2141 cl = AliAnalysisHelperJetTasks::EventClass();
2142 if(fESD) centPercent = fESD->GetCentrality()->GetCentralityPercentile("V0M"); // retrieve value 'by hand'
2143 }
2144
2145 if(cl!=fEventClass){
2146 // event not in selected event class, reject event
2147 if (fDebug > 1) Printf("%s:%d event not in selected event class: event REJECTED ...",(char*)__FILE__,__LINE__);
2148 fh1EvtSelection->Fill(2.);
2149 PostData(1, fCommonHistList);
2150 return;
2151 }
2152 }
2153
2154 // *** vertex cut ***
2155 AliAODVertex* primVtx = fAOD->GetPrimaryVertex();
2156 Int_t nTracksPrim = primVtx->GetNContributors();
2157 fh1VertexNContributors->Fill(nTracksPrim);
2158
2159
2160 if (fDebug > 1) Printf("%s:%d primary vertex selection: %d", (char*)__FILE__,__LINE__,nTracksPrim);
2161 if(!nTracksPrim){
2162 if (fDebug > 1) Printf("%s:%d primary vertex selection: event REJECTED...",(char*)__FILE__,__LINE__);
2163 fh1EvtSelection->Fill(3.);
2164 PostData(1, fCommonHistList);
2165 return;
2166 }
2167
2168 fh1VertexZ->Fill(primVtx->GetZ());
2169
2170 if(TMath::Abs(primVtx->GetZ())>fMaxVertexZ){
2171 if (fDebug > 1) Printf("%s:%d primary vertex z = %f: event REJECTED...",(char*)__FILE__,__LINE__,primVtx->GetZ());
2172 fh1EvtSelection->Fill(4.);
2173 PostData(1, fCommonHistList);
2174 return;
2175 }
2176
2177 TString primVtxName(primVtx->GetName());
2178
2179 if(primVtxName.CompareTo("TPCVertex",TString::kIgnoreCase) == 1){
2180 if (fDebug > 1) Printf("%s:%d primary vertex selection: TPC vertex, event REJECTED...",(char*)__FILE__,__LINE__);
2181 fh1EvtSelection->Fill(5.);
2182 PostData(1, fCommonHistList);
2183 return;
2184 }
2185
2186 if (fDebug > 1) Printf("%s:%d event ACCEPTED ...",(char*)__FILE__,__LINE__);
2187 fh1EvtSelection->Fill(0.);
2188 fh1EvtCent->Fill(centPercent);
2189
2190
2191 //___ get MC information __________________________________________________________________
2192
2193 fh1Trials->Fill("#sum{ntrials}",fAvgTrials);
2194
2195 Double_t ptHard = 0.;
2196 Double_t nTrials = 1; // trials for MC trigger weight for real data
2197
2198 if(fMCEvent){
2199 AliGenEventHeader* genHeader = fMCEvent->GenEventHeader();
2200
2201 if(genHeader){
2202
2203 AliGenPythiaEventHeader* pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
2204 AliGenHijingEventHeader* hijingGenHeader = 0x0;
2205
2206 if(pythiaGenHeader){
2207 if(fDebug>3) Printf("%s:%d pythiaGenHeader found", (char*)__FILE__,__LINE__);
2208 nTrials = pythiaGenHeader->Trials();
2209 ptHard = pythiaGenHeader->GetPtHard();
2210
2211 fh1PtHard->Fill(ptHard);
2212 fh1PtHardTrials->Fill(ptHard,nTrials);
2213
2214 } else { // no pythia, hijing?
2215
2216 if(fDebug>3) Printf("%s:%d no pythiaGenHeader found", (char*)__FILE__,__LINE__);
2217
2218 hijingGenHeader = dynamic_cast<AliGenHijingEventHeader*>(genHeader);
2219 if(!hijingGenHeader){
2220 Printf("%s:%d no pythiaGenHeader or hjingGenHeader found", (char*)__FILE__,__LINE__);
2221 } else {
2222 if(fDebug>3) Printf("%s:%d hijingGenHeader found", (char*)__FILE__,__LINE__);
2223 }
2224 }
2225
2226 //fh1Trials->Fill("#sum{ntrials}",fAvgTrials);
2227 }
2228 }
2229
2230 //___ fetch jets __________________________________________________________________________
2231
2232 Int_t nJ = GetListOfJets(fJetsRec, kJetsRec);
2233 Int_t nRecJets = 0;
2234 if(nJ>=0) nRecJets = fJetsRec->GetEntries();
2235 if(fDebug>2)Printf("%s:%d Selected Rec jets: %d %d",(char*)__FILE__,__LINE__,nJ,nRecJets);
2236 if(nJ != nRecJets) Printf("%s:%d Mismatch Selected Rec Jets: %d %d",(char*)__FILE__,__LINE__,nJ,nRecJets);
2237
2238 Int_t nJCuts = GetListOfJets(fJetsRecCuts, kJetsRecAcceptance);
2239 Int_t nRecJetsCuts = 0;
2240 if(nJCuts>=0) nRecJetsCuts = fJetsRecCuts->GetEntries();
2241 if(fDebug>2)Printf("%s:%d Selected Rec jets after cuts: %d %d",(char*)__FILE__,__LINE__,nJCuts,nRecJetsCuts);
2242 if(nRecJetsCuts != nJCuts) Printf("%s:%d Mismatch selected Rec jets after cuts: %d %d",(char*)__FILE__,__LINE__,nJCuts,nRecJetsCuts);
2243 fh1nRecJetsCuts->Fill(nRecJetsCuts);
2244
2245 if(fJetTypeGen==kJetsKine || fJetTypeGen == kJetsKineAcceptance) fJetsGen->SetOwner(kTRUE); // kine aod jets allocated on heap, delete them with TList::Clear()
2246
2247 Int_t nJGen = GetListOfJets(fJetsGen, fJetTypeGen);
2248 Int_t nGenJets = 0;
2249 if(nJGen>=0) nGenJets = fJetsGen->GetEntries();
2250 if(fDebug>2)Printf("%s:%d Selected Gen jets: %d %d",(char*)__FILE__,__LINE__,nJGen,nGenJets);
2251
2252 if(nJGen != nGenJets) Printf("%s:%d Mismatch selected Gen jets: %d %d",(char*)__FILE__,__LINE__,nJGen,nGenJets);
2253 fh1nGenJets->Fill(nGenJets);
2254
2255
2256 if(fJetTypeRecEff==kJetsKine || fJetTypeRecEff == kJetsKineAcceptance) fJetsRecEff->SetOwner(kTRUE); // kine aod jets allocated on heap, delete them with TList::Clear()
2257 Int_t nJRecEff = GetListOfJets(fJetsRecEff, fJetTypeRecEff);
2258 Int_t nRecEffJets = 0;
2259 if(nJRecEff>=0) nRecEffJets = fJetsRecEff->GetEntries();
2260 if(fDebug>2)Printf("%s:%d Selected RecEff jets: %d %d",(char*)__FILE__,__LINE__,nJRecEff,nRecEffJets);
2261 if(nJRecEff != nRecEffJets) Printf("%s:%d Mismatch selected RecEff jets: %d %d",(char*)__FILE__,__LINE__,nJRecEff,nRecEffJets);
2262 fh1nRecEffJets->Fill(nRecEffJets);
2263
2264
2265 Int_t nEmbeddedJets = 0;
2266 TArrayI iEmbeddedMatchIndex;
2267 TArrayF fEmbeddedPtFraction;
2268
2269
2270 if(fBranchEmbeddedJets.Length()){
2271 Int_t nJEmbedded = GetListOfJets(fJetsEmbedded, kJetsEmbedded);
2272 if(nJEmbedded>=0) nEmbeddedJets = fJetsEmbedded->GetEntries();
2273 if(fDebug>2)Printf("%s:%d Selected Embedded jets: %d %d",(char*)__FILE__,__LINE__,nJEmbedded,nEmbeddedJets);
2274 if(nJEmbedded != nEmbeddedJets) Printf("%s:%d Mismatch Selected Embedded Jets: %d %d",(char*)__FILE__,__LINE__,nJEmbedded,nEmbeddedJets);
2275 fh1nEmbeddedJets->Fill(nEmbeddedJets);
2276
2277 Float_t maxDist = 0.3;
2278
2279 iEmbeddedMatchIndex.Set(nEmbeddedJets);
2280 fEmbeddedPtFraction.Set(nEmbeddedJets);
2281
2282 iEmbeddedMatchIndex.Reset(-1);
2283 fEmbeddedPtFraction.Reset(0);
2284
2285 AliAnalysisHelperJetTasks::GetJetMatching(fJetsEmbedded, nEmbeddedJets,
2286 fJetsRecCuts, nRecJetsCuts,
2287 iEmbeddedMatchIndex, fEmbeddedPtFraction,
2288 fDebug, maxDist);
2289
2290 }
2291
2292 //____ fetch background clusters ___________________________________________________
2293 if(fBckgMode &&
2294 (fBckgType[0]==kBckgClusters || fBckgType[1]==kBckgClusters || fBckgType[2]==kBckgClusters || fBckgType[3]==kBckgClusters || fBckgType[4]==kBckgClusters ||
2295 fBckgType[0]==kBckgClustersOutLeading || fBckgType[1]==kBckgClustersOutLeading || fBckgType[2]==kBckgClustersOutLeading ||
2296 fBckgType[3]==kBckgClustersOutLeading || fBckgType[4]==kBckgClustersOutLeading)){
2297
2298 Int_t nBJ = GetListOfBckgJets(fBckgJetsRec, kJetsRec);
2299 Int_t nRecBckgJets = 0;
2300 if(nBJ>=0) nRecBckgJets = fBckgJetsRec->GetEntries();
2301 if(fDebug>2)Printf("%s:%d Selected Rec background jets: %d %d",(char*)__FILE__,__LINE__,nBJ,nRecBckgJets);
2302 if(nBJ != nRecBckgJets) Printf("%s:%d Mismatch Selected Rec background jets: %d %d",(char*)__FILE__,__LINE__,nBJ,nRecBckgJets);
2303
2304 Int_t nBJCuts = GetListOfBckgJets(fBckgJetsRecCuts, kJetsRecAcceptance);
2305 Int_t nRecBckgJetsCuts = 0;
2306 if(nBJCuts>=0) nRecBckgJetsCuts = fBckgJetsRecCuts->GetEntries();
2307 if(fDebug>2)Printf("%s:%d Selected Rec background jets after cuts: %d %d",(char*)__FILE__,__LINE__,nJCuts,nRecJetsCuts);
2308 if(nRecBckgJetsCuts != nBJCuts) Printf("%s:%d Mismatch selected Rec background jets after cuts: %d %d",(char*)__FILE__,__LINE__,nBJCuts,nRecBckgJetsCuts);
2309 fh1nRecBckgJetsCuts->Fill(nRecBckgJetsCuts);
2310
2311 if(0){ // protection OB - not yet implemented
2312 if(fJetTypeGen==kJetsKine || fJetTypeGen == kJetsKineAcceptance) fBckgJetsGen->SetOwner(kTRUE); // kine aod jets allocated on heap, delete them with TList::Clear()
2313 Int_t nBJGen = GetListOfBckgJets(fBckgJetsGen, fJetTypeGen);
2314 Int_t nGenBckgJets = 0;
2315 if(nBJGen>=0) nGenBckgJets = fBckgJetsGen->GetEntries();
2316 if(fDebug>2)Printf("%s:%d Selected Gen background jets: %d %d",(char*)__FILE__,__LINE__,nBJGen,nGenBckgJets);
2317 if(nBJGen != nGenBckgJets) Printf("%s:%d Mismatch selected Gen background jets: %d %d",(char*)__FILE__,__LINE__,nBJGen,nGenBckgJets);
2318 fh1nGenBckgJets->Fill(nGenBckgJets);
2319 }
2320 }
2321
2322
2323 //____ fetch particles __________________________________________________________
2324
2325 Int_t nTCuts;
2326 if(fUseExtraTracks == 1) nTCuts = GetListOfTracks(fTracksRecCuts, kTrackAODExtraCuts);
2327 else if(fUseExtraTracks == -1) nTCuts = GetListOfTracks(fTracksRecCuts, kTrackAODExtraonlyCuts);
2328 else nTCuts = GetListOfTracks(fTracksRecCuts, kTrackAODCuts);
2329
2330 Int_t nRecPartCuts = 0;
2331 if(nTCuts>=0) nRecPartCuts = fTracksRecCuts->GetEntries();
2332 if(fDebug>2)Printf("%s:%d Selected Rec tracks after cuts: %d %d",(char*)__FILE__,__LINE__,nTCuts,nRecPartCuts);
2333 if(nRecPartCuts != nTCuts) Printf("%s:%d Mismatch selected Rec tracks after cuts: %d %d",(char*)__FILE__,__LINE__,nTCuts,nRecPartCuts);
2334 fh1EvtMult->Fill(nRecPartCuts);
2335
2336
2337 Int_t nTGen = GetListOfTracks(fTracksGen,fTrackTypeGen);
2338 Int_t nGenPart = 0;
2339 if(nTGen>=0) nGenPart = fTracksGen->GetEntries();
2340 if(fDebug>2)Printf("%s:%d Selected Gen tracks: %d %d",(char*)__FILE__,__LINE__,nTGen,nGenPart);
2341 if(nGenPart != nTGen) Printf("%s:%d Mismatch selected Gen tracks: %d %d",(char*)__FILE__,__LINE__,nTGen,nGenPart);
2342
2343
2344 //____ analysis, fill histos ___________________________________________________
2345
2346 if(fQAMode){
2347 // loop over tracks
2348 if(fQAMode&1){
2349 for(Int_t it=0; it<nRecPartCuts; ++it){
2350 AliVParticle *part = dynamic_cast<AliVParticle*>(fTracksRecCuts->At(it));
2351 if(part)fQATrackHistosRecCuts->FillTrackQA( part->Eta(), TVector2::Phi_0_2pi(part->Phi()), part->Pt() );
2352 }
2353 for(Int_t it=0; it<nGenPart; ++it){
2354 AliVParticle *part = dynamic_cast<AliVParticle*>(fTracksGen->At(it));
2355 if(part)fQATrackHistosGen->FillTrackQA( part->Eta(), TVector2::Phi_0_2pi(part->Phi()), part->Pt());
2356 }
2357 }
2358
2359 // loop over jets
2360 if(fQAMode&2){
2361 for(Int_t ij=0; ij<nRecJets; ++ij){
2362 AliAODJet* jet = dynamic_cast<AliAODJet*>(fJetsRec->At(ij));
2363 if(jet)fQAJetHistosRec->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt());
2364 }
2365 }
2366 }
2367
2368 if(fQAMode || fFFMode){
2369 for(Int_t ij=0; ij<nRecJetsCuts; ++ij){
2370
2371 AliAODJet* jet = (AliAODJet*)(fJetsRecCuts->At(ij));
2372 if(fQAMode&2) fQAJetHistosRecCuts->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt());
37f35567 2373
2374 if(fQAMode&2 && ij==0) fQAJetHistosRecCutsLeading->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt() );
ffa175ff 2375
37f35567 2376 Double_t ptFractionEmbedded = 0;
2377 AliAODJet* embeddedJet = 0;
ffa175ff 2378
37f35567 2379 if(fBranchEmbeddedJets.Length()){ // find embedded jet
ffa175ff 2380
37f35567 2381 Int_t indexEmbedded = -1;
2382 for(Int_t i=0; i<nEmbeddedJets; i++){
2383 if(iEmbeddedMatchIndex[i] == ij){
2384 indexEmbedded = i;
2385 ptFractionEmbedded = fEmbeddedPtFraction[i];
ffa175ff 2386 }
2387 }
2388
37f35567 2389 fh1IndexEmbedded->Fill(indexEmbedded);
2390 fh1FractionPtEmbedded->Fill(ptFractionEmbedded);
ffa175ff 2391
37f35567 2392 if(indexEmbedded>-1){
ffa175ff 2393
37f35567 2394 embeddedJet = dynamic_cast<AliAODJet*>(fJetsEmbedded->At(indexEmbedded));
2395 if(!embeddedJet) continue;
ffa175ff 2396
37f35567 2397 Double_t deltaPt = jet->Pt() - embeddedJet->Pt();
2398 Double_t deltaR = jet->DeltaR((AliVParticle*) (embeddedJet));
ffa175ff 2399
37f35567 2400 fh2DeltaPtVsJetPtEmbedded->Fill(embeddedJet->Pt(),deltaPt);
2401 fh2DeltaPtVsRecJetPtEmbedded->Fill(jet->Pt(),deltaPt);
2402 fh1DeltaREmbedded->Fill(deltaR);
ffa175ff 2403 }
37f35567 2404 }
2405
2406 // get tracks in jet
2407 TList* jettracklist = new TList();
2408 Double_t sumPt = 0.;
2409 Bool_t isBadJet = kFALSE;
2410
2411 if(GetFFRadius()<=0){
2412 GetJetTracksTrackrefs(jettracklist, jet, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJet);
2413 } else {
2414 if(fUseEmbeddedJetAxis){
2415 if(embeddedJet) GetJetTracksPointing(fTracksRecCuts, jettracklist, embeddedJet,
2416 GetFFRadius(), sumPt, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJet);
2417 }
2418 else GetJetTracksPointing(fTracksRecCuts, jettracklist, jet,
2419 GetFFRadius(), sumPt, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJet);
2420 }
ffa175ff 2421
37f35567 2422 if(GetFFMinNTracks()>0 && jettracklist->GetSize()<=GetFFMinNTracks()) isBadJet = kTRUE;
2423
2424 if(isBadJet) continue;
ffa175ff 2425
37f35567 2426 if(ptFractionEmbedded>=fCutFractionPtEmbedded){ // if no embedding: ptFraction = cutFraction = 0
2427
ffa175ff 2428 for(Int_t it=0; it<jettracklist->GetSize(); ++it){
2429
2430 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(jettracklist->At(it));
2431 if(!trackVP)continue;
37f35567 2432
2433 AliAODTrack * aodtrack = dynamic_cast<AliAODTrack*>(jettracklist->At(it));
2434 if(!aodtrack) continue;
2435
ffa175ff 2436 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
2437
2438 Float_t jetPt = jet->Pt();
37f35567 2439 if(fUseEmbeddedJetPt){
2440 if(embeddedJet) jetPt = embeddedJet->Pt();
2441 else jetPt = 0;
2442 }
ffa175ff 2443 Float_t trackPt = trackV->Pt();
37f35567 2444
ffa175ff 2445 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
37f35567 2446
2447 if(fFFMode && (ij==0)) fFFHistosRecCuts->FillFF(trackPt, jetPt, incrementJetPt);
2448 if(fFFMode) fFFHistosRecCutsInc->FillFF(trackPt, jetPt, incrementJetPt);
ffa175ff 2449
37f35567 2450 if(it==0){ // leading track
2451 if(fFFMode) fFFHistosRecLeadingTrack->FillFF( trackPt, jetPt, kTRUE);
2452 }
ffa175ff 2453
2454 delete trackV;
2455 }
37f35567 2456
2457 // background ff
2458 if(fBckgMode && (ij==0)){
ffa175ff 2459 if(fBckgType[0]!=kBckgNone)
37f35567 2460 FillBckgHistos(fBckgType[0], fTracksRecCuts, fJetsRecCuts, jet,
2461 fFFBckgHisto0RecCuts,fQABckgHisto0RecCuts, fh1BckgMult0);
ffa175ff 2462 if(fBckgType[1]!=kBckgNone)
37f35567 2463 FillBckgHistos(fBckgType[1], fTracksRecCuts, fJetsRecCuts, jet,
2464 fFFBckgHisto1RecCuts,fQABckgHisto1RecCuts, fh1BckgMult1);
ffa175ff 2465 if(fBckgType[2]!=kBckgNone)
37f35567 2466 FillBckgHistos(fBckgType[2], fTracksRecCuts, fJetsRecCuts, jet,
2467 fFFBckgHisto2RecCuts,fQABckgHisto2RecCuts, fh1BckgMult2);
ffa175ff 2468 if(fBckgType[3]!=kBckgNone)
37f35567 2469 FillBckgHistos(fBckgType[3], fTracksRecCuts, fJetsRecCuts, jet,
2470 fFFBckgHisto3RecCuts,fQABckgHisto3RecCuts, fh1BckgMult3);
ffa175ff 2471 if(fBckgType[4]!=kBckgNone)
37f35567 2472 FillBckgHistos(fBckgType[4], fTracksRecCuts, fJetsRecCuts, jet,
2473 fFFBckgHisto4RecCuts,fQABckgHisto4RecCuts, fh1BckgMult4);
ffa175ff 2474 } // end if(fBckgMode)
2475
2476
37f35567 2477 if(fJSMode && (ij==0)) FillJetShape(jet, jettracklist, fProNtracksLeadingJet, fProDelRPtSum, fProDelR80pcPt);
2478
2479 delete jettracklist;
2480
2481 } // end: cut embedded ratio
2482 } // end: rec. jets after cuts
2483
2484 // generated jets
2485 for(Int_t ij=0; ij<nGenJets; ++ij){
2486
2487 AliAODJet* jet = dynamic_cast<AliAODJet*>(fJetsGen->At(ij));
2488 if(!jet)continue;
ffa175ff 2489
37f35567 2490 if(fQAMode&2) fQAJetHistosGen->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt());
2491
2492 if(fQAMode&2 && (ij==0)) fQAJetHistosGenLeading->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt());
2493
2494 TList* jettracklist = new TList();
2495 Double_t sumPt = 0.;
2496 Bool_t isBadJet = kFALSE;
2497
2498 if(GetFFRadius()<=0){
2499 GetJetTracksTrackrefs(jettracklist, jet, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJet);
2500 } else {
2501 GetJetTracksPointing(fTracksGen, jettracklist, jet, GetFFRadius(), sumPt, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJet);
ffa175ff 2502 }
37f35567 2503
2504 if(GetFFMinNTracks()>0 && jettracklist->GetSize()<=GetFFMinNTracks()) isBadJet = kTRUE;;
2505
2506 if(isBadJet) continue;
ffa175ff 2507
37f35567 2508 for(Int_t it=0; it<jettracklist->GetSize(); ++it){
2509
2510 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(jettracklist->At(it));
2511 if(!trackVP)continue;
2512 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
2513
2514 Float_t jetPt = jet->Pt();
2515 Float_t trackPt = trackV->Pt();
2516
2517 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
2518
2519 if(fFFMode && (ij==0)) fFFHistosGen->FillFF( trackPt, jetPt, incrementJetPt );
2520 if(fFFMode) fFFHistosGenInc->FillFF( trackPt, jetPt, incrementJetPt );
2521
2522 if(it==0){ // leading track
2523 if(fFFMode) fFFHistosGenLeadingTrack->FillFF( trackPt, jetPt, kTRUE );
2524 }
2525
2526 delete trackV;
2527 }
2528
2529 if(fBckgMode && (ij==0)){
2530 if(fBckgType[0]!=kBckgNone)
2531 FillBckgHistos(fBckgType[0], fTracksGen, fJetsGen, jet,
2532 fFFBckgHisto0Gen, fQABckgHisto0Gen);
2533 if(fBckgType[1]!=kBckgNone)
2534 FillBckgHistos(fBckgType[1], fTracksGen, fJetsGen, jet,
2535 fFFBckgHisto1Gen, fQABckgHisto1Gen);
2536 if(fBckgType[2]!=kBckgNone)
2537 FillBckgHistos(fBckgType[2], fTracksGen, fJetsGen, jet,
2538 fFFBckgHisto2Gen, fQABckgHisto2Gen);
2539 if(fBckgType[3]!=kBckgNone)
2540 FillBckgHistos(fBckgType[3], fTracksGen, fJetsGen, jet,
2541 fFFBckgHisto3Gen, fQABckgHisto3Gen);
2542 if(fBckgType[4]!=kBckgNone)
2543 FillBckgHistos(fBckgType[4], fTracksGen, fJetsGen, jet,
2544 fFFBckgHisto4Gen, fQABckgHisto4Gen);
2545 } // end if(fBckgMode)
2546
2547 if(fJSMode && (ij==0)) FillJetShape(jet, jettracklist, fProNtracksLeadingJetGen, fProDelRPtSumGen, fProDelR80pcPtGen);
ffa175ff 2548
37f35567 2549 delete jettracklist;
2550 }
2551 } // end: QA, FF and intra-jet
2552
2553
ffa175ff 2554 // ____ efficiency _______________________________
2555
2556 if(fEffMode && (fJetTypeRecEff != kJetsUndef)){
2557
2558 // arrays holding for each generated particle the reconstructed AOD track index & isPrimary flag, are initialized in AssociateGenRec(...) function
2559 TArrayI indexAODTr;
2560 TArrayS isGenPrim;
2561
2562 // array holding for each reconstructed AOD track generated particle index, initialized in AssociateGenRec(...) function
2563 TArrayI indexMCTr;
2564
2565 // ... and another set for secondaries from strange/non strange mothers (secondary MC tracks are stored in different lists)
2566 TArrayI indexAODTrSecNS;
2567 TArrayS isGenSecNS;
2568 TArrayI indexMCTrSecNS;
2569
2570 TArrayI indexAODTrSecS;
2571 TArrayS isGenSecS;
2572 TArrayI indexMCTrSecS;
2573
2574 Int_t nTracksAODMCCharged = GetListOfTracks(fTracksAODMCCharged, kTrackAODMCCharged);
2575 if(fDebug>2)Printf("%s:%d selected AODMC tracks: %d ",(char*)__FILE__,__LINE__,nTracksAODMCCharged);
2576
2577 Int_t nTracksAODMCChargedSecNS = GetListOfTracks(fTracksAODMCChargedSecNS, kTrackAODMCChargedSecNS);
2578 if(fDebug>2)Printf("%s:%d selected AODMC secondary tracks NS: %d ",(char*)__FILE__,__LINE__,nTracksAODMCChargedSecNS);
2579
2580 Int_t nTracksAODMCChargedSecS = GetListOfTracks(fTracksAODMCChargedSecS, kTrackAODMCChargedSecS);
2581 if(fDebug>2)Printf("%s:%d selected AODMC secondary tracks S: %d ",(char*)__FILE__,__LINE__,nTracksAODMCChargedSecS);
2582
2583 Int_t nTracksRecQualityCuts = GetListOfTracks(fTracksRecQualityCuts, kTrackAODQualityCuts);
2584 if(fDebug>2)Printf("%s:%d selected rec tracks quality after cuts, full acceptance/pt : %d ",(char*)__FILE__,__LINE__,nTracksRecQualityCuts);
2585
2586 // associate gen and rec tracks, store indices in TArrays
2587 AssociateGenRec(fTracksAODMCCharged,fTracksRecQualityCuts,indexAODTr,indexMCTr,isGenPrim,fh2PtRecVsGenPrim);
2588 AssociateGenRec(fTracksAODMCChargedSecNS,fTracksRecQualityCuts,indexAODTrSecNS,indexMCTrSecNS,isGenSecNS,fh2PtRecVsGenSec);
2589 AssociateGenRec(fTracksAODMCChargedSecS,fTracksRecQualityCuts,indexAODTrSecS,indexMCTrSecS,isGenSecS,fh2PtRecVsGenSec);
2590
2591 // single track eff
2592 if(fQAMode&1) FillSingleTrackHistosRecGen(fQATrackHistosRecEffGen,fQATrackHistosRecEffRec,fTracksAODMCCharged,indexAODTr,isGenPrim);
2593
2594 // secondaries
2595 if(fQAMode&1) FillSingleTrackHistosRecGen(0x0,fQATrackHistosSecRecNS,fTracksAODMCChargedSecNS,indexAODTrSecNS,isGenSecNS);
2596 if(fQAMode&1) FillSingleTrackHistosRecGen(0x0,fQATrackHistosSecRecS,fTracksAODMCChargedSecS,indexAODTrSecS,isGenSecS);
2597 if(fQAMode&1) FillSingleTrackHistosRecGen(0x0,fQATrackHistosSecRecSsc,fTracksAODMCChargedSecS,indexAODTrSecS,isGenSecS,kTRUE);
2598
2599
2600 // jet track eff
2601 Double_t sumPtGenLeadingJetRecEff = 0;
2602 Double_t sumPtGenLeadingJetSec = 0;
2603 Double_t sumPtRecLeadingJetRecEff = 0;
2604
2605 for(Int_t ij=0; ij<nRecEffJets; ++ij){ // jet loop
2606
2607 AliAODJet* jet = (AliAODJet*)(fJetsRecEff->At(ij));
2608
2609 Bool_t isBadJetGenPrim = kFALSE;
2610 Bool_t isBadJetGenSec = kFALSE;
2611 Bool_t isBadJetRec = kFALSE;
2612
2613
2614 if(ij==0){ // leading jet
2615
2616 // for efficiency: gen tracks from pointing with gen/rec jet
2617 TList* jettracklistGenPrim = new TList();
2618
2619 // if radius<0 -> trackRefs: collect gen tracks in wide radius + fill FF recEff rec histos with tracks contained in track refs
2620 // note : FF recEff gen histos will be somewhat useless in this approach
2621
2622 if(GetFFRadius() >0)
2623 GetJetTracksPointing(fTracksAODMCCharged, jettracklistGenPrim, jet, GetFFRadius(), sumPtGenLeadingJetRecEff, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJetGenPrim);
2624 else
2625 GetJetTracksPointing(fTracksAODMCCharged, jettracklistGenPrim, jet, TMath::Abs(GetFFRadius())+0.2, sumPtGenLeadingJetRecEff, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJetGenPrim);
2626
2627 TList* jettracklistGenSecNS = new TList();
2628 if(GetFFRadius() >0)
2629 GetJetTracksPointing(fTracksAODMCChargedSecNS, jettracklistGenSecNS, jet, GetFFRadius(), sumPtGenLeadingJetSec, GetFFMinLTrackPt() , GetFFMaxTrackPt(), isBadJetGenSec);
2630 else
2631 GetJetTracksPointing(fTracksAODMCChargedSecNS, jettracklistGenSecNS, jet, TMath::Abs(GetFFRadius())+0.2, sumPtGenLeadingJetSec, GetFFMinLTrackPt() , GetFFMaxTrackPt(), isBadJetGenSec);
2632
2633 TList* jettracklistGenSecS = new TList();
2634 if(GetFFRadius() >0)
2635 GetJetTracksPointing(fTracksAODMCChargedSecS, jettracklistGenSecS, jet, GetFFRadius(), sumPtGenLeadingJetSec, GetFFMinLTrackPt() , GetFFMaxTrackPt(), isBadJetGenSec);
2636 else
2637 GetJetTracksPointing(fTracksAODMCChargedSecS, jettracklistGenSecS, jet, TMath::Abs(GetFFRadius())+0.2, sumPtGenLeadingJetSec, GetFFMinLTrackPt() , GetFFMaxTrackPt(), isBadJetGenSec);
2638
2639
2640 // bin efficiency in jet pt bins using rec tracks
2641 TList* jettracklistRec = new TList();
2642 if(GetFFRadius() >0) GetJetTracksPointing(fTracksRecCuts,jettracklistRec, jet, GetFFRadius(), sumPtRecLeadingJetRecEff, GetFFMinLTrackPt() , GetFFMaxTrackPt(), isBadJetRec);
2643 else GetJetTracksTrackrefs(jettracklistRec, jet, GetFFMinLTrackPt() , GetFFMaxTrackPt(), isBadJetRec);
2644
2645
2646 Double_t jetEta = jet->Eta();
2647 Double_t jetPhi = TVector2::Phi_0_2pi(jet->Phi());
2648
2649 if(GetFFMinNTracks()>0 && jettracklistGenPrim->GetSize()<=GetFFMinNTracks()) isBadJetGenPrim = kTRUE;
2650 if(GetFFMinNTracks()>0 && jettracklistGenSecNS->GetSize()<=GetFFMinNTracks()) isBadJetGenSec = kTRUE;
2651 if(GetFFMinNTracks()>0 && jettracklistRec->GetSize()<=GetFFMinNTracks()) isBadJetRec = kTRUE;
2652
2653 if(isBadJetRec) continue;
2654
2655 if(fQAMode&2) fQAJetHistosRecEffLeading->FillJetQA( jetEta, jetPhi, sumPtGenLeadingJetRecEff );
2656
2657 if(fFFMode){
2658
2659 if(GetFFRadius()>0) FillJetTrackHistosRec(fFFHistosRecEffRec,jet,
2660 jettracklistGenPrim,fTracksAODMCCharged,fTracksRecQualityCuts,indexAODTr,isGenPrim,
2661 0,kFALSE,fJSMode,fProNtracksLeadingJetRecPrim,fProDelRPtSumRecPrim,fProDelR80pcPtRecPrim);
2662
2663 else FillJetTrackHistosRec(fFFHistosRecEffRec,jet,
2664 jettracklistGenPrim,fTracksAODMCCharged,fTracksRecQualityCuts,indexAODTr,isGenPrim,
2665 jettracklistRec,kFALSE,fJSMode,fProNtracksLeadingJetRecPrim,fProDelRPtSumRecPrim,fProDelR80pcPtRecPrim);
2666
2667
2668 // secondaries: use jet pt from primaries
2669 if(GetFFRadius()>0) FillJetTrackHistosRec(fFFHistosSecRecNS,jet,
2670 jettracklistGenSecNS,fTracksAODMCChargedSecNS,fTracksRecQualityCuts, indexAODTrSecNS,isGenSecNS,
2671 0,kFALSE,fJSMode,fProNtracksLeadingJetRecSecNS,fProDelRPtSumRecSecNS);
2672
2673 else FillJetTrackHistosRec(fFFHistosSecRecNS,jet,
2674 jettracklistGenSecNS,fTracksAODMCChargedSecNS,fTracksRecQualityCuts,indexAODTrSecNS,isGenSecNS,
2675 jettracklistRec,kFALSE,fJSMode,fProNtracksLeadingJetRecSecNS,fProDelRPtSumRecSecNS);
2676
2677 if(GetFFRadius()>0) FillJetTrackHistosRec(fFFHistosSecRecS,jet,
2678 jettracklistGenSecS,fTracksAODMCChargedSecS,fTracksRecQualityCuts,indexAODTrSecS,isGenSecS,
2679 0,kFALSE,fJSMode,fProNtracksLeadingJetRecSecS,fProDelRPtSumRecSecS);
2680
2681 else FillJetTrackHistosRec(fFFHistosSecRecS,jet,
2682 jettracklistGenSecS,fTracksAODMCChargedSecS,fTracksRecQualityCuts,indexAODTrSecS,isGenSecS,
2683 jettracklistRec,kFALSE,fJSMode,fProNtracksLeadingJetRecSecS,fProDelRPtSumRecSecS);
2684
2685 if(GetFFRadius()>0) FillJetTrackHistosRec(fFFHistosSecRecSsc,jet,
2686 jettracklistGenSecS,fTracksAODMCChargedSecS,fTracksRecQualityCuts,indexAODTrSecS,isGenSecS,
2687 0,kTRUE,fJSMode,fProNtracksLeadingJetRecSecSsc,fProDelRPtSumRecSecSsc);
2688
2689 else FillJetTrackHistosRec(fFFHistosSecRecSsc,jet,
2690 jettracklistGenSecS,fTracksAODMCChargedSecS,fTracksRecQualityCuts,indexAODTrSecS,isGenSecS,
2691 jettracklistRec,kTRUE,fJSMode,fProNtracksLeadingJetRecSecSsc,fProDelRPtSumRecSecSsc);
2692 }
2693
2694 delete jettracklistGenPrim;
2695 delete jettracklistGenSecNS;
2696 delete jettracklistGenSecS;
2697 delete jettracklistRec;
2698
2699
2700 if(fBckgMode && fFFMode){
2701
2702 TList* perpjettracklistGen = new TList();
2703 TList* perpjettracklistGen1 = new TList();
2704 TList* perpjettracklistGen2 = new TList();
2705
2706 Double_t sumPtGenPerp = 0.;
2707 Double_t sumPtGenPerp1 = 0.;
2708 Double_t sumPtGenPerp2 = 0.;
2709 GetTracksTiltedwrpJetAxis(TMath::Pi()/2.,fTracksAODMCCharged, perpjettracklistGen1, jet, TMath::Abs(GetFFRadius()) , sumPtGenPerp1);
2710 GetTracksTiltedwrpJetAxis(-1*TMath::Pi()/2.,fTracksAODMCCharged, perpjettracklistGen2, jet, TMath::Abs(GetFFRadius()) , sumPtGenPerp2);
2711
2712 perpjettracklistGen->AddAll(perpjettracklistGen1);
2713 perpjettracklistGen->AddAll(perpjettracklistGen2);
2714 sumPtGenPerp = 0.5*(sumPtGenPerp1+sumPtGenPerp2);
2715
2716 TList* perpjettracklistGenSecNS = new TList();
2717 TList* perpjettracklistGenSecNS1 = new TList();
2718 TList* perpjettracklistGenSecNS2 = new TList();
2719
2720 Double_t sumPtGenPerpNS;
2721 Double_t sumPtGenPerpNS1;
2722 Double_t sumPtGenPerpNS2;
2723 GetTracksTiltedwrpJetAxis(TMath::Pi()/2.,fTracksAODMCChargedSecNS, perpjettracklistGenSecNS1, jet, TMath::Abs(GetFFRadius()) , sumPtGenPerpNS1);
2724 GetTracksTiltedwrpJetAxis(-1*TMath::Pi()/2.,fTracksAODMCChargedSecNS, perpjettracklistGenSecNS2, jet, TMath::Abs(GetFFRadius()) , sumPtGenPerpNS2);
2725
2726 perpjettracklistGenSecNS->AddAll(perpjettracklistGenSecNS1);
2727 perpjettracklistGenSecNS->AddAll(perpjettracklistGenSecNS2);
2728 sumPtGenPerpNS = 0.5*(sumPtGenPerpNS1+sumPtGenPerpNS2);
2729
2730
2731 TList* perpjettracklistGenSecS = new TList();
2732 TList* perpjettracklistGenSecS1 = new TList();
2733 TList* perpjettracklistGenSecS2 = new TList();
2734
2735 Double_t sumPtGenPerpS;
2736 Double_t sumPtGenPerpS1;
2737 Double_t sumPtGenPerpS2;
2738 GetTracksTiltedwrpJetAxis(TMath::Pi()/2.,fTracksAODMCChargedSecS, perpjettracklistGenSecS1, jet, TMath::Abs(GetFFRadius()) , sumPtGenPerpS1);
2739 GetTracksTiltedwrpJetAxis(-1*TMath::Pi()/2.,fTracksAODMCChargedSecS, perpjettracklistGenSecS2, jet, TMath::Abs(GetFFRadius()) , sumPtGenPerpS2);
2740
2741 perpjettracklistGenSecS->AddAll(perpjettracklistGenSecS1);
2742 perpjettracklistGenSecS->AddAll(perpjettracklistGenSecS2);
2743 sumPtGenPerpS = 0.5*(sumPtGenPerpS1+sumPtGenPerpS2);
2744
2745
2746 if(perpjettracklistGen->GetSize() != perpjettracklistGen1->GetSize() + perpjettracklistGen2->GetSize()){
2747 cout<<" ERROR: perpjettracklistGen size "<<perpjettracklistGen->GetSize()<<" perp1 "<<perpjettracklistGen1->GetSize()
2748 <<" perp2 "<<perpjettracklistGen2->GetSize()<<endl;
2749 exit(0);
2750 }
2751
2752 if(perpjettracklistGenSecNS->GetSize() != perpjettracklistGenSecNS1->GetSize() + perpjettracklistGenSecNS2->GetSize()){
2753 cout<<" ERROR: perpjettracklistGenSecNS size "<<perpjettracklistGenSecNS->GetSize()<<" perp1 "<<perpjettracklistGenSecNS1->GetSize()
2754 <<" perp2 "<<perpjettracklistGenSecNS2->GetSize()<<endl;
2755 exit(0);
2756 }
2757
2758 if(perpjettracklistGenSecS->GetSize() != perpjettracklistGenSecS1->GetSize() + perpjettracklistGenSecS2->GetSize()){
2759 cout<<" ERROR: perpjettracklistGenSecS size "<<perpjettracklistGenSecS->GetSize()<<" perp1 "<<perpjettracklistGenSecS1->GetSize()
2760 <<" perp2 "<<perpjettracklistGenSecS2->GetSize()<<endl;
2761 exit(0);
2762 }
2763
2764
2765 FillJetTrackHistosRec(fFFBckgHisto0RecEffRec,jet,
2766 perpjettracklistGen,fTracksAODMCCharged,fTracksRecQualityCuts,indexAODTr,isGenPrim);
2767
2768 FillJetTrackHistosRec(fFFBckgHisto0SecRecNS,jet,
2769 perpjettracklistGenSecNS,fTracksAODMCChargedSecNS,fTracksRecQualityCuts,indexAODTrSecNS,isGenSecNS);
2770
2771 FillJetTrackHistosRec(fFFBckgHisto0SecRecS,jet,
2772 perpjettracklistGenSecS,fTracksAODMCChargedSecS,fTracksRecQualityCuts,indexAODTrSecS,isGenSecS);
2773
2774 FillJetTrackHistosRec(fFFBckgHisto0SecRecSsc,jet,
2775 perpjettracklistGenSecS,fTracksAODMCChargedSecS,fTracksRecQualityCuts,indexAODTrSecS,isGenSecS,0,kTRUE);
2776
2777
2778 delete perpjettracklistGen;
2779 delete perpjettracklistGen1;
2780 delete perpjettracklistGen2;
2781
2782 delete perpjettracklistGenSecNS;
2783 delete perpjettracklistGenSecNS1;
2784 delete perpjettracklistGenSecNS2;
2785
2786 delete perpjettracklistGenSecS;
2787 delete perpjettracklistGenSecS1;
2788 delete perpjettracklistGenSecS2;
2789
2790 }
2791 }
2792 }
2793 }
2794
2795 //___________________
2796
2797 fTracksRecCuts->Clear();
2798 fTracksGen->Clear();
2799 fTracksAODMCCharged->Clear();
2800 fTracksAODMCChargedSecNS->Clear();
2801 fTracksAODMCChargedSecS->Clear();
2802 fTracksRecQualityCuts->Clear();
2803
2804 fJetsRec->Clear();
2805 fJetsRecCuts->Clear();
2806 fJetsGen->Clear();
2807 fJetsRecEff->Clear();
2808 fJetsEmbedded->Clear();
2809
2810
2811 if(fBckgMode &&
2812 (fBckgType[0]==kBckgClusters || fBckgType[1]==kBckgClusters || fBckgType[2]==kBckgClusters || fBckgType[3]==kBckgClusters || fBckgType[4]==kBckgClusters ||
2813 fBckgType[0]==kBckgClustersOutLeading || fBckgType[1]==kBckgClustersOutLeading || fBckgType[2]==kBckgClustersOutLeading ||
2814 fBckgType[3]==kBckgClustersOutLeading || fBckgType[4]==kBckgClustersOutLeading)){
2815
2816 fBckgJetsRec->Clear();
2817 fBckgJetsRecCuts->Clear();
2818 fBckgJetsGen->Clear();
2819 }
2820
2821
2822 //Post output data.
2823 PostData(1, fCommonHistList);
2824}
2825
2826//______________________________________________________________
2827void AliAnalysisTaskFragmentationFunction::Terminate(Option_t *)
2828{
2829 // terminated
2830
2831 if(fDebug > 1) printf("AliAnalysisTaskFragmentationFunction::Terminate() \n");
2832}
2833
2834//_________________________________________________________________________________
2835Int_t AliAnalysisTaskFragmentationFunction::GetListOfTracks(TList *list, Int_t type)
2836{
2837 // fill list of tracks selected according to type
2838
2839 if(fDebug > 2) Printf("%s:%d Selecting tracks with %d", (char*)__FILE__,__LINE__,type);
2840
2841 if(!list){
2842 if(fDebug>1) Printf("%s:%d no input list", (char*)__FILE__,__LINE__);
2843 return -1;
2844 }
2845
2846 if(!fAOD) return -1;
2847
2848 if(!fAOD->GetTracks()) return 0;
2849
2850 if(type==kTrackUndef) return 0;
2851
2852 Int_t iCount = 0;
2853
2854 if(type==kTrackAODExtraCuts || type==kTrackAODExtraonlyCuts || type==kTrackAODExtra || type==kTrackAODExtraonly){
2855
2856 TClonesArray *aodExtraTracks = dynamic_cast<TClonesArray*>(fAOD->FindListObject("aodExtraTracks"));
2857 if(!aodExtraTracks)return iCount;
2858 for(int it =0; it<aodExtraTracks->GetEntries(); it++) {
2859 AliVParticle *track = dynamic_cast<AliVParticle*> ((*aodExtraTracks)[it]);
2860 if (!track) continue;
2861
2862 AliAODTrack *tr = dynamic_cast<AliAODTrack*> (track);
2863 if(!tr)continue;
2864
2865 if(type==kTrackAODExtraCuts || type==kTrackAODExtraonlyCuts){
2866
2867 if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask))) continue;
2868
2869 if(tr->Eta() < fTrackEtaMin || tr->Eta() > fTrackEtaMax) continue;
2870 if(tr->Phi() < fTrackPhiMin || tr->Phi() > fTrackPhiMax) continue;
2871 if(tr->Pt() < fTrackPtCut) continue;
2872 }
2873
2874 list->Add(tr);
2875 iCount++;
2876 }
2877 }
2878
2879 if(type==kTrackAODCuts || type==kTrackAODQualityCuts || type==kTrackAOD || type==kTrackAODExtraCuts || type==kTrackAODExtra){
2880
2881 // all rec. tracks, esd filter mask, eta range
2882
2883 for(Int_t it=0; it<fAOD->GetNumberOfTracks(); ++it){
2884 AliAODTrack *tr = fAOD->GetTrack(it);
2885
2886 if(type == kTrackAODCuts || type==kTrackAODQualityCuts || type==kTrackAODExtraCuts){
2887
2888 if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask))) continue;
2889 if(type == kTrackAODCuts){
2890 if(tr->Eta() < fTrackEtaMin || tr->Eta() > fTrackEtaMax) continue;
2891 if(tr->Phi() < fTrackPhiMin || tr->Phi() > fTrackPhiMax) continue;
2892 if(tr->Pt() < fTrackPtCut) continue;
2893 }
2894 }
2895 list->Add(tr);
2896 iCount++;
2897 }
2898 }
2899 else if (type==kTrackKineAll || type==kTrackKineCharged || type==kTrackKineChargedAcceptance){
2900 // kine particles, all or rather charged
2901 if(!fMCEvent) return iCount;
2902
2903 for(Int_t it=0; it<fMCEvent->GetNumberOfTracks(); ++it){
2904 AliMCParticle* part = (AliMCParticle*) fMCEvent->GetTrack(it);
2905
2906 if(type == kTrackKineCharged || type == kTrackKineChargedAcceptance){
2907 if(part->Charge()==0) continue;
2908
2909 if(type == kTrackKineChargedAcceptance &&
2910 ( part->Eta() < fTrackEtaMin
2911 || part->Eta() > fTrackEtaMax
2912 || part->Phi() < fTrackPhiMin
2913 || part->Phi() > fTrackPhiMax
2914 || part->Pt() < fTrackPtCut)) continue;
2915 }
2916
2917 list->Add(part);
2918 iCount++;
2919 }
2920 }
2921 else if (type==kTrackAODMCCharged || type==kTrackAODMCAll || type==kTrackAODMCChargedAcceptance || type==kTrackAODMCChargedSecNS || type==kTrackAODMCChargedSecS) {
2922 // MC particles (from AOD), physical primaries, all or rather charged or rather charged within acceptance
2923 if(!fAOD) return -1;
2924
2925 TClonesArray *tca = dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName()));
2926 if(!tca)return iCount;
2927
2928 for(int it=0; it<tca->GetEntriesFast(); ++it){
2929 AliAODMCParticle *part = dynamic_cast<AliAODMCParticle*>(tca->At(it));
2930 if(!part)continue;
2931 if(type != kTrackAODMCChargedSecNS && type != kTrackAODMCChargedSecS && !part->IsPhysicalPrimary())continue;
2932 if((type == kTrackAODMCChargedSecNS || type == kTrackAODMCChargedSecS) && part->IsPhysicalPrimary())continue;
2933
2934 if (type==kTrackAODMCCharged || type==kTrackAODMCChargedAcceptance || type==kTrackAODMCChargedSecNS || type==kTrackAODMCChargedSecS){
2935 if(part->Charge()==0) continue;
2936
2937 if(type==kTrackAODMCChargedSecNS || type==kTrackAODMCChargedSecS){
2938 Bool_t isFromStrange = kFALSE;
2939 Int_t iMother = part->GetMother();
ffa175ff 2940
47830b3f 2941 if(iMother < 0) continue; // throw out PYTHIA stack partons + incoming protons
2942
2943 AliAODMCParticle *partM = dynamic_cast<AliAODMCParticle*>(tca->At(iMother));
2944 if(!partM) continue;
2945
2946 Int_t codeM = TMath::Abs(partM->GetPdgCode());
2947 Int_t mfl = Int_t (codeM/ TMath::Power(10, Int_t(TMath::Log10(codeM))));
2948 if (mfl == 3 && codeM != 3) isFromStrange = kTRUE;
ffa175ff 2949
47830b3f 2950 // if(mfl ==3){
2951 // cout<<" mfl "<<mfl<<" codeM "<<partM->GetPdgCode()<<" code this track "<<part->GetPdgCode()<<endl;
2952 // cout<<" index this track "<<it<<" index daughter 0 "<<partM->GetDaughter(0)<<" 1 "<<partM->GetDaughter(1)<<endl;
2953 // }
ffa175ff 2954
47830b3f 2955 if(type==kTrackAODMCChargedSecNS && isFromStrange) continue;
2956 if(type==kTrackAODMCChargedSecS && !isFromStrange) continue;
ffa175ff 2957 }
47830b3f 2958
ffa175ff 2959
2960 if(type==kTrackAODMCChargedAcceptance &&
2961 ( part->Eta() > fTrackEtaMax
2962 || part->Eta() < fTrackEtaMin
2963 || part->Phi() > fTrackPhiMax
2964 || part->Phi() < fTrackPhiMin
2965 || part->Pt() < fTrackPtCut)) continue;
2966 }
2967
2968 list->Add(part);
2969 iCount++;
2970 }
2971 }
2972
2973 list->Sort();
2974 return iCount;
2975
2976}
2977// _______________________________________________________________________________
2978Int_t AliAnalysisTaskFragmentationFunction::GetListOfJets(TList *list, Int_t type)
2979{
2980 // fill list of jets selected according to type
2981
2982 if(!list){
2983 if(fDebug>1) Printf("%s:%d no input list", (char*)__FILE__,__LINE__);
2984 return -1;
2985 }
2986
2987 if(type == kJetsRec || type == kJetsRecAcceptance){ // reconstructed jets
2988
2989 if(fBranchRecJets.Length()==0){
2990 Printf("%s:%d no rec jet branch specified", (char*)__FILE__,__LINE__);
2991 if(fDebug>1)fAOD->Print();
2992 return 0;
2993 }
2994
2995 TClonesArray *aodRecJets = 0;
2996 if(fBranchRecJets.Length()) aodRecJets = dynamic_cast<TClonesArray*>(fAODJets->FindListObject(fBranchRecJets.Data()));
2997 if(!aodRecJets) aodRecJets = dynamic_cast<TClonesArray*>(fAODJets->GetList()->FindObject(fBranchRecJets.Data()));
2998 if(fAODExtension&&!aodRecJets) aodRecJets = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fBranchRecJets.Data()));
2999
3000 if(!aodRecJets){
3001 if(fBranchRecJets.Length()) Printf("%s:%d no reconstructed jet array with name %s in AOD", (char*)__FILE__,__LINE__,fBranchRecJets.Data());
3002 if(fDebug>1)fAOD->Print();
3003 return 0;
3004 }
3005
3006 // Reorder jet pt and fill new temporary AliAODJet objects
3007 Int_t nRecJets = 0;
3008
3009 for(Int_t ij=0; ij<aodRecJets->GetEntries(); ++ij){
3010
3011 AliAODJet *tmp = dynamic_cast<AliAODJet*>(aodRecJets->At(ij));
3012 if(!tmp) continue;
3013
3014 if( tmp->Pt() < fJetPtCut ) continue;
3015 if( type == kJetsRecAcceptance &&
3016 ( tmp->Eta() < fJetEtaMin
3017 || tmp->Eta() > fJetEtaMax
3018 || tmp->Phi() < fJetPhiMin
3019 || tmp->Phi() > fJetPhiMax )) continue;
3020
3021
3022 list->Add(tmp);
3023 nRecJets++;
3024 }
3025
3026 list->Sort();
3027
3028 return nRecJets;
3029 }
3030 else if(type == kJetsKine || type == kJetsKineAcceptance){
3031
3032 // generated jets
3033 Int_t nGenJets = 0;
3034
3035 if(!fMCEvent){
3036 if(fDebug>1) Printf("%s:%d no mcEvent",(char*)__FILE__,__LINE__);
3037 return 0;
3038 }
3039
3040 AliGenEventHeader* genHeader = fMCEvent->GenEventHeader();
3041 AliGenPythiaEventHeader* pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
3042 AliGenHijingEventHeader* hijingGenHeader = 0x0;
3043
3044 if(!pythiaGenHeader){
3045 hijingGenHeader = dynamic_cast<AliGenHijingEventHeader*>(genHeader);
3046
3047 if(!hijingGenHeader){
3048 Printf("%s:%d no pythiaGenHeader or hijingGenHeader found", (char*)__FILE__,__LINE__);
3049 return 0;
3050 }else{
3051 TLorentzVector mom[4];
3052 AliAODJet* jet[4];
3053 hijingGenHeader->GetJets(mom[0], mom[1], mom[2], mom[3]);
3054
3055 for(Int_t i=0; i<2; ++i){
3056 if(!mom[i].Pt()) continue;
3057 jet[i] = new AliAODJet(mom[i]);
3058
3059 if( type == kJetsKineAcceptance &&
3060 ( jet[i]->Eta() < fJetEtaMin
3061 || jet[i]->Eta() > fJetEtaMax
3062 || jet[i]->Phi() < fJetPhiMin
3063 || jet[i]->Phi() > fJetPhiMax )) continue;
3064
3065 list->Add(jet[i]);
3066 nGenJets++;
3067 }
3068 list->Sort();
3069 return nGenJets;
3070 }
3071 }
3072
3073 // fetch the pythia generated jets
3074 for(int ip=0; ip<pythiaGenHeader->NTriggerJets(); ++ip){
3075
3076 Float_t p[4];
3077 AliAODJet *jet = new AliAODJet();
3078 pythiaGenHeader->TriggerJet(ip, p);
3079 jet->SetPxPyPzE(p[0], p[1], p[2], p[3]);
3080
3081 if( type == kJetsKineAcceptance &&
3082 ( jet->Eta() < fJetEtaMin
3083 || jet->Eta() > fJetEtaMax
3084 || jet->Phi() < fJetPhiMin
3085 || jet->Phi() > fJetPhiMax )) continue;
3086
3087 list->Add(jet);
3088 nGenJets++;
3089 }
3090 list->Sort();
3091 return nGenJets;
3092 }
3093 else if(type == kJetsGen || type == kJetsGenAcceptance ){
3094
3095 if(fBranchGenJets.Length()==0){
3096 if(fDebug>1) Printf("%s:%d no gen jet branch specified", (char*)__FILE__,__LINE__);
3097 return 0;
3098 }
3099
3100 TClonesArray *aodGenJets = 0;
3101 if(fBranchGenJets.Length()) aodGenJets = dynamic_cast<TClonesArray*>(fAODJets->FindListObject(fBranchGenJets.Data()));
3102 if(!aodGenJets) aodGenJets = dynamic_cast<TClonesArray*>(fAODJets->GetList()->FindObject(fBranchGenJets.Data()));
3103 if(fAODExtension&&!aodGenJets) aodGenJets = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fBranchGenJets.Data()));
3104
3105 if(!aodGenJets){
3106 if(fDebug>0){
3107 if(fBranchGenJets.Length()) Printf("%s:%d Generated jet branch %s not found",(char*)__FILE__,__LINE__,fBranchGenJets.Data());
3108 }
3109 if(fDebug>1)fAOD->Print();
3110 return 0;
3111 }
3112
3113 Int_t nGenJets = 0;
3114
3115 for(Int_t ig=0; ig<aodGenJets->GetEntries(); ++ig){
3116
3117 AliAODJet *tmp = dynamic_cast<AliAODJet*>(aodGenJets->At(ig));
3118 if(!tmp) continue;
3119
3120 if( tmp->Pt() < fJetPtCut ) continue;
3121 if( type == kJetsGenAcceptance &&
3122 ( tmp->Eta() < fJetEtaMin
3123 || tmp->Eta() > fJetEtaMax
3124 || tmp->Phi() < fJetPhiMin
3125 || tmp->Phi() > fJetPhiMax )) continue;
3126
3127 list->Add(tmp);
3128 nGenJets++;
3129 }
3130 list->Sort();
3131 return nGenJets;
3132 }
3133 else if(type == kJetsEmbedded){ // embedded jets
3134
3135 if(fBranchEmbeddedJets.Length()==0){
3136 Printf("%s:%d no embedded jet branch specified", (char*)__FILE__,__LINE__);
3137 if(fDebug>1)fAOD->Print();
3138 return 0;
3139 }
3140
3141 TClonesArray *aodEmbeddedJets = 0;
3142 if(fBranchEmbeddedJets.Length()) aodEmbeddedJets = dynamic_cast<TClonesArray*>(fAODJets->FindListObject(fBranchEmbeddedJets.Data()));
3143 if(!aodEmbeddedJets) aodEmbeddedJets = dynamic_cast<TClonesArray*>(fAODJets->GetList()->FindObject(fBranchEmbeddedJets.Data()));
3144 if(fAODExtension&&!aodEmbeddedJets) aodEmbeddedJets = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fBranchEmbeddedJets.Data()));
3145
3146 if(!aodEmbeddedJets){
3147 if(fBranchEmbeddedJets.Length()) Printf("%s:%d no reconstructed jet array with name %s in AOD", (char*)__FILE__,__LINE__,fBranchEmbeddedJets.Data());
3148 if(fDebug>1)fAOD->Print();
3149 return 0;
3150 }
3151
3152 // Reorder jet pt and fill new temporary AliAODJet objects
3153 Int_t nEmbeddedJets = 0;
3154
3155 for(Int_t ij=0; ij<aodEmbeddedJets->GetEntries(); ++ij){
3156
3157 AliAODJet *tmp = dynamic_cast<AliAODJet*>(aodEmbeddedJets->At(ij));
3158 if(!tmp) continue;
3159
3160 if( tmp->Pt() < fJetPtCut ) continue;
3161 if( tmp->Eta() < fJetEtaMin
3162 || tmp->Eta() > fJetEtaMax
3163 || tmp->Phi() < fJetPhiMin
3164 || tmp->Phi() > fJetPhiMax ) continue;
3165
3166 list->Add(tmp);
3167 nEmbeddedJets++;
3168 }
3169
3170 list->Sort();
3171
3172 return nEmbeddedJets;
3173 }
3174 else{
3175 if(fDebug>0)Printf("%s:%d no such type %d",(char*)__FILE__,__LINE__,type);
3176 return 0;
3177 }
3178}
3179
3180// ___________________________________________________________________________________
3181Int_t AliAnalysisTaskFragmentationFunction::GetListOfBckgJets(TList *list, Int_t type)
3182{
3183 // fill list of bgr clusters selected according to type
3184
3185 if(type == kJetsRec || type == kJetsRecAcceptance){ // reconstructed jets
3186
3187 if(fBranchRecBckgClusters.Length()==0){
3188 Printf("%s:%d no rec jet branch specified", (char*)__FILE__,__LINE__);
3189 if(fDebug>1)fAOD->Print();
3190 return 0;
3191 }
3192
3193 TClonesArray *aodRecJets = 0;
3194 if(fBranchRecBckgClusters.Length()) aodRecJets = dynamic_cast<TClonesArray*>(fAODJets->FindListObject(fBranchRecBckgClusters.Data()));
3195 if(!aodRecJets) aodRecJets = dynamic_cast<TClonesArray*>(fAODJets->GetList()->FindObject(fBranchRecBckgClusters.Data()));
3196 if(fAODExtension&&!aodRecJets) aodRecJets = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fBranchRecBckgClusters.Data()));
3197
3198 if(!aodRecJets){
3199 if(fBranchRecBckgClusters.Length()) Printf("%s:%d no reconstructed jet array with name %s in AOD", (char*)__FILE__,__LINE__,fBranchRecBckgClusters.Data());
3200 if(fDebug>1)fAOD->Print();
3201 return 0;
3202 }
3203
3204 // Reorder jet pt and fill new temporary AliAODJet objects
3205 Int_t nRecJets = 0;
3206
3207 for(Int_t ij=0; ij<aodRecJets->GetEntries(); ++ij){
3208
3209 AliAODJet *tmp = dynamic_cast<AliAODJet*>(aodRecJets->At(ij));
3210 if(!tmp) continue;
3211
3212 // if( tmp->Pt() < fJetPtCut ) continue; // no pt cut on bckg clusters !
3213 if( type == kJetsRecAcceptance &&
3214 ( tmp->Eta() < fJetEtaMin
3215 || tmp->Eta() > fJetEtaMax
3216 || tmp->Phi() < fJetPhiMin
3217 || tmp->Phi() > fJetPhiMax )) continue;
3218
3219 list->Add(tmp);
3220
3221 nRecJets++;
3222
3223 }
3224
3225 list->Sort();
3226
3227 return nRecJets;
3228 }
3229
3230 // /*
3231 // MC clusters still Under construction
3232 // */
3233
3234 return 0;
3235}
3236
3237// _________________________________________________________________________________________________________
3238void AliAnalysisTaskFragmentationFunction::SetProperties(THnSparse* h,const Int_t dim, const char** labels)
3239{
3240 // Set properties of THnSparse
3241
3242 for(Int_t i=0; i<dim; i++){
3243 h->GetAxis(i)->SetTitle(labels[i]);
3244 h->GetAxis(i)->SetTitleColor(1);
3245 }
3246}
3247
3248// __________________________________________________________________________________________
3249void AliAnalysisTaskFragmentationFunction::SetProperties(TH1* h,const char* x, const char* y)
3250{
3251 //Set properties of histos (x and y title)
3252
3253 h->SetXTitle(x);
3254 h->SetYTitle(y);
3255 h->GetXaxis()->SetTitleColor(1);
3256 h->GetYaxis()->SetTitleColor(1);
3257}
3258
3259// _________________________________________________________________________________________________________
3260void AliAnalysisTaskFragmentationFunction::SetProperties(TH1* h,const char* x, const char* y, const char* z)
3261{
3262 //Set properties of histos (x,y and z title)
3263
3264 h->SetXTitle(x);
3265 h->SetYTitle(y);
3266 h->SetZTitle(z);
3267 h->GetXaxis()->SetTitleColor(1);
3268 h->GetYaxis()->SetTitleColor(1);
3269 h->GetZaxis()->SetTitleColor(1);
3270}
3271
3272// ________________________________________________________________________________________________________________________________________________________
3273void AliAnalysisTaskFragmentationFunction::GetJetTracksPointing(TList* inputlist, TList* outputlist, const AliAODJet* jet,
3274 const Double_t radius, Double_t& sumPt, const Double_t minPtL, const Double_t maxPt, Bool_t& isBadPt)
3275{
3276 // fill list of tracks in cone around jet axis
3277
3278 sumPt = 0;
3279 Bool_t isBadMaxPt = kFALSE;
3280 Bool_t isBadMinPt = kTRUE;
3281
3282 Double_t jetMom[3];
3283 jet->PxPyPz(jetMom);
3284 TVector3 jet3mom(jetMom);
3285
3286 for (Int_t itrack=0; itrack<inputlist->GetSize(); itrack++){
3287
3288 AliVParticle* track = dynamic_cast<AliVParticle*>(inputlist->At(itrack));
3289 if(!track)continue;
3290 Double_t trackMom[3];
3291 track->PxPyPz(trackMom);
3292 TVector3 track3mom(trackMom);
3293
3294 Double_t dR = jet3mom.DeltaR(track3mom);
3295
3296 if(dR<radius){
3297
3298 outputlist->Add(track);
3299
3300 sumPt += track->Pt();
3301
3302 if(maxPt>0 && track->Pt()>maxPt) isBadMaxPt = kTRUE;
3303 if(minPtL>0 && track->Pt()>minPtL) isBadMinPt = kFALSE;
3304 }
3305 }
3306
3307 isBadPt = kFALSE;
3308 if(minPtL>0 && isBadMinPt) isBadPt = kTRUE;
3309 if(maxPt>0 && isBadMaxPt) isBadPt = kTRUE;
3310
3311 outputlist->Sort();
3312}
3313
3314// _________________________________________________________________________________________________________________________________________________________________
3315void AliAnalysisTaskFragmentationFunction::GetJetTracksTrackrefs(TList* list, const AliAODJet* jet, const Double_t minPtL, const Double_t maxPt, Bool_t& isBadPt)
3316{
3317 // list of jet tracks from trackrefs
3318
3319 Int_t nTracks = jet->GetRefTracks()->GetEntriesFast();
3320
3321 Bool_t isBadMaxPt = kFALSE;
3322 Bool_t isBadMinPt = kTRUE;
3323
3324 for(Int_t itrack=0; itrack<nTracks; itrack++) {
3325
3326 AliVParticle* track = dynamic_cast<AliVParticle*>(jet->GetRefTracks()->At(itrack));
3327 if(!track){
3328 AliError("expected ref track not found ");
3329 continue;
3330 }
3331
3332 if(track->Pt() < fTrackPtCut) continue; // track refs may contain low pt cut (bug in AliFastJetInput)
3333 if(maxPt>0 && track->Pt()>maxPt) isBadMaxPt = kTRUE;
3334 if(minPtL>0 && track->Pt()>minPtL) isBadMinPt = kFALSE;
3335
3336 list->Add(track);
3337 }
3338
3339 isBadPt = kFALSE;
3340 if(minPtL>0 && isBadMinPt) isBadPt = kTRUE;
3341 if(maxPt>0 && isBadMaxPt) isBadPt = kTRUE;
3342
3343 list->Sort();
3344}
3345
3346// _ ________________________________________________________________________________________________________________________________
3347void AliAnalysisTaskFragmentationFunction::AssociateGenRec(TList* tracksAODMCCharged,TList* tracksRec, TArrayI& indexAODTr,TArrayI& indexMCTr,
3348 TArrayS& isRefGen,TH2F* fh2PtRecVsGen)
3349{
3350 // associate generated and reconstructed tracks, fill TArrays of list indices
3351
3352 Int_t nTracksRec = tracksRec->GetSize();
3353 Int_t nTracksGen = tracksAODMCCharged->GetSize();
3354 TClonesArray *tca = dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName()));
3355
3356
3357 if(!nTracksGen) return;
3358 if(!tca) return;
3359
3360 // set size
3361 indexAODTr.Set(nTracksGen);
3362 indexMCTr.Set(nTracksRec);
3363 isRefGen.Set(nTracksGen);
3364
3365 indexAODTr.Reset(-1);
3366 indexMCTr.Reset(-1);
3367 isRefGen.Reset(0);
3368
3369 // loop over reconstructed tracks, get generated track
3370
3371 for(Int_t iRec=0; iRec<nTracksRec; iRec++){
3372
3373 AliAODTrack* rectrack = dynamic_cast<AliAODTrack*>(tracksRec->At(iRec));
3374 if(!rectrack)continue;
3375 Int_t label = TMath::Abs(rectrack->GetLabel());
3376
3377 // find MC track in our list
3378 AliAODMCParticle* gentrack = dynamic_cast<AliAODMCParticle*> (tca->At(label));
3379
3380 Int_t listIndex = -1;
3381 if(gentrack) listIndex = tracksAODMCCharged->IndexOf(gentrack);
3382
3383 if(listIndex>=0){
3384
3385 indexAODTr[listIndex] = iRec;
3386 indexMCTr[iRec] = listIndex;
3387 }
3388 }
3389
3390
3391 // define reference sample of primaries/secondaries (for reconstruction efficiency / contamination)
3392
3393 for(Int_t iGen=0; iGen<nTracksGen; iGen++){
3394
3395 AliAODMCParticle* gentrack = dynamic_cast<AliAODMCParticle*> (tracksAODMCCharged->At(iGen));
3396 if(!gentrack)continue;
3397 Int_t pdg = gentrack->GetPdgCode();
3398
3399 // 211 - pi, 2212 - proton, 321 - Kaon, 11 - electron, 13 - muon
3400 if(TMath::Abs(pdg) == 211 || TMath::Abs(pdg) == 2212 || TMath::Abs(pdg) == 321 ||
3401 TMath::Abs(pdg) == 11 || TMath::Abs(pdg) == 13){
3402
3403 isRefGen[iGen] = kTRUE;
3404
3405 Int_t iRec = indexAODTr[iGen]; // can be -1 if no good reconstructed track
3406
3407 if(iRec>=0){
3408 Float_t genPt = gentrack->Pt();
3409 AliAODTrack* vt = dynamic_cast<AliAODTrack*>(tracksRec->At(iRec));
3410 if(vt){
3411 Float_t recPt = vt->Pt();
3412 fh2PtRecVsGen->Fill(genPt,recPt);
3413 }
3414 }
3415 }
3416 }
3417}
3418
3419// _____________________________________________________________________________________________________________________________________________
3420void AliAnalysisTaskFragmentationFunction::FillSingleTrackHistosRecGen(AliFragFuncQATrackHistos* trackQAGen, AliFragFuncQATrackHistos* trackQARec, TList* tracksGen,
3421 const TArrayI& indexAODTr, const TArrayS& isRefGen, const Bool_t scaleStrangeness){
3422
3423 // fill QA for single track reconstruction efficiency
3424
3425 Int_t nTracksGen = tracksGen->GetSize();
3426
3427 if(!nTracksGen) return;
3428
3429 for(Int_t iGen=0; iGen<nTracksGen; iGen++){
3430
3431 if(isRefGen[iGen] != 1) continue; // select primaries
3432
3433 AliAODMCParticle* gentrack = dynamic_cast<AliAODMCParticle*> (tracksGen->At(iGen));
3434 if(!gentrack) continue;
3435 Double_t ptGen = gentrack->Pt();
3436 Double_t etaGen = gentrack->Eta();
3437 Double_t phiGen = TVector2::Phi_0_2pi(gentrack->Phi());
3438
3439 // apply same acc & pt cuts as for FF
3440
3441 if(etaGen < fTrackEtaMin || etaGen > fTrackEtaMax) continue;
3442 if(phiGen < fTrackPhiMin || phiGen > fTrackPhiMax) continue;
3443 if(ptGen < fTrackPtCut) continue;
3444
3445 if(trackQAGen) trackQAGen->FillTrackQA(etaGen, phiGen, ptGen);
3446
3447 Int_t iRec = indexAODTr[iGen]; // can be -1 if no good reconstructed track
3448
3449 if(iRec>=0 && trackQARec){
3450 if(scaleStrangeness){
2b3fb225 3451 //Double_t weight = GetMCStrangenessFactor(ptGen);
3452 Double_t weight = GetMCStrangenessFactorCMS(gentrack);
ffa175ff 3453 trackQARec->FillTrackQA(etaGen, phiGen, ptGen, kFALSE, 0, kTRUE, weight);
3454 }
3455 else trackQARec->FillTrackQA(etaGen, phiGen, ptGen);
3456 }
3457 }
3458}
3459
3460// ______________________________________________________________________________________________________________________________________________________
3461
3462void AliAnalysisTaskFragmentationFunction::FillJetTrackHistosRec(AliFragFuncHistos* ffhistRec, AliAODJet* jet,
3463 TList* jetTrackList, const TList* tracksGen, const TList* tracksRec, const TArrayI& indexAODTr,
3464 const TArrayS& isRefGen, TList* jetTrackListTR, const Bool_t scaleStrangeness,
3465 Bool_t fillJS, TProfile* hProNtracksLeadingJet, TProfile** hProDelRPtSum, TProfile* hProDelR80pcPt)
3466{
3467 // fill objects for jet track reconstruction efficiency or secondaries contamination
3468 // arguments histGen/histRec can be of different type: AliFragFuncHistos*, AliFragFuncIntraJetHistos*, ...
3469 // jetTrackListTR pointer: track refs if not NULL
3470
3471
3472 // ensure proper normalization, even for secondaries
3473 Double_t jetPtRec = jet->Pt();
3474 ffhistRec->FillFF(-1, jetPtRec, kTRUE);
3475
3476 Int_t nTracksJet = jetTrackList->GetSize(); // list with AODMC tracks
3477 if(nTracksJet == 0) return;
3478
3479 TList* listRecTracks = new TList();
3480 listRecTracks->Clear();
3481
3482 for(Int_t iTr=0; iTr<nTracksJet; iTr++){ // jet tracks loop
3483
3484 AliAODMCParticle* gentrack = dynamic_cast<AliAODMCParticle*> (jetTrackList->At(iTr));
3485 if(!gentrack)continue;
3486 // find jet track in gen tracks list
3487 Int_t iGen = tracksGen->IndexOf(gentrack);
3488
3489 if(iGen<0){
3490 if(fDebug>0) Printf("%s:%d gen jet track not found ",(char*)__FILE__,__LINE__);
3491 continue;
3492 }
3493
3494 if(isRefGen[iGen] != 1) continue; // select primaries
3495
3496 Double_t ptGen = gentrack->Pt();
3497 Double_t etaGen = gentrack->Eta();
3498 Double_t phiGen = TVector2::Phi_0_2pi(gentrack->Phi());
3499
3500 // gen level acc & pt cuts - skip in case of track refs
3501 if(!jetTrackListTR && (etaGen < fTrackEtaMin || etaGen > fTrackEtaMax)) continue;
3502 if(!jetTrackListTR && (phiGen < fTrackPhiMin || phiGen > fTrackPhiMax)) continue;
3503 if(!jetTrackListTR && ptGen < fTrackPtCut) continue;
3504
3505
3506 Double_t ptRec = -1;
3507
3508 Int_t iRec = indexAODTr[iGen]; // can be -1 if no good reconstructed track
3509 Bool_t isRec = (iRec>=0) ? kTRUE : kFALSE;
3510
3511 Bool_t isJetTrack = kFALSE;
3512 if(!jetTrackListTR) isJetTrack = kTRUE; // skip trackRefs check for tracks in ideal cone
3513
3514 if(isRec){
3515
3516 AliAODTrack* rectrack = dynamic_cast<AliAODTrack*> (tracksRec->At(iRec));
3517 if(!rectrack) continue;
3518
3519 ptRec = rectrack->Pt();
3520
3521 if(jetTrackListTR){
3522 Int_t iRecTR = jetTrackListTR->IndexOf(rectrack);
3523 if(iRecTR >=0 ) isJetTrack = kTRUE; // rec tracks assigned to jet
3524 }
3525
3526 if(isJetTrack){
3527
3528 Double_t trackPt = ptRec;
3529 Bool_t incrementJetPt = kFALSE;
3530
3531 if(scaleStrangeness){
2b3fb225 3532 //Double_t weight = GetMCStrangenessFactor(ptGen);
3533 Double_t weight = GetMCStrangenessFactorCMS(gentrack);
3534
ffa175ff 3535 ffhistRec->FillFF( trackPt, jetPtRec, incrementJetPt, 0, kTRUE, weight );
3536 }
3537 else{
3538 ffhistRec->FillFF( trackPt, jetPtRec, incrementJetPt );
3539 }
3540
3541 listRecTracks->Add(rectrack);
3542
3543 }
3544 }
3545 }
3546
3547
3548 if(fillJS) FillJetShape(jet,listRecTracks,hProNtracksLeadingJet, hProDelRPtSum, hProDelR80pcPt,0,0,scaleStrangeness);
3549
3550 delete listRecTracks;
3551
3552}
3553
3554// _____________________________________________________________________________________________________________________________________________________________________
3555void AliAnalysisTaskFragmentationFunction::GetTracksTiltedwrpJetAxis(Float_t alpha, TList* inputlist, TList* outputlist, const AliAODJet* jet, Double_t radius,Double_t& sumPt)
3556{
3557 // List of tracks in cone perpendicular to the jet azimuthal direction
3558
3559 Double_t jetMom[3];
3560 jet->PxPyPz(jetMom);
3561
3562 TVector3 jet3mom(jetMom);
3563 // Rotate phi and keep eta unchanged
3564 Double_t etaTilted = jet3mom.Eta();
3565 Double_t phiTilted = TVector2::Phi_0_2pi(jet3mom.Phi()) + alpha;
3566 if(phiTilted > 2*TMath::Pi()) phiTilted = phiTilted - 2*TMath::Pi();
3567
3568 for (Int_t itrack=0; itrack<inputlist->GetSize(); itrack++){
3569
3570 // embedded tracks
3571 if( fUseExtraTracksBgr != 1){
3572 if(AliAODTrack* trackAOD = dynamic_cast<AliAODTrack*> (inputlist->At(itrack))){
3573 if(fUseExtraTracksBgr == 0 && ((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
3574 if(fUseExtraTracksBgr == -1 && !((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
3575 }
3576 }
3577
3578 AliVParticle* track = dynamic_cast<AliVParticle*>(inputlist->At(itrack));
3579 if(!track)continue;
3580 Double_t trackMom[3];
3581 track->PxPyPz(trackMom);
3582 TVector3 track3mom(trackMom);
3583
3584 Double_t deta = track3mom.Eta() - etaTilted;
3585 Double_t dphi = TMath::Abs(track3mom.Phi() - phiTilted);
3586 if (dphi > TMath::Pi()) dphi = 2. * TMath::Pi() - dphi;
3587 Double_t dR = TMath::Sqrt(deta * deta + dphi * dphi);
3588
3589
3590 if(dR<=radius){
3591 outputlist->Add(track);
3592 sumPt += track->Pt();
3593 }
3594 }
3595
3596}
3597
3598// ________________________________________________________________________________________________________________________________________________________
3599void AliAnalysisTaskFragmentationFunction::GetTracksTiltedwrpJetAxisWindow(Float_t alpha, TList* inputlist, TList* outputlist, const AliAODJet* jet, Double_t radius,Double_t& sumPt,Double_t &normFactor)
3600{
3601 // List of tracks in cone perpendicular to the jet azimuthal direction
3602
3603 Double_t jetMom[3];
3604 jet->PxPyPz(jetMom);
3605
3606 TVector3 jet3mom(jetMom);
3607 // Rotate phi and keep eta unchanged
3608 Double_t etaTilted = jet3mom.Eta();
3609 Double_t phiTilted = TVector2::Phi_0_2pi(jet3mom.Phi()) + alpha;
3610 if(phiTilted > 2*TMath::Pi()) phiTilted = phiTilted - 2*TMath::Pi();
3611
3612 for (Int_t itrack=0; itrack<inputlist->GetSize(); itrack++)
3613 {
3614
3615 // embedded tracks
3616 if( fUseExtraTracksBgr != 1){
3617 if(AliAODTrack* trackAOD = dynamic_cast<AliAODTrack*> (inputlist->At(itrack))){
3618 if(fUseExtraTracksBgr == 0 && ((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
3619 if(fUseExtraTracksBgr == -1 && !((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
3620 }
3621 }
3622
3623 AliVParticle* track = dynamic_cast<AliVParticle*>(inputlist->At(itrack));
3624 if(!track)continue;
3625 Float_t trackEta = track->Eta();
3626 Float_t trackPhi = track->Phi();
3627
3628 if( ( phiTilted-radius >= 0 ) && ( phiTilted+radius <= 2*TMath::Pi()))
3629 {
3630 if((trackPhi<=phiTilted+radius) &&
3631 (trackPhi>=phiTilted-radius) &&
3632 (trackEta<=fTrackEtaMax)&&(trackEta>=fTrackEtaMin)) // 0.9 and - 0.9
3633 {
3634 outputlist->Add(track);
3635 sumPt += track->Pt();
3636 }
3637 }
3638 else if( phiTilted-radius < 0 )
3639 {
3640 if((( trackPhi < phiTilted+radius ) ||
3641 ( trackPhi > 2*TMath::Pi()-(radius-phiTilted) )) &&
3642 (( trackEta <= fTrackEtaMax ) && ( trackEta >= fTrackEtaMin )))
3643 {
3644 outputlist->Add(track);
3645 sumPt += track->Pt();
3646 }
3647 }
3648 else if( phiTilted+radius > 2*TMath::Pi() )
3649 {
3650 if((( trackPhi > phiTilted-radius ) ||
3651 ( trackPhi < phiTilted+radius-2*TMath::Pi() )) &&
3652 (( trackEta <= fTrackEtaMax ) && ( trackEta >= fTrackEtaMin )))
3653 {
3654 outputlist->Add(track);
3655 sumPt += track->Pt();
3656 }
3657 }
3658 }
3659
3660 // Jet area - Temporarily added should be modified with the proper jet area value
3661 Float_t areaJet = CalcJetArea(etaTilted,radius);
3662 Float_t areaTilted = 2*radius*(fTrackEtaMax-fTrackEtaMin);
3663
3664 normFactor = (Float_t) 1. / (areaJet / areaTilted);
3665
3666}
3667
3668
3669// ________________________________________________________________________________________________________________________________________________________
3670void AliAnalysisTaskFragmentationFunction::GetTracksOutOfNJets(Int_t nCases, TList* inputlist, TList* outputlist, TList* jetlist, Double_t& sumPt)
3671{
3672 // List of tracks outside cone around N jet axis
3673 // Particles taken randomly
3674
3675 sumPt = 0;
3676 // Int_t nj = jetlist->GetSize();
3677 Float_t rc = TMath::Abs(GetFFRadius());
3678 Float_t rcl = GetFFBckgRadius();
3679
3680 // Estimate jet and background areas
3681 Float_t* areaJet = new Float_t[nCases];
3682 memset(areaJet, 0, sizeof(Float_t) * nCases);
3683 Float_t* areaJetLarge = new Float_t[nCases];
3684 memset(areaJetLarge, 0, sizeof(Float_t) * nCases);
3685 Float_t areaFull = (fTrackEtaMax-fTrackEtaMin)*(fTrackPhiMax-fTrackPhiMin);
3686 Float_t areaOut = areaFull;
3687
3688 //estimate jets and background areas
3689 Int_t nOut = 0;
3690 Int_t ijet = 0;
3691 TList* templist = new TList();
3692 TClonesArray *vect3Jet = new TClonesArray("TVector3",nCases);
3693
3694 for(Int_t ij=0; ij<nCases; ++ij)
3695 {
3696 // Get jet information
3697 AliAODJet* jet = dynamic_cast<AliAODJet*>(jetlist->At(ij));
3698 if(!jet)continue;
3699 TVector3 jet3mom;
3700 jet3mom.SetPtEtaPhi(jet->Pt(),jet->Eta(),jet->Phi());
3701 new((*vect3Jet)[ijet]) TVector3((TVector3)jet3mom);
3702 Float_t etaJet = (Float_t)((TVector3*) vect3Jet->At(ij))->Eta();
3703
3704 // Jet area
3705 areaJet[ij] = CalcJetArea(etaJet,rc);
3706
3707 // Area jet larger angle
3708 areaJetLarge[ij] = CalcJetArea(etaJet,rcl);
3709
3710 // Outside jet area
3711 areaOut = areaOut - areaJetLarge[ij];
3712 ijet++;
3713 }
3714
3715 // List of all tracks outside jet areas
3716 for (Int_t itrack=0; itrack<inputlist->GetSize(); itrack++){
3717
3718 // embedded tracks
3719 if( fUseExtraTracksBgr != 1){
3720 if(AliAODTrack* trackAOD = dynamic_cast<AliAODTrack*> (inputlist->At(itrack))){
3721 if(fUseExtraTracksBgr == 0 && ((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
3722 if(fUseExtraTracksBgr == -1 && !((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
3723 }
3724 }
3725
3726 AliVParticle* track = dynamic_cast<AliVParticle*>(inputlist->At(itrack));
3727
3728 if(!track)continue;
3729 Double_t trackMom[3];
3730 track->PxPyPz(trackMom);
3731 TVector3 track3mom(trackMom);
3732
3733 Double_t *dR = new Double_t[nCases];
3734 for(Int_t ij=0; ij<nCases; ij++)
3735 dR[ij] = (Double_t)((TVector3*) vect3Jet->At(ij))->DeltaR(track3mom);
3736
3737 if((nCases==1 && (dR[0]>rcl)) ||
3738 (nCases==2 && (dR[0]>rcl && dR[1]>rcl)) ||
3739 (nCases==3 && (dR[0]>rcl && dR[1]>rcl && dR[2]>rcl)))
3740 {
3741 templist->Add(track);
3742 nOut++;
3743 }
3744 delete [] dR;
3745 }
3746
3747 // Take tracks randomly
3748 Int_t nScaled = (Int_t) (nOut * areaJet[0] / areaOut + 0.5);
3749 TArrayI* ar = new TArrayI(nOut);
3750
3751 for(Int_t init=0; init<nOut; init++)
3752 (*ar)[init] = init;
3753
3754 Int_t *randIndex = new Int_t[nScaled];
3755 for(Int_t init2=0; init2<nScaled; init2++)
3756 randIndex[init2] = -1;
3757
3758 // Select nScaled different random numbers in nOut
3759 for(Int_t i=0; i<nScaled; i++)
3760 {
3761 Int_t* tmpArr = new Int_t[nOut-i];
3762 Int_t temp = fRandom->Integer(nOut-i);
3763 for(Int_t ind = 0; ind< ar->GetSize()-1; ind++)
3764 {
3765 if(ind<temp) tmpArr[ind] = (*ar)[ind];
3766 else tmpArr[ind] = (*ar)[ind+1];
3767 }
3768 randIndex[i] = (*ar)[temp];
3769
3770 ar->Set(nOut-i-1,tmpArr);
3771
3772 delete [] tmpArr;
3773
3774 }
3775
3776 for(Int_t ipart=0; ipart<nScaled; ipart++)
3777 {
3778 AliVParticle* track = (AliVParticle*)(templist->At(randIndex[ipart]));
3779 outputlist->Add(track);
3780 sumPt += track->Pt();
3781 }
3782
3783 outputlist->Sort();
3784
3785 delete vect3Jet;
3786 delete templist;
3787 delete [] areaJetLarge;
3788 delete [] areaJet;
3789 delete ar;
3790 delete [] randIndex;
3791
3792}
3793
3794// ________________________________________________________________________________________________________________________________________________________
3795void AliAnalysisTaskFragmentationFunction::GetTracksOutOfNJetsStat(Int_t nCases, TList* inputlist, TList* outputlist, TList* jetlist, Double_t& sumPt, Double_t &normFactor)
3796{
3797 // List of tracks outside cone around N jet axis
3798 // All particles taken + final scaling factor
3799
3800 sumPt = 0;
3801 Float_t rc = TMath::Abs(GetFFRadius());
3802 Float_t rcl = GetFFBckgRadius();
3803
3804 // Estimate jet and background areas
3805 Float_t* areaJet = new Float_t[nCases];
3806 memset(areaJet, 0, sizeof(Float_t) * nCases);
3807 Float_t* areaJetLarge = new Float_t[nCases];
3808 memset(areaJetLarge, 0, sizeof(Float_t) * nCases);
3809 Float_t areaFull = (fTrackEtaMax-fTrackEtaMin)*(fTrackPhiMax-fTrackPhiMin);
3810 Float_t areaOut = areaFull;
3811
3812 //estimate jets and background areas
3813 Int_t nOut = 0;
3814 Int_t ijet = 0;
3815 TClonesArray *vect3Jet = new TClonesArray("TVector3",nCases);
3816
3817 for(Int_t ij=0; ij<nCases; ++ij)
3818 {
3819 // Get jet information
3820 AliAODJet* jet = dynamic_cast<AliAODJet*>(jetlist->At(ij));
3821 if(!jet)continue;
3822 TVector3 jet3mom;
3823 jet3mom.SetPtEtaPhi(jet->Pt(),jet->Eta(),jet->Phi());
3824 new((*vect3Jet)[ijet]) TVector3((TVector3)jet3mom);
3825 Float_t etaJet = (Float_t)((TVector3*) vect3Jet->At(ij))->Eta();
3826
3827 // Jet area
3828 areaJet[ij] = CalcJetArea(etaJet,rc);
3829
3830 // Area jet larger angle
3831 areaJetLarge[ij] = CalcJetArea(etaJet,rcl);
3832
3833 // Outside jets area
3834 areaOut = areaOut - areaJetLarge[ij];
3835 ijet++;
3836 }
3837
3838 for (Int_t itrack=0; itrack<inputlist->GetSize(); itrack++){
3839
3840 // embedded tracks
3841 if( fUseExtraTracksBgr != 1){
3842 if(AliAODTrack* trackAOD = dynamic_cast<AliAODTrack*> (inputlist->At(itrack))){
3843 if(fUseExtraTracksBgr == 0 && ((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
3844 if(fUseExtraTracksBgr == -1 && !((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
3845 }
3846 }
3847
3848 AliVParticle* track = dynamic_cast<AliVParticle*>(inputlist->At(itrack));
3849 if(!track)continue;
3850 Double_t trackMom[3];
3851 track->PxPyPz(trackMom);
3852 TVector3 track3mom(trackMom);
3853
3854 Double_t *dR = new Double_t[nCases];
3855 for(Int_t ij=0; ij<nCases; ij++)
3856 dR[ij] = (Double_t)((TVector3*) vect3Jet->At(ij))->DeltaR(track3mom);
3857
3858 if((nCases==0) ||
3859 (nCases==1 && (dR[0]>rcl)) ||
3860 (nCases==2 && (dR[0]>rcl && dR[1]>rcl)) ||
3861 (nCases==3 && (dR[0]>rcl && dR[1]>rcl && dR[2]>rcl)))
3862 {
3863 outputlist->Add(track);
3864 sumPt += track->Pt();
3865 nOut++;
3866 }
3867 delete [] dR;
3868 }
3869
3870 if(nCases==0) areaJet[0] = TMath::Pi()*rc*rc;
3871 normFactor = (Float_t) 1./(areaJet[0] / areaOut);
3872
3873 outputlist->Sort();
3874
3875 delete vect3Jet;
3876 delete [] areaJetLarge;
3877 delete [] areaJet;
3878
3879}
3880
3881// ______________________________________________________________________________________________________________________________________________________
3882Float_t AliAnalysisTaskFragmentationFunction::CalcJetArea(const Float_t etaJet, const Float_t rc) const
3883{
3884 // calculate area of jet with eta etaJet and radius rc
3885
3886 Float_t detamax = etaJet + rc;
3887 Float_t detamin = etaJet - rc;
3888 Float_t accmax = 0.0; Float_t accmin = 0.0;
3889 if(detamax > fTrackEtaMax){ // sector outside etamax
3890 Float_t h = fTrackEtaMax - etaJet;
3891 accmax = rc*rc*TMath::ACos(h/rc) - h*TMath::Sqrt(rc*rc - h*h);
3892 }
3893 if(detamin < fTrackEtaMin){ // sector outside etamin
3894 Float_t h = fTrackEtaMax + etaJet;
3895 accmin = rc*rc*TMath::ACos(h/rc) - h*TMath::Sqrt(rc*rc - h*h);
3896 }
3897 Float_t areaJet = rc*rc*TMath::Pi() - accmax - accmin;
3898
3899 return areaJet;
3900
3901}
3902
3903// ___________________________________________________________________________________________________________________________
3904void AliAnalysisTaskFragmentationFunction::GetClusterTracksOutOf1Jet(AliAODJet* jet, TList* outputlist, Double_t &normFactor)
3905{
3906 // fill tracks from bckgCluster branch in list,
3907 // for all clusters outside jet cone
3908 // sum up total area of clusters
3909
3910 Double_t rc = GetFFRadius();
3911 Double_t rcl = GetFFBckgRadius();
3912
3913 Double_t areaTotal = 0;
3914 Double_t sumPtTotal = 0;
3915
3916 for(Int_t ij=0; ij<fBckgJetsRec->GetEntries(); ++ij){
3917
3918 AliAODJet* bgrCluster = (AliAODJet*)(fBckgJetsRec->At(ij)); // not 'recCuts': use all clusters in full eta range
3919
3920 Double_t dR = jet->DeltaR(bgrCluster);
3921
3922 if(dR<rcl) continue;
3923
3924 Double_t clusterPt = bgrCluster->Pt();
3925 Double_t area = bgrCluster->EffectiveAreaCharged();
3926 areaTotal += area;
3927 sumPtTotal += clusterPt;
3928
3929 Int_t nTracksJet = bgrCluster->GetRefTracks()->GetEntries();
3930
3931 for(Int_t it = 0; it<nTracksJet; it++){
3932
3933 // embedded tracks - note: using ref tracks here, fBranchRecBckgClusters has to be consistent
3934 if( fUseExtraTracksBgr != 1){
3935 if(AliAODTrack* trackAOD = dynamic_cast<AliAODTrack*> (bgrCluster->GetTrack(it))){
3936 if(fUseExtraTracksBgr == 0 && ((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
3937 if(fUseExtraTracksBgr == -1 && !((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
3938 }
3939 }
3940
3941 AliVParticle* track = dynamic_cast<AliVParticle*>(bgrCluster->GetTrack(it));
3942 if(!track) continue;
3943
3944 Float_t trackPt = track->Pt();
3945 Float_t trackEta = track->Eta();
3946 Float_t trackPhi = TVector2::Phi_0_2pi(track->Phi());
3947
3948 if(trackEta < fTrackEtaMin || trackEta > fTrackEtaMax) continue;
3949 if(trackPhi < fTrackPhiMin || trackPhi > fTrackPhiMax) continue;
3950 if(trackPt < fTrackPtCut) continue;
3951
3952 outputlist->Add(track);
3953 }
3954 }
3955
3956 Double_t areaJet = TMath::Pi()*rc*rc;
3957 if(areaTotal) normFactor = (Float_t) 1./(areaJet / areaTotal);
3958
3959 outputlist->Sort();
3960}
3961
3962// _______________________________________________________________________________________________________________________
3963void AliAnalysisTaskFragmentationFunction::GetClusterTracksMedian(TList* outputlist, Double_t &normFactor)
3964{
3965 // fill tracks from bckgCluster branch,
3966 // using cluster with median density (odd number of clusters)
3967 // or picking randomly one of the two closest to median (even number)
3968
3969 normFactor = 0;
3970
3971 Int_t nBckgClusters = fBckgJetsRec->GetEntries(); // not 'recCuts': use all clusters in full eta range
3972
3973 if(nBckgClusters<3) return; // need at least 3 clusters (skipping 2 highest)
3974
3975 Double_t* bgrDensity = new Double_t[nBckgClusters];
3976 Int_t* indices = new Int_t[nBckgClusters];
3977
3978 for(Int_t ij=0; ij<nBckgClusters; ++ij){
3979
3980 AliAODJet* bgrCluster = (AliAODJet*)(fBckgJetsRec->At(ij));
3981 Double_t clusterPt = bgrCluster->Pt();
3982 Double_t area = bgrCluster->EffectiveAreaCharged();
3983
3984 Double_t density = 0;
3985 if(area>0) density = clusterPt/area;
3986
3987 bgrDensity[ij] = density;
3988 indices[ij] = ij;
3989 }
3990
3991 TMath::Sort(nBckgClusters, bgrDensity, indices);
3992
3993 // get median cluster
3994
3995 AliAODJet* medianCluster = 0;
3996 Double_t medianDensity = 0;
3997
3998 if(TMath::Odd(nBckgClusters)){
3999
4000 Int_t medianIndex = indices[(Int_t) (0.5*(nBckgClusters-1))];
4001 medianCluster = (AliAODJet*)(fBckgJetsRec->At(medianIndex));
4002
4003 Double_t clusterPt = medianCluster->Pt();
4004 Double_t area = medianCluster->EffectiveAreaCharged();
4005
4006 if(area>0) medianDensity = clusterPt/area;
4007 }
4008 else{
4009
4010 Int_t medianIndex1 = indices[(Int_t) (0.5*nBckgClusters-1)];
4011 Int_t medianIndex2 = indices[(Int_t) (0.5*nBckgClusters)];
4012
4013 AliAODJet* medianCluster1 = (AliAODJet*)(fBckgJetsRec->At(medianIndex1));
4014 AliAODJet* medianCluster2 = (AliAODJet*)(fBckgJetsRec->At(medianIndex2));
4015
4016 Double_t density1 = 0;
4017 Double_t clusterPt1 = medianCluster1->Pt();
4018 Double_t area1 = medianCluster1->EffectiveAreaCharged();
4019 if(area1>0) density1 = clusterPt1/area1;
4020
4021 Double_t density2 = 0;
4022 Double_t clusterPt2 = medianCluster2->Pt();
4023 Double_t area2 = medianCluster2->EffectiveAreaCharged();
4024 if(area2>0) density2 = clusterPt2/area2;
4025
4026 medianDensity = 0.5*(density1+density2);
4027
4028 medianCluster = ( (fRandom->Rndm()>0.5) ? medianCluster1 : medianCluster2 ); // select one randomly to avoid adding areas
4029 }
4030
4031 Int_t nTracksJet = medianCluster->GetRefTracks()->GetEntries();
4032
4033 for(Int_t it = 0; it<nTracksJet; it++){
4034
4035 // embedded tracks - note: using ref tracks here, fBranchRecBckgClusters has to be consistent
4036 if( fUseExtraTracksBgr != 1){
4037 if(AliAODTrack* trackAOD = dynamic_cast<AliAODTrack*> (medianCluster->GetTrack(it))){
4038 if(fUseExtraTracksBgr == 0 && ((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
4039 if(fUseExtraTracksBgr == -1 && !((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
4040 }
4041 }
4042
4043 AliVParticle* track = dynamic_cast<AliVParticle*>(medianCluster->GetTrack(it));
4044 if(!track) continue;
4045
4046 Float_t trackPt = track->Pt();
4047 Float_t trackEta = track->Eta();
4048 Float_t trackPhi = TVector2::Phi_0_2pi(track->Phi());
4049
4050 if(trackEta < fTrackEtaMin || trackEta > fTrackEtaMax) continue;
4051 if(trackPhi < fTrackPhiMin || trackPhi > fTrackPhiMax) continue;
4052 if(trackPt < fTrackPtCut) continue;
4053
4054 outputlist->Add(track);
4055 }
4056
4057 Double_t areaMedian = medianCluster->EffectiveAreaCharged();
4058 Double_t areaJet = TMath::Pi()*GetFFRadius()*GetFFRadius();
4059
4060 if(areaMedian) normFactor = (Float_t) 1./(areaJet / areaMedian);
4061
4062 outputlist->Sort();
4063
4064 delete[] bgrDensity;
4065 delete[] indices;
4066}
4067
4068// ______________________________________________________________________________________________________________________________________________________
4069void AliAnalysisTaskFragmentationFunction::FillBckgHistos(Int_t type, TList* inputtracklist, TList* inputjetlist, AliAODJet* jet,
4070 AliFragFuncHistos* ffbckghistocuts, AliFragFuncQATrackHistos* qabckghistocuts, TH1F* fh1Mult){
4071
4072 // List of tracks outside jets for background study
4073 TList* tracklistout2jets = new TList();
4074 TList* tracklistout3jets = new TList();
4075 TList* tracklistout2jetsStat = new TList();
4076 TList* tracklistout3jetsStat = new TList();
4077 Double_t sumPtOut2Jets = 0.;
4078 Double_t sumPtOut3Jets = 0.;
4079 Double_t sumPtOut2JetsStat = 0.;
4080 Double_t sumPtOut3JetsStat = 0.;
4081 Double_t normFactor2Jets = 0.;
4082 Double_t normFactor3Jets = 0.;
4083
4084 Int_t nRecJetsCuts = inputjetlist->GetEntries();
4085
4086 if(nRecJetsCuts>1) {
4087 GetTracksOutOfNJets(2,inputtracklist, tracklistout2jets, inputjetlist, sumPtOut2Jets);
4088 GetTracksOutOfNJetsStat(2,inputtracklist, tracklistout2jetsStat, inputjetlist,sumPtOut2JetsStat, normFactor2Jets);
4089
4090 }
4091 if(nRecJetsCuts>2) {
4092 GetTracksOutOfNJets(3,inputtracklist, tracklistout3jets, inputjetlist, sumPtOut3Jets);
4093 GetTracksOutOfNJetsStat(3,inputtracklist, tracklistout3jetsStat, inputjetlist, sumPtOut3JetsStat, normFactor3Jets);
4094 }
4095
4096 if(type==kBckgOutLJ || type==kBckgOutAJ)
4097 {
4098 TList* tracklistoutleading = new TList();
4099 Double_t sumPtOutLeading = 0.;
4100 GetTracksOutOfNJets(1,inputtracklist, tracklistoutleading, inputjetlist, sumPtOutLeading);
4101 if(type==kBckgOutLJ && fh1Mult) fh1Mult->Fill(tracklistoutleading->GetSize());
4102
4103 for(Int_t it=0; it<tracklistoutleading->GetSize(); ++it){
4104
4105 AliVParticle* trackVP = (AliVParticle*)(tracklistoutleading->At(it));
4106 if(!trackVP) continue;
4107 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
4108
4109 Float_t jetPt = jet->Pt();
4110 Float_t trackPt = trackV->Pt();
4111
4112 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
4113
4114 if(type==kBckgOutLJ)
4115 {
4116 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt);
4117
4118 // Fill track QA for background
4119 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt);
4120 }
4121
4122 // All cases included
4123 if(nRecJetsCuts==1 && type==kBckgOutAJ)
4124 {
4125 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt );
4126 }
4127 delete trackV;
4128 }
4129 // Increment jet pt with one entry in case #tracks outside jets = 0
4130 if(tracklistoutleading->GetSize()==0) {
4131 Float_t jetPt = jet->Pt();
4132 Bool_t incrementJetPt = kTRUE;
4133 if(type==kBckgOutLJ)
4134 {
4135 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt );
4136 }
4137 // All cases included
4138 if(nRecJetsCuts==1 && type==kBckgOutAJ)
4139 {
4140 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt );
4141 }
4142 }
4143 delete tracklistoutleading;
4144 }
4145 if(type==kBckgOutLJStat || type==kBckgOutAJStat)
4146 {
4147 TList* tracklistoutleadingStat = new TList();
4148 Double_t sumPtOutLeadingStat = 0.;
4149 Double_t normFactorLeading = 0.;
4150
4151 GetTracksOutOfNJetsStat(1,inputtracklist, tracklistoutleadingStat, inputjetlist, sumPtOutLeadingStat, normFactorLeading);
4152 if(type==kBckgOutLJStat && fh1Mult) fh1Mult->Fill(tracklistoutleadingStat->GetSize());
4153
4154 for(Int_t it=0; it<tracklistoutleadingStat->GetSize(); ++it){
4155
4156 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistoutleadingStat->At(it));
4157 if(!trackVP) continue;
4158 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
4159
4160 Float_t jetPt = jet->Pt();
4161 Float_t trackPt = trackV->Pt();
4162 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
4163
4164 // Stat plots
4165 if(type==kBckgOutLJStat)
4166 {
4167 if(fFFMode)ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactorLeading);
4168
4169 // Fill track QA for background
4170 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt); // OB added bgr QA
4171 }
4172
4173 // All cases included
4174 if(nRecJetsCuts==1 && type==kBckgOutAJStat)
4175 {
4176 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactorLeading);
4177 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt ); // OB added bgr QA
4178
4179 }
4180 delete trackV;
4181 }
4182 // Increment jet pt with one entry in case #tracks outside jets = 0
4183 if(tracklistoutleadingStat->GetSize()==0) {
4184 Float_t jetPt = jet->Pt();
4185 Bool_t incrementJetPt = kTRUE;
4186 if(type==kBckgOutLJStat)
4187 {
4188 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt, normFactorLeading);
4189 }
4190 // All cases included
4191 if(nRecJetsCuts==1 && type==kBckgOutLJStat)
4192 {
4193 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt, normFactorLeading);
4194 }
4195 }
4196
4197 delete tracklistoutleadingStat;
4198 }
4199
4200 if(type==kBckgPerp || type==kBckgPerp2 || type==kBckgPerp2Area)
4201 {
4202 Double_t sumPtPerp1 = 0.;
4203 Double_t sumPtPerp2 = 0.;
4204 TList* tracklistperp = new TList();
4205 TList* tracklistperp1 = new TList();
4206 TList* tracklistperp2 = new TList();
4207
4208 Double_t norm = 1;
4209 if(type == kBckgPerp2) norm = 2; // in FillFF() scaleFac = 1/norm = 0.5 - account for double area;
4210 if(type == kBckgPerp2Area) norm = 2*TMath::Pi()*TMath::Abs(GetFFRadius())*TMath::Abs(GetFFRadius()) / jet->EffectiveAreaCharged(); // in FillFF() scaleFac = 1/norm;
4211
4212 GetTracksTiltedwrpJetAxis(TMath::Pi()/2., inputtracklist,tracklistperp1,jet,TMath::Abs(GetFFRadius()),sumPtPerp1);
4213 if(type==kBckgPerp2 || type==kBckgPerp2Area) GetTracksTiltedwrpJetAxis(-1*TMath::Pi()/2., inputtracklist,tracklistperp2,jet,TMath::Abs(GetFFRadius()),sumPtPerp2);
4214
4215 tracklistperp->AddAll(tracklistperp1);
4216 tracklistperp->AddAll(tracklistperp2);
4217
4218 if(tracklistperp->GetSize() != tracklistperp1->GetSize() + tracklistperp2->GetSize()){
4219 cout<<" ERROR: tracklistperp size "<<tracklistperp->GetSize()<<" perp1 "<<tracklistperp1->GetSize()<<" perp2 "<<tracklistperp2->GetSize()<<endl;
4220 exit(0);
4221 }
4222
4223 if(fh1Mult) fh1Mult->Fill(tracklistperp->GetSize());
4224
4225 for(Int_t it=0; it<tracklistperp->GetSize(); ++it){
4226
4227 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistperp->At(it));
4228 if(!trackVP)continue;
4229 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
4230
4231 Float_t jetPt = jet->Pt();
4232 Float_t trackPt = trackV->Pt();
4233
4234 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
4235
4236 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, norm );
4237
4238 // Fill track QA for background
4239 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt);
4240
4241 delete trackV;
4242 }
4243 // Increment jet pt with one entry in case #tracks outside jets = 0
4244 if(tracklistperp->GetSize()==0) {
4245 Float_t jetPt = jet->Pt();
4246 Bool_t incrementJetPt = kTRUE;
4247 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt );
4248 }
4249
4250
4251 if(fJSMode){
4252 // fill for tracklistperp1/2 separately, divide norm by 2
4253 if(type==kBckgPerp){
4254 FillJetShape(jet, tracklistperp, fProNtracksLeadingJetBgrPerp2, fProDelRPtSumBgrPerp2, 0, TMath::Pi()/2., 0., kFALSE);
4255 }
4256 if(type==kBckgPerp2){
4257 FillJetShape(jet, tracklistperp1, fProNtracksLeadingJetBgrPerp2, fProDelRPtSumBgrPerp2, 0, TMath::Pi()/2., 0., kFALSE);
4258 FillJetShape(jet, tracklistperp2, fProNtracksLeadingJetBgrPerp2, fProDelRPtSumBgrPerp2, 0, -1*TMath::Pi()/2., 0., kFALSE);
4259 }
4260 if(type==kBckgPerp2Area){ // divide norm by 2: listperp1/2 filled separately
4261 FillJetShape(jet, tracklistperp1, fProNtracksLeadingJetBgrPerp2, fProDelRPtSumBgrPerp2, 0, TMath::Pi()/2., 0.5*norm, kFALSE);
4262 FillJetShape(jet, tracklistperp2, fProNtracksLeadingJetBgrPerp2, fProDelRPtSumBgrPerp2, 0, -1*TMath::Pi()/2., 0.5*norm, kFALSE);
4263 }
4264 }
4265
4266 delete tracklistperp;
4267 delete tracklistperp1;
4268 delete tracklistperp2;
4269 }
4270
4271 if(type==kBckgASide)
4272 {
4273 Double_t sumPtASide = 0.;
4274 TList* tracklistaside = new TList();
4275 GetTracksTiltedwrpJetAxis(TMath::Pi(),inputtracklist,tracklistaside,jet,TMath::Abs(GetFFRadius()),sumPtASide);
4276 if(fh1Mult) fh1Mult->Fill(tracklistaside->GetSize());
4277
4278 for(Int_t it=0; it<tracklistaside->GetSize(); ++it){
4279
4280 AliVParticle* trackVP = (AliVParticle*)(tracklistaside->At(it));
4281 if(!trackVP) continue;
4282 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
4283
4284 Float_t jetPt = jet->Pt();
4285 Float_t trackPt = trackV->Pt();
4286
4287 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
4288
4289 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt );
4290
4291 // Fill track QA for background
4292 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt);
4293
4294 delete trackV;
4295 }
4296 if(tracklistaside->GetSize()==0) {
4297 Float_t jetPt = jet->Pt();
4298 Bool_t incrementJetPt = kTRUE;
4299 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt );
4300 }
4301
4302 delete tracklistaside;
4303 }
4304
4305 if(type==kBckgASideWindow)
4306 {
4307 Double_t normFactorASide = 0.;
4308 Double_t sumPtASideW = 0.;
4309 TList* tracklistasidew = new TList();
4310 GetTracksTiltedwrpJetAxisWindow(TMath::Pi(),inputtracklist,tracklistasidew,jet,TMath::Abs(GetFFRadius()),sumPtASideW,normFactorASide);
4311 if(fh1Mult) fh1Mult->Fill(tracklistasidew->GetSize());
4312
4313 for(Int_t it=0; it<tracklistasidew->GetSize(); ++it){
4314
4315 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistasidew->At(it));
4316 if(!trackVP) continue;
4317 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
4318
4319 Float_t jetPt = jet->Pt();
4320 Float_t trackPt = trackV->Pt();
4321 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
4322
4323 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactorASide);
4324
4325 // Fill track QA for background
4326 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt, kFALSE, normFactorASide);
4327
4328 delete trackV;
4329 }
4330 if(tracklistasidew->GetSize()==0) {
4331 Float_t jetPt = jet->Pt();
4332 Bool_t incrementJetPt = kTRUE;
4333 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt, normFactorASide);
4334 }
4335
4336 delete tracklistasidew;
4337 }
4338
4339 if(type==kBckgPerpWindow)
4340 {
4341 Double_t normFactorPerp = 0.;
4342 Double_t sumPtPerpW = 0.;
4343 TList* tracklistperpw = new TList();
4344 GetTracksTiltedwrpJetAxisWindow(TMath::Pi()/2.,inputtracklist,tracklistperpw,jet,TMath::Abs(GetFFRadius()),sumPtPerpW,normFactorPerp);
4345 if(fh1Mult) fh1Mult->Fill(tracklistperpw->GetSize());
4346
4347 for(Int_t it=0; it<tracklistperpw->GetSize(); ++it){
4348
4349 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistperpw->At(it));
4350 if(!trackVP) continue;
4351 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
4352
4353 Float_t jetPt = jet->Pt();
4354 Float_t trackPt = trackV->Pt();
4355 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
4356
4357 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactorPerp);
4358
4359 // Fill track QA for background
4360 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt, kFALSE, normFactorPerp);
4361
4362 delete trackV;
4363 }
4364 if(tracklistperpw->GetSize()==0) {
4365 Float_t jetPt = jet->Pt();
4366 Bool_t incrementJetPt = kTRUE;
4367 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt, normFactorPerp);
4368 }
4369
4370 delete tracklistperpw;
4371 }
4372
4373
4374 if(type==kBckgOut2J || type==kBckgOutAJ)
4375 {
4376 if(type==kBckgOut2J && fh1Mult) fh1Mult->Fill(tracklistout2jets->GetSize());
4377 for(Int_t it=0; it<tracklistout2jets->GetSize(); ++it){
4378
4379 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistout2jets->At(it));
4380 if(!trackVP) continue;
4381 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
4382
4383 Float_t jetPt = jet->Pt();
4384 Float_t trackPt = trackV->Pt();
4385
4386 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
4387
4388 if(type==kBckgOut2J)
4389 {
4390 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt );
4391 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt);
4392 }
4393
4394 // All cases included
4395 if(nRecJetsCuts==2 && type==kBckgOutAJ)
4396 {
4397 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt );
4398
4399 }
4400 delete trackV;
4401 }
4402 // Increment jet pt with one entry in case #tracks outside jets = 0
4403 if(tracklistout2jets->GetSize()==0) {
4404 Float_t jetPt = jet->Pt();
4405 Bool_t incrementJetPt = kTRUE;
4406 if(type==kBckgOut2J)
4407 {
4408 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt );
4409 }
4410 // All cases included
4411 if(nRecJetsCuts==2 && type==kBckgOutAJ)
4412 {
4413 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt );
4414 }
4415 }
4416 }
4417
4418 if(type==kBckgOut2JStat || type==kBckgOutAJStat)
4419 {
4420 for(Int_t it=0; it<tracklistout2jetsStat->GetSize(); ++it){
4421
4422 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistout2jetsStat->At(it));
4423 if(!trackVP) continue;
4424 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
4425
4426 Float_t jetPt = jet->Pt();
4427 Float_t trackPt = trackV->Pt();
4428 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
4429
4430 if(type==kBckgOut2JStat)
4431 {
4432 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactor2Jets);
4433
4434 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt ); // OB added bgr QA
4435 }
4436
4437 // All cases included
4438 if(nRecJetsCuts==2 && type==kBckgOutAJStat)
4439 {
4440 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactor2Jets);
4441
4442 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt ); // OB added bgr QA
4443 }
4444 delete trackV;
4445 }
4446 // Increment jet pt with one entry in case #tracks outside jets = 0
4447 if(tracklistout2jetsStat->GetSize()==0) {
4448 Float_t jetPt = jet->Pt();
4449 Bool_t incrementJetPt = kTRUE;
4450 if(type==kBckgOut2JStat)
4451 {
4452 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt, normFactor2Jets);
4453 }
4454 // All cases included
4455 if(nRecJetsCuts==2 && type==kBckgOutAJStat)
4456 {
4457 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt, normFactor2Jets);
4458 }
4459 }
4460
4461 }
4462
4463 if(type==kBckgOut3J || type==kBckgOutAJ)
4464 {
4465 if(type==kBckgOut3J && fh1Mult) fh1Mult->Fill(tracklistout3jets->GetSize());
4466
4467 for(Int_t it=0; it<tracklistout3jets->GetSize(); ++it){
4468
4469 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistout3jets->At(it));
4470 if(!trackVP) continue;
4471 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
4472
4473 Float_t jetPt = jet->Pt();
4474 Float_t trackPt = trackV->Pt();
4475
4476 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
4477
4478 if(type==kBckgOut3J)
4479 {
4480 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt );
4481
4482 qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt);
4483 }
4484
4485 // All cases included
4486 if(nRecJetsCuts==3 && type==kBckgOutAJ)
4487 {
4488 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt );
4489
4490 }
4491 delete trackV;
4492 }
4493 // Increment jet pt with one entry in case #tracks outside jets = 0
4494 if(tracklistout3jets->GetSize()==0) {
4495 Float_t jetPt = jet->Pt();
4496 Bool_t incrementJetPt = kTRUE;
4497 if(type==kBckgOut3J)
4498 {
4499 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt );
4500 }
4501 // All cases included
4502 if(nRecJetsCuts==3 && type==kBckgOutAJ)
4503 {
4504 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt );
4505 }
4506 }
4507 }
4508
4509 if(type==kBckgOut3JStat || type==kBckgOutAJStat)
4510 {
4511 for(Int_t it=0; it<tracklistout3jetsStat->GetSize(); ++it){
4512
4513 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistout3jetsStat->At(it));
4514 if(!trackVP) continue;
4515 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
4516
4517 Float_t jetPt = jet->Pt();
4518 Float_t trackPt = trackV->Pt();
4519 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
4520
4521 if(type==kBckgOut3JStat)
4522 {
4523 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactor3Jets);
4524
4525 //if(fQAMode&1) qabckghistocuts->FillTrackQA( trackEta, TVector2::Phi_0_2pi(trackPhi), trackPt);
4526 }
4527
4528 // All cases included
4529 if(nRecJetsCuts==3 && type==kBckgOutAJStat)
4530 {
4531 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactor3Jets );
4532
4533 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt );
4534
4535 }
4536 delete trackV;
4537 }
4538 // Increment jet pt with one entry in case #tracks outside jets = 0
4539 if(tracklistout3jetsStat->GetSize()==0) {
4540 Float_t jetPt = jet->Pt();
4541 Bool_t incrementJetPt = kTRUE;
4542 if(type==kBckgOut3JStat)
4543 {
4544 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt, normFactor3Jets);
4545 }
4546 // All cases included
4547 if(nRecJetsCuts==3 && type==kBckgOutAJStat)
4548 {
4549 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt, normFactor3Jets);
4550 }
4551 }
4552
4553 }
4554
4555 if(type==kBckgClustersOutLeading){ // clusters bgr: all tracks in clusters out of leading jet
4556
4557 TList* tracklistClustersOutLeading = new TList();
4558 Double_t normFactorClusters = 0;
4559 Float_t jetPt = jet->Pt();
4560
4561 GetClusterTracksOutOf1Jet(jet, tracklistClustersOutLeading, normFactorClusters);
4562 if(fh1Mult) fh1Mult->Fill(tracklistClustersOutLeading->GetSize());
4563
4564 for(Int_t it=0; it<tracklistClustersOutLeading->GetSize(); ++it){
4565
4566 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistClustersOutLeading->At(it));
4567 if(!trackVP) continue;
4568 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
4569
4570 Float_t trackPt = trackVP->Pt();
4571
4572 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
4573
4574 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactorClusters );
4575 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt );
4576
4577 delete trackV;
4578 }
4579
4580 delete tracklistClustersOutLeading;
4581
4582 }
4583
4584 if(type == kBckgClusters){ // clusters bgr: all tracks in 'median cluster'
4585
4586 TList* tracklistClustersMedian = new TList();
4587 Double_t normFactorClusters = 0;
4588 Float_t jetPt = jet->Pt();
4589
4590 GetClusterTracksMedian(tracklistClustersMedian, normFactorClusters);
4591 if(fh1Mult) fh1Mult->Fill(tracklistClustersMedian->GetSize());
4592
4593 for(Int_t it=0; it<tracklistClustersMedian->GetSize(); ++it){
4594
4595 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistClustersMedian->At(it));
4596 if(!trackVP) continue;
4597 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
4598
4599 Float_t trackPt = trackVP->Pt();
4600
4601 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
4602
4603 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactorClusters );
4604 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt );
4605
4606 delete trackV;
4607 }
4608
4609 delete tracklistClustersMedian;
4610 }
4611
4612 delete tracklistout2jets;
4613 delete tracklistout3jets;
4614 delete tracklistout2jetsStat;
4615 delete tracklistout3jetsStat;
4616}
4617
2b3fb225 4618//_____________________________________________________________________________________
ffa175ff 4619Double_t AliAnalysisTaskFragmentationFunction::GetMCStrangenessFactor(const Double_t pt)
4620{
4621 // factor strangeness data/MC as function of pt from UE analysis (Sara Vallero)
4622
4623 Double_t alpha = 1;
4624
4625 if(0.150<pt && pt<0.200) alpha = 3.639;
4626 if(0.200<pt && pt<0.250) alpha = 2.097;
4627 if(0.250<pt && pt<0.300) alpha = 1.930;
4628 if(0.300<pt && pt<0.350) alpha = 1.932;
4629 if(0.350<pt && pt<0.400) alpha = 1.943;
4630 if(0.400<pt && pt<0.450) alpha = 1.993;
4631 if(0.450<pt && pt<0.500) alpha = 1.989;
4632 if(0.500<pt && pt<0.600) alpha = 1.963;
4633 if(0.600<pt && pt<0.700) alpha = 1.917;
4634 if(0.700<pt && pt<0.800) alpha = 1.861;
4635 if(0.800<pt && pt<0.900) alpha = 1.820;
4636 if(0.900<pt && pt<1.000) alpha = 1.741;
4637 if(1.000<pt && pt<1.500) alpha = 0.878;
4638
4639 return alpha;
4640}
4641
2b3fb225 4642//__________________________________________________________________________________________________
4643Double_t AliAnalysisTaskFragmentationFunction::GetMCStrangenessFactorCMS(AliAODMCParticle* daughter)
4644{
4645 // strangeness ratio MC/data as function of mother pt from CMS data in |eta|<2.0
4646
4647 TClonesArray *tca = dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName()));
4648 if(!tca) return 1;
4649
4650 AliAODMCParticle* currentMother = daughter;
4651 AliAODMCParticle* currentDaughter = daughter;
4652
4653
4654 // find first primary mother K0s, Lambda or Xi
4655 while(1){
4656
4657 Int_t daughterPDG = currentDaughter->GetPdgCode();
4658
4659 Int_t motherLabel = currentDaughter->GetMother();
4660 if(motherLabel >= tca->GetEntriesFast()){ // protection
4661 currentMother = currentDaughter;
4662 break;
4663 }
4664
4665 currentMother = (AliAODMCParticle*) tca->At(motherLabel);
4666
4667 if(!currentMother){
4668 currentMother = currentDaughter;
4669 break;
4670 }
4671
4672 Int_t motherPDG = currentMother->GetPdgCode();
4673
4674 // phys. primary found ?
4675 if(currentMother->IsPhysicalPrimary()) break;
4676
4677 if(TMath::Abs(daughterPDG) == 321){ // K+/K- e.g. from phi (ref data not feeddown corrected)
4678 currentMother = currentDaughter; break;
4679 }
4680 if(TMath::Abs(motherPDG) == 310 ){ // K0s e.g. from phi (ref data not feeddown corrected)
4681 break;
4682 }
4683 if(TMath::Abs(motherPDG) == 3212 && TMath::Abs(daughterPDG) == 3122){ // mother Sigma0, daughter Lambda (this case not included in feeddown corr.)
4684 currentMother = currentDaughter; break;
4685 }
4686
4687 currentDaughter = currentMother;
4688 }
4689
4690
4691 Int_t motherPDG = currentMother->GetPdgCode();
4692 Double_t motherPt = currentMother->Pt();
4693
4694 Double_t fac = 1;
4695
47830b3f 4696 if(TMath::Abs(motherPDG) == 310 || TMath::Abs(motherPDG)==321){ // K0s / K+ / K-
2b3fb225 4697
4698 if(0.00 <= motherPt && motherPt < 0.20) fac = 0.768049;
4699 else if(0.20 <= motherPt && motherPt < 0.40) fac = 0.732933;
4700 else if(0.40 <= motherPt && motherPt < 0.60) fac = 0.650298;
4701 else if(0.60 <= motherPt && motherPt < 0.80) fac = 0.571332;
4702 else if(0.80 <= motherPt && motherPt < 1.00) fac = 0.518734;
4703 else if(1.00 <= motherPt && motherPt < 1.20) fac = 0.492543;
4704 else if(1.20 <= motherPt && motherPt < 1.40) fac = 0.482704;
4705 else if(1.40 <= motherPt && motherPt < 1.60) fac = 0.488056;
4706 else if(1.60 <= motherPt && motherPt < 1.80) fac = 0.488861;
4707 else if(1.80 <= motherPt && motherPt < 2.00) fac = 0.492862;
4708 else if(2.00 <= motherPt && motherPt < 2.20) fac = 0.504332;
4709 else if(2.20 <= motherPt && motherPt < 2.40) fac = 0.501858;
4710 else if(2.40 <= motherPt && motherPt < 2.60) fac = 0.512970;
4711 else if(2.60 <= motherPt && motherPt < 2.80) fac = 0.524131;
4712 else if(2.80 <= motherPt && motherPt < 3.00) fac = 0.539130;
4713 else if(3.00 <= motherPt && motherPt < 3.20) fac = 0.554101;
4714 else if(3.20 <= motherPt && motherPt < 3.40) fac = 0.560348;
4715 else if(3.40 <= motherPt && motherPt < 3.60) fac = 0.568869;
4716 else if(3.60 <= motherPt && motherPt < 3.80) fac = 0.583310;
4717 else if(3.80 <= motherPt && motherPt < 4.00) fac = 0.604818;
4718 else if(4.00 <= motherPt && motherPt < 5.00) fac = 0.632630;
4719 else if(5.00 <= motherPt && motherPt < 6.00) fac = 0.710070;
4720 else if(6.00 <= motherPt && motherPt < 8.00) fac = 0.736365;
4721 else if(8.00 <= motherPt && motherPt < 10.00) fac = 0.835865;
4722 }
4723
47830b3f 4724 if(TMath::Abs(motherPDG) == 3122){ // Lambda
2b3fb225 4725
4726 if(0.00 <= motherPt && motherPt < 0.20) fac = 0.645162;
4727 else if(0.20 <= motherPt && motherPt < 0.40) fac = 0.627431;
4728 else if(0.40 <= motherPt && motherPt < 0.60) fac = 0.457136;
4729 else if(0.60 <= motherPt && motherPt < 0.80) fac = 0.384369;
4730 else if(0.80 <= motherPt && motherPt < 1.00) fac = 0.330597;
4731 else if(1.00 <= motherPt && motherPt < 1.20) fac = 0.309571;
4732 else if(1.20 <= motherPt && motherPt < 1.40) fac = 0.293620;
4733 else if(1.40 <= motherPt && motherPt < 1.60) fac = 0.283709;
4734 else if(1.60 <= motherPt && motherPt < 1.80) fac = 0.282047;
4735 else if(1.80 <= motherPt && motherPt < 2.00) fac = 0.277261;
4736 else if(2.00 <= motherPt && motherPt < 2.20) fac = 0.275772;
4737 else if(2.20 <= motherPt && motherPt < 2.40) fac = 0.280726;
4738 else if(2.40 <= motherPt && motherPt < 2.60) fac = 0.288540;
4739 else if(2.60 <= motherPt && motherPt < 2.80) fac = 0.288315;
4740 else if(2.80 <= motherPt && motherPt < 3.00) fac = 0.296619;
4741 else if(3.00 <= motherPt && motherPt < 3.20) fac = 0.302993;
4742 else if(3.20 <= motherPt && motherPt < 3.40) fac = 0.338121;
4743 else if(3.40 <= motherPt && motherPt < 3.60) fac = 0.349800;
4744 else if(3.60 <= motherPt && motherPt < 3.80) fac = 0.356802;
4745 else if(3.80 <= motherPt && motherPt < 4.00) fac = 0.391202;
4746 else if(4.00 <= motherPt && motherPt < 5.00) fac = 0.422573;
4747 else if(5.00 <= motherPt && motherPt < 6.00) fac = 0.573815;
4748 else if(6.00 <= motherPt && motherPt < 8.00) fac = 0.786984;
4749 else if(8.00 <= motherPt && motherPt < 10.00) fac = 1.020021;
4750 }
4751
47830b3f 4752 if(TMath::Abs(motherPDG) == 3312 || TMath::Abs(motherPDG) == 3322){ // xi
2b3fb225 4753
4754 if(0.00 <= motherPt && motherPt < 0.20) fac = 0.666620;
4755 else if(0.20 <= motherPt && motherPt < 0.40) fac = 0.575908;
4756 else if(0.40 <= motherPt && motherPt < 0.60) fac = 0.433198;
4757 else if(0.60 <= motherPt && motherPt < 0.80) fac = 0.340901;
4758 else if(0.80 <= motherPt && motherPt < 1.00) fac = 0.290896;
4759 else if(1.00 <= motherPt && motherPt < 1.20) fac = 0.236074;
4760 else if(1.20 <= motherPt && motherPt < 1.40) fac = 0.218681;
4761 else if(1.40 <= motherPt && motherPt < 1.60) fac = 0.207763;
4762 else if(1.60 <= motherPt && motherPt < 1.80) fac = 0.222848;
4763 else if(1.80 <= motherPt && motherPt < 2.00) fac = 0.208806;
4764 else if(2.00 <= motherPt && motherPt < 2.20) fac = 0.197275;
4765 else if(2.20 <= motherPt && motherPt < 2.40) fac = 0.183645;
4766 else if(2.40 <= motherPt && motherPt < 2.60) fac = 0.188788;
4767 else if(2.60 <= motherPt && motherPt < 2.80) fac = 0.188282;
4768 else if(2.80 <= motherPt && motherPt < 3.00) fac = 0.207442;
4769 else if(3.00 <= motherPt && motherPt < 3.20) fac = 0.240388;
4770 else if(3.20 <= motherPt && motherPt < 3.40) fac = 0.241916;
4771 else if(3.40 <= motherPt && motherPt < 3.60) fac = 0.208276;
4772 else if(3.60 <= motherPt && motherPt < 3.80) fac = 0.234550;
4773 else if(3.80 <= motherPt && motherPt < 4.00) fac = 0.251689;
4774 else if(4.00 <= motherPt && motherPt < 5.00) fac = 0.310204;
4775 else if(5.00 <= motherPt && motherPt < 6.00) fac = 0.343492;
4776 }
4777
4778 Double_t weight = 1;
4779 if(fac > 0) weight = 1/fac;
4780
4781 return weight;
4782}
4783
4784// _________________________________________________________________________________
ffa175ff 4785void AliAnalysisTaskFragmentationFunction::FillJetShape(AliAODJet* jet, TList* list,
4786 TProfile* hProNtracksLeadingJet, TProfile** hProDelRPtSum, TProfile* hProDelR80pcPt,
4787 Double_t dPhiUE, Double_t normUE, Bool_t scaleStrangeness){
4788
4789 const Int_t kNbinsR = 50;
4790 const Float_t kBinWidthR = 0.02;
4791
4792 Int_t nJetTracks = list->GetEntries();
4793
4794 Float_t PtSumA[kNbinsR] = {0.0};
ffa175ff 4795
4796 Float_t *delRA = new Float_t[nJetTracks];
4797 Float_t *trackPtA = new Float_t[nJetTracks];
4798 Int_t *index = new Int_t[nJetTracks];
4799
4800 for(Int_t i=0; i<nJetTracks; i++){
4801 delRA[i] = 0;
4802 trackPtA[i] = 0;
4803 index[i] = 0;
4804 }
4805
4806 Double_t jetMom[3];
4807 jet->PxPyPz(jetMom);
4808 TVector3 jet3mom(jetMom);
4809
4810 if(TMath::Abs(dPhiUE)>0){
4811 Double_t phiTilted = jet3mom.Phi();
4812 phiTilted += dPhiUE;
4813 phiTilted = TVector2::Phi_0_2pi(phiTilted);
4814 jet3mom.SetPhi(phiTilted);
4815 }
4816
4817 Double_t jetPt = jet->Pt();
4818 Double_t sumWeights = 0;
4819
4820 for (Int_t j =0; j<nJetTracks; j++){
4821
4822 AliVParticle* track = dynamic_cast<AliVParticle*>(list->At(j));
4823 if(!track)continue;
4824
4825 Double_t trackMom[3];
4826 track->PxPyPz(trackMom);
4827 TVector3 track3mom(trackMom);
4828
4829 Double_t dR = jet3mom.DeltaR(track3mom);
4830
4831 delRA[j] = dR;
4832 trackPtA[j] = track->Pt();
4833
4834 Double_t weight = GetMCStrangenessFactor(track->Pt()); // more correctly should be gen pt
4835 sumWeights += weight;
4836
4837 for(Int_t ibin=1; ibin<=kNbinsR; ibin++){
4838 Float_t xlow = kBinWidthR*(ibin-1);
4839 Float_t xup = kBinWidthR*ibin;
4840 if(xlow <= dR && dR < xup){
fd5c68ba 4841
4842 if(scaleStrangeness) PtSumA[ibin-1] += track->Pt()*weight;
4843 else PtSumA[ibin-1] += track->Pt();
ffa175ff 4844 }
4845 }
4846 } // track loop
4847
4848 Float_t jetPtMin=0;
4849 Float_t jetPtMax=0;
4850
4851 for(Int_t ibin=0; ibin<kNbinsR; ibin++){
4852 Float_t fR = kBinWidthR*(ibin+0.5);
4853
4854 for(Int_t k=0; k<5; k++){
4855 if(k==0){jetPtMin=20.0;jetPtMax=30.0;}
4856 if(k==1){jetPtMin=30.0;jetPtMax=40.0;}
4857 if(k==2){jetPtMin=40.0;jetPtMax=60.0;}
4858 if(k==3){jetPtMin=60.0;jetPtMax=80.0;}
4859 if(k==4){jetPtMin=80.0;jetPtMax=100.0;}
4860 if(jetPt>jetPtMin && jetPt<jetPtMax){
4861
fd5c68ba 4862 hProDelRPtSum[k]->Fill(fR,PtSumA[ibin]);
4863
ffa175ff 4864 }
4865 }
4866 }
4867
4868 if(scaleStrangeness) hProNtracksLeadingJet->Fill(jetPt,sumWeights);
4869 else hProNtracksLeadingJet->Fill(jetPt,nJetTracks);
4870
4871 if(normUE) hProNtracksLeadingJet->Fill(jetPt,nJetTracks/normUE);
4872
4873 if(hProDelR80pcPt){
4874
4875 Float_t PtSum = 0;
4876 Float_t delRPtSum80pc = 0;
4877
4878 TMath::Sort(nJetTracks,delRA,index,0);
4879
4880 for(Int_t ii=0; ii<nJetTracks; ii++){
4881
4882 if(scaleStrangeness){
4883 Double_t weight = GetMCStrangenessFactor(trackPtA[index[ii]]); // more correctly should be gen pt
4884 PtSum += weight*trackPtA[index[ii]];
4885 }
4886 else PtSum += trackPtA[index[ii]];
4887
4888
4889 if(PtSum/jetPt >= 0.8000){
4890 delRPtSum80pc = delRA[index[ii]];
4891 break;
4892 }
4893 }
4894 hProDelR80pcPt->Fill(jetPt,delRPtSum80pc);
4895 }
4896
4897 delete[] delRA;
4898 delete[] trackPtA;
4899 delete[] index;
4900}