Updated Fragmentation function task (added Di jet part and intra-jet radiation) ...
[u/mrichter/AliRoot.git] / PWG4 / JetTasks / AliAnalysisTaskFragmentationFunction.cxx
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 "TString.h"
29 #include "THnSparse.h"
30
31 #include "AliAODInputHandler.h" 
32 #include "AliAODHandler.h" 
33 #include "AliESDEvent.h"
34 #include "AliAODMCParticle.h"
35 #include "AliAODJet.h"
36 #include "AliGenPythiaEventHeader.h"
37 #include "AliInputEventHandler.h"
38
39 #include "AliAnalysisHelperJetTasks.h"
40 #include "AliAnalysisManager.h"
41 #include "AliAnalysisTaskSE.h"
42
43 #include "AliAnalysisTaskFragmentationFunction.h"
44
45
46 ClassImp(AliAnalysisTaskFragmentationFunction)
47
48 //____________________________________________________________________________
49 AliAnalysisTaskFragmentationFunction::AliAnalysisTaskFragmentationFunction()
50    : AliAnalysisTaskSE()
51    ,fESD(0)
52    ,fAOD(0)
53    ,fMCEvent(0)
54    ,fBranchRecJets("jets")
55    ,fBranchGenJets("")
56    ,fTrackTypeGen(0)
57    ,fJetTypeGen(0)
58    ,fFilterMask(0)
59    ,fTrackPtCut(0)
60    ,fTrackEtaMin(0)
61    ,fTrackEtaMax(0)
62    ,fTrackPhiMin(0)
63    ,fTrackPhiMax(0)
64    ,fJetPtCut(0)
65    ,fJetEtaMin(0)
66    ,fJetEtaMax(0)
67    ,fJetPhiMin(0)
68    ,fJetPhiMax(0)
69    ,fDiJetCut(0)
70    ,fDiJetDeltaPhiCut(0)
71    ,fDiJetPtFractionCut(0)
72    ,fDiJetCDFCut(0)
73    ,fDiJetKindBins(0)
74    ,fFFRadius(0)
75    ,fTracksRec(0)
76    ,fTracksRecCuts(0)
77    ,fTracksGen(0)
78    ,fJetsRec(0)
79    ,fJetsRecCuts(0)
80    ,fJetsGen(0)
81    ,fQATrackHistosRec(0)
82    ,fQATrackHistosRecCuts(0)
83    ,fQATrackHistosGen(0)
84    ,fQAJetHistosRec(0)
85    ,fQAJetHistosRecCuts(0)
86    ,fQAJetHistosRecCutsLeading(0)
87    ,fQAJetHistosGen(0)
88    ,fQAJetHistosGenLeading(0)
89    ,fFFHistosRecCuts(0)
90    ,fFFHistosRecLeading(0)
91    ,fFFHistosRecLeadingTrack(0)
92    ,fFFHistosGen(0)
93    ,fFFHistosGenLeading(0)
94    ,fFFHistosGenLeadingTrack(0)
95    ,fIJHistosRecCuts(0)
96    ,fIJHistosRecLeading(0)
97    ,fIJHistosRecLeadingTrack(0)
98    ,fIJHistosGen(0)
99    ,fIJHistosGenLeading(0)
100    ,fIJHistosGenLeadingTrack(0)
101    ,fFFDiJetHistosRecCuts(0)
102    ,fFFDiJetHistosRecLeading(0)
103    ,fFFDiJetHistosRecLeadingTrack(0)
104    ,fFFDiJetHistosGen(0)
105    ,fFFDiJetHistosGenLeading(0)
106    ,fFFDiJetHistosGenLeadingTrack(0)
107    ,fQADiJetHistosRecCuts(0)
108    ,fQADiJetHistosGen(0)
109    ,fQATrackHighPtThreshold(0)
110    ,fFFNBinsJetPt(0)    
111    ,fFFJetPtMin(0) 
112    ,fFFJetPtMax(0)
113    ,fFFNBinsPt(0)      
114    ,fFFPtMin(0)        
115    ,fFFPtMax(0)        
116    ,fFFNBinsXi(0)      
117    ,fFFXiMin(0)        
118    ,fFFXiMax(0)        
119    ,fFFNBinsZ(0)       
120    ,fFFZMin(0)         
121    ,fFFZMax(0)         
122    ,fQAJetNBinsPt(0)   
123    ,fQAJetPtMin(0)     
124    ,fQAJetPtMax(0)     
125    ,fQAJetNBinsEta(0)  
126    ,fQAJetEtaMin(0)    
127    ,fQAJetEtaMax(0)    
128    ,fQAJetNBinsPhi(0)  
129    ,fQAJetPhiMin(0)    
130    ,fQAJetPhiMax(0)    
131    ,fQATrackNBinsPt(0) 
132    ,fQATrackPtMin(0)   
133    ,fQATrackPtMax(0)   
134    ,fQATrackNBinsEta(0)
135    ,fQATrackEtaMin(0)  
136    ,fQATrackEtaMax(0)  
137    ,fQATrackNBinsPhi(0)
138    ,fQATrackPhiMin(0)  
139    ,fQATrackPhiMax(0)
140    ,fIJNBinsJetPt(0)
141    ,fIJJetPtMin(0)
142    ,fIJJetPtMax(0)
143    ,fIJNBinsPt(0)
144    ,fIJPtMin(0)
145    ,fIJPtMax(0)
146    ,fIJNBinsZ(0)
147    ,fIJZMin(0)
148    ,fIJZMax(0)
149    ,fIJNBinsCosTheta(0)
150    ,fIJCosThetaMin(0)
151    ,fIJCosThetaMax(0)
152    ,fIJNBinsTheta(0)
153    ,fIJThetaMin(0)
154    ,fIJThetaMax(0)
155    ,fIJNBinsJt(0)
156    ,fIJJtMin(0)
157    ,fIJJtMax(0)
158    ,fDiJetNBinsJetInvMass(0)
159    ,fDiJetJetInvMassMin(0)
160    ,fDiJetJetInvMassMax(0)
161    ,fDiJetNBinsJetPt(0)
162    ,fDiJetJetPtMin(0)
163    ,fDiJetJetPtMax(0)
164    ,fDiJetNBinsPt(0)
165    ,fDiJetPtMin(0)
166    ,fDiJetPtMax(0)
167    ,fDiJetNBinsXi(0)
168    ,fDiJetXiMin(0)
169    ,fDiJetXiMax(0)
170    ,fDiJetNBinsZ(0)
171    ,fDiJetZMin(0)
172    ,fDiJetZMax(0)
173    ,fQADiJetNBinsInvMass(0)
174    ,fQADiJetInvMassMin(0)
175    ,fQADiJetInvMassMax(0)
176    ,fQADiJetNBinsJetPt(0)
177    ,fQADiJetJetPtMin(0)
178    ,fQADiJetJetPtMax(0)
179    ,fQADiJetNBinsDeltaPhi(0)
180    ,fQADiJetDeltaPhiMin(0)
181    ,fQADiJetDeltaPhiMax(0)
182    ,fQADiJetNBinsDeltaEta(0)
183    ,fQADiJetDeltaEtaMin(0)
184    ,fQADiJetDeltaEtaMax(0)
185    ,fQADiJetNBinsDeltaPt(0)
186    ,fQADiJetDeltaPtMin(0)
187    ,fQADiJetDeltaPtMax(0)
188    ,fCommonHistList(0)
189    ,fh1EvtSelection(0)
190    ,fh1VertexNContributors(0)
191    ,fh1VertexZ(0)
192    ,fh1EvtMult(0)
193    ,fh1nRecJetsCuts(0)
194    ,fh1nGenJets(0)
195 {
196    // default constructor
197 }
198
199 //__________________________________________________________________________________________
200 AliAnalysisTaskFragmentationFunction::AliAnalysisTaskFragmentationFunction(const char *name) 
201   : AliAnalysisTaskSE(name)
202   ,fESD(0)
203   ,fAOD(0)
204   ,fMCEvent(0)
205   ,fBranchRecJets("jets")
206   ,fBranchGenJets("")
207   ,fTrackTypeGen(0)
208   ,fJetTypeGen(0)
209   ,fFilterMask(0)
210   ,fTrackPtCut(0)
211   ,fTrackEtaMin(0)
212   ,fTrackEtaMax(0)
213   ,fTrackPhiMin(0)
214   ,fTrackPhiMax(0)
215   ,fJetPtCut(0)
216   ,fJetEtaMin(0)
217   ,fJetEtaMax(0)
218   ,fJetPhiMin(0)
219   ,fJetPhiMax(0)
220   ,fDiJetCut(0)
221   ,fDiJetDeltaPhiCut(0)
222   ,fDiJetPtFractionCut(0)
223   ,fDiJetCDFCut(0)
224   ,fDiJetKindBins(0)
225   ,fFFRadius(0)
226   ,fTracksRec(0)
227   ,fTracksRecCuts(0)
228   ,fTracksGen(0)
229   ,fJetsRec(0)
230   ,fJetsRecCuts(0)
231   ,fJetsGen(0)
232   ,fQATrackHistosRec(0)
233   ,fQATrackHistosRecCuts(0)
234   ,fQATrackHistosGen(0)
235   ,fQAJetHistosRec(0)
236   ,fQAJetHistosRecCuts(0)
237   ,fQAJetHistosRecCutsLeading(0)
238   ,fQAJetHistosGen(0)
239   ,fQAJetHistosGenLeading(0)
240   ,fFFHistosRecCuts(0)
241   ,fFFHistosRecLeading(0)
242   ,fFFHistosRecLeadingTrack(0)
243   ,fFFHistosGen(0)
244   ,fFFHistosGenLeading(0)
245   ,fFFHistosGenLeadingTrack(0)
246   ,fIJHistosRecCuts(0)
247   ,fIJHistosRecLeading(0)
248   ,fIJHistosRecLeadingTrack(0)
249   ,fIJHistosGen(0)
250   ,fIJHistosGenLeading(0)
251   ,fIJHistosGenLeadingTrack(0)
252   ,fFFDiJetHistosRecCuts(0)
253   ,fFFDiJetHistosRecLeading(0)
254   ,fFFDiJetHistosRecLeadingTrack(0)
255   ,fFFDiJetHistosGen(0)
256   ,fFFDiJetHistosGenLeading(0)
257   ,fFFDiJetHistosGenLeadingTrack(0)
258   ,fQADiJetHistosRecCuts(0)
259   ,fQADiJetHistosGen(0)
260   ,fQATrackHighPtThreshold(0) 
261   ,fFFNBinsJetPt(0)    
262   ,fFFJetPtMin(0) 
263   ,fFFJetPtMax(0)
264   ,fFFNBinsPt(0)      
265   ,fFFPtMin(0)        
266   ,fFFPtMax(0)        
267   ,fFFNBinsXi(0)      
268   ,fFFXiMin(0)        
269   ,fFFXiMax(0)        
270   ,fFFNBinsZ(0)       
271   ,fFFZMin(0)         
272   ,fFFZMax(0)         
273   ,fQAJetNBinsPt(0)   
274   ,fQAJetPtMin(0)     
275   ,fQAJetPtMax(0)     
276   ,fQAJetNBinsEta(0)  
277   ,fQAJetEtaMin(0)    
278   ,fQAJetEtaMax(0)    
279   ,fQAJetNBinsPhi(0)  
280   ,fQAJetPhiMin(0)    
281   ,fQAJetPhiMax(0)    
282   ,fQATrackNBinsPt(0) 
283   ,fQATrackPtMin(0)   
284   ,fQATrackPtMax(0)   
285   ,fQATrackNBinsEta(0)
286   ,fQATrackEtaMin(0)  
287   ,fQATrackEtaMax(0)  
288   ,fQATrackNBinsPhi(0)
289   ,fQATrackPhiMin(0)  
290   ,fQATrackPhiMax(0)  
291   ,fIJNBinsJetPt(0)
292   ,fIJJetPtMin(0)
293   ,fIJJetPtMax(0)
294   ,fIJNBinsPt(0)
295   ,fIJPtMin(0)
296   ,fIJPtMax(0)
297   ,fIJNBinsZ(0)
298   ,fIJZMin(0)
299   ,fIJZMax(0)
300   ,fIJNBinsCosTheta(0)
301   ,fIJCosThetaMin(0)
302   ,fIJCosThetaMax(0)
303   ,fIJNBinsTheta(0)
304   ,fIJThetaMin(0)
305   ,fIJThetaMax(0)
306   ,fIJNBinsJt(0)
307   ,fIJJtMin(0)
308   ,fIJJtMax(0)
309   ,fDiJetNBinsJetInvMass(0)
310   ,fDiJetJetInvMassMin(0)
311   ,fDiJetJetInvMassMax(0)
312   ,fDiJetNBinsJetPt(0)
313   ,fDiJetJetPtMin(0)
314   ,fDiJetJetPtMax(0)
315   ,fDiJetNBinsPt(0)
316   ,fDiJetPtMin(0)
317   ,fDiJetPtMax(0)
318   ,fDiJetNBinsXi(0)
319   ,fDiJetXiMin(0)
320   ,fDiJetXiMax(0)
321   ,fDiJetNBinsZ(0)
322   ,fDiJetZMin(0)
323   ,fDiJetZMax(0)
324   ,fQADiJetNBinsInvMass(0)
325   ,fQADiJetInvMassMin(0)
326   ,fQADiJetInvMassMax(0)
327   ,fQADiJetNBinsJetPt(0)
328   ,fQADiJetJetPtMin(0)
329   ,fQADiJetJetPtMax(0)
330   ,fQADiJetNBinsDeltaPhi(0)
331   ,fQADiJetDeltaPhiMin(0)
332   ,fQADiJetDeltaPhiMax(0)
333   ,fQADiJetNBinsDeltaEta(0)
334   ,fQADiJetDeltaEtaMin(0)
335   ,fQADiJetDeltaEtaMax(0)
336   ,fQADiJetNBinsDeltaPt(0)
337   ,fQADiJetDeltaPtMin(0)
338   ,fQADiJetDeltaPtMax(0)
339   ,fCommonHistList(0)
340   ,fh1EvtSelection(0)
341   ,fh1VertexNContributors(0)
342   ,fh1VertexZ(0)
343   ,fh1EvtMult(0)
344   ,fh1nRecJetsCuts(0)
345   ,fh1nGenJets(0)
346 {
347   // constructor
348   
349   DefineOutput(1,TList::Class());
350   
351
352 }
353
354 //__________________________________________________________________________________________________________________________
355 AliAnalysisTaskFragmentationFunction::AliAnalysisTaskFragmentationFunction(const  AliAnalysisTaskFragmentationFunction &copy)
356   : AliAnalysisTaskSE()
357   ,fESD(copy.fESD)
358   ,fAOD(copy.fAOD)
359   ,fMCEvent(copy.fMCEvent)
360   ,fBranchRecJets(copy.fBranchRecJets)
361   ,fBranchGenJets(copy.fBranchGenJets)
362   ,fTrackTypeGen(copy.fTrackTypeGen)
363   ,fJetTypeGen(copy.fJetTypeGen)
364   ,fFilterMask(copy.fFilterMask)
365   ,fTrackPtCut(copy.fTrackPtCut)
366   ,fTrackEtaMin(copy.fTrackEtaMin)
367   ,fTrackEtaMax(copy.fTrackEtaMax)
368   ,fTrackPhiMin(copy.fTrackPhiMin)
369   ,fTrackPhiMax(copy.fTrackPhiMax)
370   ,fJetPtCut(copy.fJetPtCut)
371   ,fJetEtaMin(copy.fJetEtaMin)
372   ,fJetEtaMax(copy.fJetEtaMax)
373   ,fJetPhiMin(copy.fJetPhiMin)
374   ,fJetPhiMax(copy.fJetPhiMin)
375   ,fDiJetCut(copy.fDiJetCut)
376   ,fDiJetDeltaPhiCut(copy.fDiJetDeltaPhiCut)
377   ,fDiJetPtFractionCut(copy.fDiJetPtFractionCut)
378   ,fDiJetCDFCut(copy.fDiJetCDFCut)
379   ,fDiJetKindBins(copy.fDiJetKindBins)
380   ,fFFRadius(copy.fFFRadius)
381   ,fTracksRec(copy.fTracksRec)
382   ,fTracksRecCuts(copy.fTracksRecCuts)
383   ,fTracksGen(copy.fTracksGen)
384   ,fJetsRec(copy.fJetsRec)
385   ,fJetsRecCuts(copy.fJetsRecCuts)
386   ,fJetsGen(copy.fJetsGen)
387   ,fQATrackHistosRec(copy.fQATrackHistosRec)
388   ,fQATrackHistosRecCuts(copy.fQATrackHistosRecCuts)
389   ,fQATrackHistosGen(copy.fQATrackHistosGen)
390   ,fQAJetHistosRec(copy.fQAJetHistosRec)
391   ,fQAJetHistosRecCuts(copy.fQAJetHistosRecCuts)
392   ,fQAJetHistosRecCutsLeading(copy.fQAJetHistosRecCutsLeading)
393   ,fQAJetHistosGen(copy.fQAJetHistosGen)
394   ,fQAJetHistosGenLeading(copy.fQAJetHistosGenLeading)
395   ,fFFHistosRecCuts(copy.fFFHistosRecCuts)
396   ,fFFHistosRecLeading(copy.fFFHistosRecLeading)
397   ,fFFHistosRecLeadingTrack(copy.fFFHistosRecLeadingTrack)
398   ,fFFHistosGen(copy.fFFHistosGen)
399   ,fFFHistosGenLeading(copy.fFFHistosGenLeading)
400   ,fFFHistosGenLeadingTrack(copy.fFFHistosGenLeadingTrack)
401   ,fIJHistosRecCuts(copy.fIJHistosRecCuts)
402   ,fIJHistosRecLeading(copy.fIJHistosRecLeading)
403   ,fIJHistosRecLeadingTrack(copy.fIJHistosRecLeadingTrack)
404   ,fIJHistosGen(copy.fIJHistosGen)
405   ,fIJHistosGenLeading(copy.fIJHistosGenLeading)
406   ,fIJHistosGenLeadingTrack(copy.fIJHistosGenLeadingTrack)
407   ,fFFDiJetHistosRecCuts(copy.fFFDiJetHistosRecCuts)
408   ,fFFDiJetHistosRecLeading(copy.fFFDiJetHistosRecLeading)
409   ,fFFDiJetHistosRecLeadingTrack(copy.fFFDiJetHistosRecLeadingTrack)
410   ,fFFDiJetHistosGen(copy.fFFDiJetHistosGen)
411   ,fFFDiJetHistosGenLeading(copy.fFFDiJetHistosGenLeading)
412   ,fFFDiJetHistosGenLeadingTrack(copy.fFFDiJetHistosGenLeadingTrack)
413   ,fQADiJetHistosRecCuts(copy.fQADiJetHistosRecCuts)
414   ,fQADiJetHistosGen(copy.fQADiJetHistosGen)
415   ,fQATrackHighPtThreshold(copy.fQATrackHighPtThreshold) 
416   ,fFFNBinsJetPt(copy.fFFNBinsJetPt)    
417   ,fFFJetPtMin(copy.fFFJetPtMin) 
418   ,fFFJetPtMax(copy.fFFJetPtMax)
419   ,fFFNBinsPt(copy.fFFNBinsPt)      
420   ,fFFPtMin(copy.fFFPtMin)        
421   ,fFFPtMax(copy.fFFPtMax)        
422   ,fFFNBinsXi(copy.fFFNBinsXi)      
423   ,fFFXiMin(copy.fFFXiMin)        
424   ,fFFXiMax(copy.fFFXiMax)        
425   ,fFFNBinsZ(copy.fFFNBinsZ)       
426   ,fFFZMin(copy.fFFZMin)         
427   ,fFFZMax(copy.fFFZMax)         
428   ,fQAJetNBinsPt(copy.fQAJetNBinsPt)   
429   ,fQAJetPtMin(copy.fQAJetPtMin)     
430   ,fQAJetPtMax(copy.fQAJetPtMax)     
431   ,fQAJetNBinsEta(copy.fQAJetNBinsEta)  
432   ,fQAJetEtaMin(copy.fQAJetEtaMin)    
433   ,fQAJetEtaMax(copy.fQAJetEtaMax)    
434   ,fQAJetNBinsPhi(copy.fQAJetNBinsPhi)  
435   ,fQAJetPhiMin(copy.fQAJetPhiMin)    
436   ,fQAJetPhiMax(copy.fQAJetPhiMax)    
437   ,fQATrackNBinsPt(copy.fQATrackNBinsPt) 
438   ,fQATrackPtMin(copy.fQATrackPtMin)   
439   ,fQATrackPtMax(copy.fQATrackPtMax)   
440   ,fQATrackNBinsEta(copy.fQATrackNBinsEta)
441   ,fQATrackEtaMin(copy.fQATrackEtaMin)  
442   ,fQATrackEtaMax(copy.fQATrackEtaMax)  
443   ,fQATrackNBinsPhi(copy.fQATrackNBinsPhi)
444   ,fQATrackPhiMin(copy.fQATrackPhiMin)  
445   ,fQATrackPhiMax(copy.fQATrackPhiMax)
446   ,fIJNBinsJetPt(copy.fIJNBinsJetPt)
447   ,fIJJetPtMin(copy.fIJJetPtMin)
448   ,fIJJetPtMax(copy.fIJJetPtMax)
449   ,fIJNBinsPt(copy.fIJNBinsPt)
450   ,fIJPtMin(copy.fIJPtMin)
451   ,fIJPtMax(copy.fIJPtMax)
452   ,fIJNBinsZ(copy.fIJNBinsZ)
453   ,fIJZMin(copy.fIJZMin)
454   ,fIJZMax(copy.fIJZMax)
455   ,fIJNBinsCosTheta(copy.fIJNBinsCosTheta)
456   ,fIJCosThetaMin(copy.fIJCosThetaMin)
457   ,fIJCosThetaMax(copy.fIJCosThetaMax)
458   ,fIJNBinsTheta(copy.fIJNBinsTheta)
459   ,fIJThetaMin(copy.fIJThetaMin)
460   ,fIJThetaMax(copy.fIJThetaMax)
461   ,fIJNBinsJt(copy.fIJNBinsJt)
462   ,fIJJtMin(copy.fIJJtMin)
463   ,fIJJtMax(copy.fIJJtMax)
464   ,fDiJetNBinsJetInvMass(copy.fDiJetNBinsJetInvMass)
465   ,fDiJetJetInvMassMin(copy.fDiJetJetInvMassMin)
466   ,fDiJetJetInvMassMax(copy.fDiJetJetInvMassMax)
467   ,fDiJetNBinsJetPt(copy.fDiJetNBinsJetPt)
468   ,fDiJetJetPtMin(copy.fDiJetJetPtMin)
469   ,fDiJetJetPtMax(copy.fDiJetJetPtMax)
470   ,fDiJetNBinsPt(copy.fDiJetNBinsPt)
471   ,fDiJetPtMin(copy.fDiJetPtMin)
472   ,fDiJetPtMax(copy.fDiJetPtMax)
473   ,fDiJetNBinsXi(copy.fDiJetNBinsXi)
474   ,fDiJetXiMin(copy.fDiJetXiMin)
475   ,fDiJetXiMax(copy.fDiJetXiMax)
476   ,fDiJetNBinsZ(copy.fDiJetNBinsZ)
477   ,fDiJetZMin(copy.fDiJetZMin)
478   ,fDiJetZMax(copy.fDiJetZMax)
479   ,fQADiJetNBinsInvMass(copy.fQADiJetNBinsInvMass)
480   ,fQADiJetInvMassMin(copy.fQADiJetInvMassMin)
481   ,fQADiJetInvMassMax(copy.fQADiJetInvMassMax)
482   ,fQADiJetNBinsJetPt(copy.fQADiJetNBinsJetPt)
483   ,fQADiJetJetPtMin(copy.fQADiJetJetPtMin)
484   ,fQADiJetJetPtMax(copy.fQADiJetJetPtMax)
485   ,fQADiJetNBinsDeltaPhi(copy.fQADiJetNBinsDeltaPhi)
486   ,fQADiJetDeltaPhiMin(copy.fQADiJetDeltaPhiMin)
487   ,fQADiJetDeltaPhiMax(copy.fQADiJetDeltaPhiMax)
488   ,fQADiJetNBinsDeltaEta(copy.fQADiJetNBinsDeltaEta)
489   ,fQADiJetDeltaEtaMin(copy.fQADiJetDeltaEtaMin)
490   ,fQADiJetDeltaEtaMax(copy.fQADiJetDeltaEtaMax)
491   ,fQADiJetNBinsDeltaPt(copy.fQADiJetNBinsDeltaPt)
492   ,fQADiJetDeltaPtMin(copy.fQADiJetDeltaPtMin)
493   ,fQADiJetDeltaPtMax(copy.fQADiJetDeltaPtMax)
494   ,fCommonHistList(copy.fCommonHistList)
495   ,fh1EvtSelection(copy.fh1EvtSelection)
496   ,fh1VertexNContributors(copy.fh1VertexNContributors)
497   ,fh1VertexZ(copy.fh1VertexZ)
498   ,fh1EvtMult(copy.fh1EvtMult)
499   ,fh1nRecJetsCuts(copy.fh1nRecJetsCuts)
500   ,fh1nGenJets(copy.fh1nGenJets)
501 {
502   // copy constructor
503
504 }
505
506 // _________________________________________________________________________________________________________________________________
507 AliAnalysisTaskFragmentationFunction& AliAnalysisTaskFragmentationFunction::operator=(const AliAnalysisTaskFragmentationFunction& o)
508 {
509   // assignment
510   
511   if(this!=&o){
512
513     AliAnalysisTaskSE::operator=(o);
514     fESD                          = o.fESD;
515     fAOD                          = o.fAOD;
516     fMCEvent                      = o.fMCEvent;
517     fBranchRecJets                = o.fBranchRecJets;
518     fBranchGenJets                = o.fBranchGenJets;
519     fTrackTypeGen                 = o.fTrackTypeGen;
520     fJetTypeGen                   = o.fJetTypeGen;
521     fFilterMask                   = o.fFilterMask;
522     fTrackPtCut                   = o.fTrackPtCut;
523     fTrackEtaMin                  = o.fTrackEtaMin;
524     fTrackEtaMax                  = o.fTrackEtaMax;
525     fTrackPhiMin                  = o.fTrackPhiMin;
526     fTrackPhiMax                  = o.fTrackPhiMax;
527     fJetPtCut                     = o.fJetPtCut;
528     fJetEtaMin                    = o.fJetEtaMin;
529     fJetEtaMax                    = o.fJetEtaMax;
530     fJetPhiMin                    = o.fJetPhiMin;
531     fJetPhiMax                    = o.fJetPhiMin;
532     fDiJetCut                     = o.fDiJetCut;
533     fDiJetDeltaPhiCut             = o.fDiJetDeltaPhiCut;
534     fDiJetPtFractionCut           = o.fDiJetPtFractionCut;
535     fDiJetCDFCut                  = o.fDiJetCDFCut;
536     fDiJetKindBins                = o.fDiJetKindBins;
537     fFFRadius                     = o.fFFRadius;
538     fTracksRec                    = o.fTracksRec;
539     fTracksRecCuts                = o.fTracksRecCuts;
540     fTracksGen                    = o.fTracksGen;
541     fJetsRec                      = o.fJetsRec;
542     fJetsRecCuts                  = o.fJetsRecCuts;
543     fJetsGen                      = o.fJetsGen;
544     fQATrackHistosRec             = o.fQATrackHistosRec;
545     fQATrackHistosRecCuts         = o.fQATrackHistosRecCuts;
546     fQATrackHistosGen             = o.fQATrackHistosGen;
547     fQAJetHistosRec               = o.fQAJetHistosRec;
548     fQAJetHistosRecCuts           = o.fQAJetHistosRecCuts;
549     fQAJetHistosRecCutsLeading    = o.fQAJetHistosRecCutsLeading;
550     fQAJetHistosGen               = o.fQAJetHistosGen;
551     fQAJetHistosGenLeading        = o.fQAJetHistosGenLeading;
552     fFFHistosRecCuts              = o.fFFHistosRecCuts;
553     fFFHistosRecLeading           = o.fFFHistosRecLeading;
554     fFFHistosRecLeadingTrack      = o.fFFHistosRecLeadingTrack;
555     fFFHistosGen                  = o.fFFHistosGen;
556     fFFHistosGenLeading           = o.fFFHistosGenLeading;
557     fFFHistosGenLeadingTrack      = o.fFFHistosGenLeadingTrack;
558     fIJHistosRecCuts              = o.fIJHistosRecCuts;
559     fIJHistosRecLeading           = o.fIJHistosRecLeading;
560     fIJHistosRecLeadingTrack      = o.fIJHistosRecLeadingTrack;
561     fIJHistosGen                  = o.fIJHistosGen;
562     fIJHistosGenLeading           = o.fIJHistosGenLeading;
563     fIJHistosGenLeadingTrack      = o.fIJHistosGenLeadingTrack;
564     fFFDiJetHistosRecCuts         = o.fFFDiJetHistosRecCuts;
565     fFFDiJetHistosRecLeading      = o.fFFDiJetHistosRecLeading;
566     fFFDiJetHistosRecLeadingTrack = o.fFFDiJetHistosRecLeadingTrack;
567     fFFDiJetHistosGen             = o.fFFDiJetHistosGen;
568     fFFDiJetHistosGenLeading      = o.fFFDiJetHistosGenLeading;
569     fFFDiJetHistosGenLeadingTrack = o.fFFDiJetHistosGenLeadingTrack;
570     fQADiJetHistosRecCuts         = o.fQADiJetHistosRecCuts;
571     fQADiJetHistosGen             = o.fQADiJetHistosGen;
572     fQATrackHighPtThreshold       = o.fQATrackHighPtThreshold; 
573     fFFNBinsJetPt                 = o.fFFNBinsJetPt;    
574     fFFJetPtMin                   = o.fFFJetPtMin; 
575     fFFJetPtMax                   = o.fFFJetPtMax;
576     fFFNBinsPt                    = o.fFFNBinsPt;      
577     fFFPtMin                      = o.fFFPtMin;        
578     fFFPtMax                      = o.fFFPtMax;        
579     fFFNBinsXi                    = o.fFFNBinsXi;      
580     fFFXiMin                      = o.fFFXiMin;        
581     fFFXiMax                      = o.fFFXiMax;        
582     fFFNBinsZ                     = o.fFFNBinsZ;       
583     fFFZMin                       = o.fFFZMin;         
584     fFFZMax                       = o.fFFZMax;         
585     fQAJetNBinsPt                 = o.fQAJetNBinsPt;   
586     fQAJetPtMin                   = o.fQAJetPtMin;     
587     fQAJetPtMax                   = o.fQAJetPtMax;     
588     fQAJetNBinsEta                = o.fQAJetNBinsEta;  
589     fQAJetEtaMin                  = o.fQAJetEtaMin;    
590     fQAJetEtaMax                  = o.fQAJetEtaMax;    
591     fQAJetNBinsPhi                = o.fQAJetNBinsPhi;  
592     fQAJetPhiMin                  = o.fQAJetPhiMin;    
593     fQAJetPhiMax                  = o.fQAJetPhiMax;    
594     fQATrackNBinsPt               = o.fQATrackNBinsPt; 
595     fQATrackPtMin                 = o.fQATrackPtMin;   
596     fQATrackPtMax                 = o.fQATrackPtMax;   
597     fQATrackNBinsEta              = o.fQATrackNBinsEta;
598     fQATrackEtaMin                = o.fQATrackEtaMin;  
599     fQATrackEtaMax                = o.fQATrackEtaMax;  
600     fQATrackNBinsPhi              = o.fQATrackNBinsPhi;
601     fQATrackPhiMin                = o.fQATrackPhiMin;  
602     fQATrackPhiMax                = o.fQATrackPhiMax;  
603     fIJNBinsJetPt                 = o.fIJNBinsJetPt;
604     fIJJetPtMin                   = o.fIJJetPtMin;
605     fIJJetPtMax                   = o.fIJJetPtMax;
606     fIJNBinsPt                    = o.fIJNBinsPt;
607     fIJPtMin                      = o.fIJPtMin;
608     fIJPtMax                      = o.fIJPtMax;
609     fIJNBinsZ                     = o.fIJNBinsZ;
610     fIJZMin                       = o.fIJZMin;
611     fIJZMax                       = o.fIJZMax;
612     fIJNBinsCosTheta              = o.fIJNBinsCosTheta;
613     fIJCosThetaMin                = o.fIJCosThetaMin;
614     fIJCosThetaMax                = o.fIJCosThetaMax;
615     fIJNBinsTheta                 = o.fIJNBinsTheta;
616     fIJThetaMin                   = o.fIJThetaMin;
617     fIJThetaMax                   = o.fIJThetaMax;
618     fIJNBinsJt                    = o.fIJNBinsJt;
619     fIJJtMin                      = o.fIJJtMin;
620     fIJJtMax                      = o.fIJJtMax;
621     fDiJetNBinsJetInvMass         = o.fDiJetNBinsJetInvMass;
622     fDiJetJetInvMassMin           = o.fDiJetJetInvMassMin;
623     fDiJetJetInvMassMax           = o.fDiJetJetInvMassMax;
624     fDiJetNBinsJetPt              = o.fDiJetNBinsJetPt;
625     fDiJetJetPtMin                = o.fDiJetJetPtMin;
626     fDiJetJetPtMax                = o.fDiJetJetPtMax;
627     fDiJetNBinsPt                 = o.fDiJetNBinsPt;
628     fDiJetPtMin                   = o.fDiJetPtMin;
629     fDiJetPtMax                   = o.fDiJetPtMax;
630     fDiJetNBinsXi                 = o.fDiJetNBinsXi;
631     fDiJetXiMin                   = o.fDiJetXiMin;
632     fDiJetXiMax                   = o.fDiJetXiMax;
633     fDiJetNBinsZ                  = o.fDiJetNBinsZ;
634     fDiJetZMin                    = o.fDiJetZMin;
635     fDiJetZMax                    = o.fDiJetZMax;
636     fQADiJetNBinsInvMass          = o.fQADiJetNBinsInvMass;
637     fQADiJetInvMassMin            = o.fQADiJetInvMassMin;
638     fQADiJetInvMassMax            = o.fQADiJetInvMassMax;
639     fQADiJetNBinsJetPt            = o.fQADiJetNBinsJetPt;
640     fQADiJetJetPtMin              = o.fQADiJetJetPtMin;
641     fQADiJetJetPtMax              = o.fQADiJetJetPtMax;
642     fQADiJetNBinsDeltaPhi         = o.fQADiJetNBinsDeltaPhi;
643     fQADiJetDeltaPhiMin           = o.fQADiJetDeltaPhiMin;
644     fQADiJetDeltaPhiMax           = o.fQADiJetDeltaPhiMax;
645     fQADiJetNBinsDeltaEta         = o.fQADiJetNBinsDeltaEta;
646     fQADiJetDeltaEtaMin           = o.fQADiJetDeltaEtaMin;
647     fQADiJetDeltaEtaMax           = o.fQADiJetDeltaEtaMax;
648     fQADiJetNBinsDeltaPt          = o.fQADiJetNBinsDeltaPt;
649     fQADiJetDeltaPtMin            = o.fQADiJetDeltaPtMin;
650     fQADiJetDeltaPtMax            = o.fQADiJetDeltaPtMax;
651     fCommonHistList               = o.fCommonHistList;
652     fh1EvtSelection               = o.fh1EvtSelection;
653     fh1VertexNContributors        = o.fh1VertexNContributors;
654     fh1VertexZ                    = o.fh1VertexZ;
655     fh1EvtMult                    = o.fh1EvtMult;
656     fh1nRecJetsCuts               = o.fh1nRecJetsCuts;
657     fh1nGenJets                   = o.fh1nGenJets; 
658   }
659     
660   return *this;
661 }
662
663 //___________________________________________________________________________
664 AliAnalysisTaskFragmentationFunction::~AliAnalysisTaskFragmentationFunction()
665 {
666   // destructor
667   
668   if(fTracksRec)     delete fTracksRec;
669   if(fTracksRecCuts) delete fTracksRecCuts;
670   if(fTracksGen)     delete fTracksGen;
671   if(fJetsRec)       delete fJetsRec;
672   if(fJetsRecCuts)   delete fJetsRecCuts;
673   if(fJetsGen)       delete fJetsGen;
674
675   //  if(fDiJetBins)     delete fDiJetBins;
676
677 }
678
679 //______________________________________________________________________________________________________
680 AliAnalysisTaskFragmentationFunction::AliFragFuncHistos::AliFragFuncHistos(const char* name, 
681                                                          Int_t nJetPt, Float_t jetPtMin, Float_t jetPtMax,  
682                                                          Int_t nPt, Float_t ptMin, Float_t ptMax,
683                                                          Int_t nXi, Float_t xiMin, Float_t xiMax,
684                                                          Int_t nZ , Float_t zMin , Float_t zMax )
685   : TObject()
686   ,fNBinsJetPt(nJetPt)
687   ,fJetPtMin(jetPtMin)
688   ,fJetPtMax(jetPtMax)
689   ,fNBinsPt(nPt) 
690   ,fPtMin(ptMin)   
691   ,fPtMax(ptMax)   
692   ,fNBinsXi(nXi) 
693   ,fXiMin(xiMin)   
694   ,fXiMax(xiMax)   
695   ,fNBinsZ(nZ)  
696   ,fZMin(zMin)    
697   ,fZMax(zMax)    
698   ,fh2TrackPt(0)
699   ,fh2Xi(0)
700   ,fh2Z(0)
701   ,fh1JetPt(0)
702   ,fName(name)
703 {
704   // default constructor
705
706 }
707
708 //___________________________________________________________________________
709 AliAnalysisTaskFragmentationFunction::AliFragFuncHistos::AliFragFuncHistos(const AliFragFuncHistos& copy)
710   : TObject()
711   ,fNBinsJetPt(copy.fNBinsJetPt)
712   ,fJetPtMin(copy.fJetPtMin)
713   ,fJetPtMax(copy.fJetPtMax)
714   ,fNBinsPt(copy.fNBinsPt) 
715   ,fPtMin(copy.fPtMin)   
716   ,fPtMax(copy.fPtMax)   
717   ,fNBinsXi(copy.fNBinsXi) 
718   ,fXiMin(copy.fXiMin)   
719   ,fXiMax(copy.fXiMax)   
720   ,fNBinsZ(copy.fNBinsZ)  
721   ,fZMin(copy.fZMin)    
722   ,fZMax(copy.fZMax)    
723   ,fh2TrackPt(copy.fh2TrackPt)
724   ,fh2Xi(copy.fh2Xi)
725   ,fh2Z(copy.fh2Z)
726   ,fh1JetPt(copy.fh1JetPt)
727   ,fName(copy.fName)
728 {
729   // copy constructor
730 }
731
732 //_______________________________________________________________________________________________________________________________________________________________
733 AliAnalysisTaskFragmentationFunction::AliFragFuncHistos& AliAnalysisTaskFragmentationFunction::AliFragFuncHistos::operator=(const AliAnalysisTaskFragmentationFunction::AliFragFuncHistos& o)
734 {
735   // assignment
736   
737   if(this!=&o){
738     TObject::operator=(o);
739     fNBinsJetPt = o.fNBinsJetPt;
740     fJetPtMin   = o.fJetPtMin;
741     fJetPtMax   = o.fJetPtMax;
742     fNBinsPt    = o.fNBinsPt; 
743     fPtMin      = o.fPtMin;   
744     fPtMax      = o.fPtMax;   
745     fNBinsXi    = o.fNBinsXi; 
746     fXiMin      = o.fXiMin;   
747     fXiMax      = o.fXiMax;   
748     fNBinsZ     = o.fNBinsZ;  
749     fZMin       = o.fZMin;    
750     fZMax       = o.fZMax;    
751     fh2TrackPt  = o.fh2TrackPt;
752     fh2Xi       = o.fh2Xi;
753     fh2Z        = o.fh2Z;
754     fh1JetPt    = o.fh1JetPt;
755     fName       = o.fName;
756   }
757     
758   return *this;
759 }
760
761 //_________________________________________________________
762 AliAnalysisTaskFragmentationFunction::AliFragFuncHistos::~AliFragFuncHistos()
763 {
764   // destructor 
765
766   if(fh1JetPt)   delete fh1JetPt;
767   if(fh2TrackPt) delete fh2TrackPt;
768   if(fh2Xi)      delete fh2Xi;
769   if(fh2Z)       delete fh2Z;
770 }
771
772 //_________________________________________________________________
773 void AliAnalysisTaskFragmentationFunction::AliFragFuncHistos::DefineHistos()
774 {
775   // book FF histos
776
777   fh1JetPt   = new TH1F(Form("fh1FFJetPt%s", fName.Data()),"",fNBinsJetPt,fJetPtMin,fJetPtMax);
778   fh2TrackPt = new TH2F(Form("fh2FFTrackPt%s",fName.Data()),"",fNBinsJetPt, fJetPtMin, fJetPtMax,fNBinsPt, fPtMin, fPtMax);
779   fh2Xi      = new TH2F(Form("fh2FFXi%s",fName.Data()),"",fNBinsJetPt, fJetPtMin, fJetPtMax, fNBinsXi, fXiMin, fXiMax);
780   fh2Z       = new TH2F(Form("fh2FFZ%s",fName.Data()),"",fNBinsJetPt, fJetPtMin, fJetPtMax, fNBinsZ, fZMin, fZMax);
781
782   AliAnalysisTaskFragmentationFunction::SetProperties(fh1JetPt, "p_{T} [GeV/c]", "entries"); 
783   AliAnalysisTaskFragmentationFunction::SetProperties(fh2TrackPt,"jet p_{T} [GeV/c]","p_{T} [GeV/c]","entries");
784   AliAnalysisTaskFragmentationFunction::SetProperties(fh2Xi,"jet p_{T} [GeV/c]","#xi", "entries");
785   AliAnalysisTaskFragmentationFunction::SetProperties(fh2Z,"jet p_{T} [GeV/c]","z","entries");
786 }
787
788 //_______________________________________________________________________________________________________________
789 void AliAnalysisTaskFragmentationFunction::AliFragFuncHistos::FillFF(Float_t trackPt, Float_t jetPt, Bool_t incrementJetPt)
790 {
791   // fill FF
792  
793   if(incrementJetPt) fh1JetPt->Fill(jetPt);    
794   fh2TrackPt->Fill(jetPt,trackPt);
795   
796   // Added by me
797   Double_t z = 0.;
798   if(jetPt>0) z = trackPt / jetPt;
799   Double_t xi = 0;
800   if(z>0) xi = TMath::Log(1/z);
801   
802   fh2Xi->Fill(jetPt,xi);
803   fh2Z->Fill(jetPt,z);
804 }
805
806 //_________________________________________________________________________________
807 void AliAnalysisTaskFragmentationFunction::AliFragFuncHistos::AddToOutput(TList* list) const
808 {
809   // add histos to list
810
811   list->Add(fh1JetPt);
812   
813   list->Add(fh2TrackPt);
814   list->Add(fh2Xi);
815   list->Add(fh2Z);
816 }
817
818 //_________________________________________________________________________________________________________
819 AliAnalysisTaskFragmentationFunction::AliFragFuncQAJetHistos::AliFragFuncQAJetHistos(const char* name,
820                                                                Int_t nPt,   Float_t ptMin,   Float_t ptMax,
821                                                                Int_t nEta,  Float_t etaMin,  Float_t etaMax,
822                                                                Int_t nPhi,  Float_t phiMin,  Float_t phiMax)
823   : TObject()
824   ,fNBinsPt(nPt)
825   ,fPtMin(ptMin)
826   ,fPtMax(ptMax)
827   ,fNBinsEta(nEta)
828   ,fEtaMin(etaMin)
829   ,fEtaMax(etaMax)
830   ,fNBinsPhi(nPhi)
831   ,fPhiMin(phiMin)
832   ,fPhiMax(phiMax)
833   ,fh2EtaPhi(0)
834   ,fh1Pt(0)
835   ,fName(name)
836 {
837   // default constructor
838 }
839
840 //____________________________________________________________________________________
841 AliAnalysisTaskFragmentationFunction::AliFragFuncQAJetHistos::AliFragFuncQAJetHistos(const AliFragFuncQAJetHistos& copy)
842   : TObject()
843   ,fNBinsPt(copy.fNBinsPt)
844   ,fPtMin(copy.fPtMin)
845   ,fPtMax(copy.fPtMax)
846   ,fNBinsEta(copy.fNBinsEta)
847   ,fEtaMin(copy.fEtaMin)
848   ,fEtaMax(copy.fEtaMax)
849   ,fNBinsPhi(copy.fNBinsPhi)
850   ,fPhiMin(copy.fPhiMin)
851   ,fPhiMax(copy.fPhiMax)
852   ,fh2EtaPhi(copy.fh2EtaPhi)
853   ,fh1Pt(copy.fh1Pt)
854   ,fName(copy.fName)
855 {
856   // copy constructor
857 }
858
859 //________________________________________________________________________________________________________________________________________________________________________
860 AliAnalysisTaskFragmentationFunction::AliFragFuncQAJetHistos& AliAnalysisTaskFragmentationFunction::AliFragFuncQAJetHistos::operator=(const AliAnalysisTaskFragmentationFunction::AliFragFuncQAJetHistos& o)
861 {
862   // assignment
863   
864   if(this!=&o){
865     TObject::operator=(o);
866     fNBinsPt  = o.fNBinsPt;
867     fPtMin    = o.fPtMin;
868     fPtMax    = o.fPtMax;
869     fNBinsEta = o.fNBinsEta;
870     fEtaMin   = o.fEtaMin;
871     fEtaMax   = o.fEtaMax;
872     fNBinsPhi = o.fNBinsPhi;
873     fPhiMin   = o.fPhiMin;
874     fPhiMax   = o.fPhiMax;
875     fh2EtaPhi = o.fh2EtaPhi;
876     fh1Pt     = o.fh1Pt;
877     fName     = o.fName;
878   }
879   
880   return *this;
881 }
882
883 //______________________________________________________________
884 AliAnalysisTaskFragmentationFunction::AliFragFuncQAJetHistos::~AliFragFuncQAJetHistos()
885 {
886   // destructor 
887   
888   if(fh2EtaPhi) delete fh2EtaPhi;
889   if(fh1Pt)     delete fh1Pt;
890 }
891
892 //____________________________________________________________________
893 void AliAnalysisTaskFragmentationFunction::AliFragFuncQAJetHistos::DefineHistos()
894 {
895   // book jet QA histos
896
897   fh2EtaPhi  = new TH2F(Form("fh2JetQAEtaPhi%s", fName.Data()), Form("%s: #eta - #phi distribution", fName.Data()), fNBinsEta, fEtaMin, fEtaMax, fNBinsPhi, fPhiMin, fPhiMax);
898   fh1Pt      = new TH1F(Form("fh1JetQAPt%s", fName.Data()), Form("%s: p_{T} distribution", fName.Data()), fNBinsPt, fPtMin, fPtMax);
899         
900   AliAnalysisTaskFragmentationFunction::SetProperties(fh2EtaPhi, "#eta", "#phi"); 
901   AliAnalysisTaskFragmentationFunction::SetProperties(fh1Pt, "p_{T} [GeV/c]", "entries");
902 }
903
904 //____________________________________________________________________________________________________
905 void AliAnalysisTaskFragmentationFunction::AliFragFuncQAJetHistos::FillJetQA(Float_t eta, Float_t phi, Float_t pt)
906 {
907   // fill jet QA histos 
908
909   fh2EtaPhi->Fill( eta, phi);
910   fh1Pt->Fill( pt );
911 }
912
913 //____________________________________________________________________________________
914 void AliAnalysisTaskFragmentationFunction::AliFragFuncQAJetHistos::AddToOutput(TList* list) const 
915 {
916   // add histos to list
917
918   list->Add(fh2EtaPhi);
919   list->Add(fh1Pt);
920 }
921
922 //___________________________________________________________________________________________________________
923 AliAnalysisTaskFragmentationFunction::AliFragFuncQATrackHistos::AliFragFuncQATrackHistos(const char* name,
924                                                                    Int_t nPt, Float_t ptMin, Float_t ptMax,
925                                                                    Int_t nEta, Float_t etaMin, Float_t etaMax,
926                                                                    Int_t nPhi, Float_t phiMin, Float_t phiMax,
927                                                                    Float_t ptThresh) 
928   : TObject()
929   ,fNBinsPt(nPt)
930   ,fPtMin(ptMin)
931   ,fPtMax(ptMax)
932   ,fNBinsEta(nEta)
933   ,fEtaMin(etaMin)
934   ,fEtaMax(etaMax)
935   ,fNBinsPhi(nPhi)
936   ,fPhiMin(phiMin)
937   ,fPhiMax(phiMax)
938   ,fHighPtThreshold(ptThresh)
939   ,fh2EtaPhi(0)
940   ,fh1Pt(0)
941   ,fh2HighPtEtaPhi(0)
942   ,fName(name)
943 {
944   // default constructor
945 }
946
947 //__________________________________________________________________________________________
948 AliAnalysisTaskFragmentationFunction::AliFragFuncQATrackHistos::AliFragFuncQATrackHistos(const AliFragFuncQATrackHistos& copy)
949   : TObject()
950   ,fNBinsPt(copy.fNBinsPt)
951   ,fPtMin(copy.fPtMin)
952   ,fPtMax(copy.fPtMax)
953   ,fNBinsEta(copy.fNBinsEta)
954   ,fEtaMin(copy.fEtaMin)
955   ,fEtaMax(copy.fEtaMax)
956   ,fNBinsPhi(copy.fNBinsPhi)
957   ,fPhiMin(copy.fPhiMin)
958   ,fPhiMax(copy.fPhiMax)
959   ,fHighPtThreshold(copy.fHighPtThreshold)
960   ,fh2EtaPhi(copy.fh2EtaPhi)
961   ,fh1Pt(copy.fh1Pt)
962   ,fh2HighPtEtaPhi(copy.fh2HighPtEtaPhi)
963   ,fName(copy.fName)
964 {
965   // copy constructor
966 }
967
968 // _____________________________________________________________________________________________________________________________________________________________________________
969 AliAnalysisTaskFragmentationFunction::AliFragFuncQATrackHistos& AliAnalysisTaskFragmentationFunction::AliFragFuncQATrackHistos::operator=(const AliAnalysisTaskFragmentationFunction::AliFragFuncQATrackHistos& o)
970 {
971   // assignment
972   
973   if(this!=&o){
974     TObject::operator=(o);
975     fNBinsPt         = o.fNBinsPt;
976     fPtMin           = o.fPtMin;
977     fPtMax           = o.fPtMax;
978     fNBinsEta        = o.fNBinsEta;
979     fEtaMin          = o.fEtaMin;
980     fEtaMax          = o.fEtaMax;
981     fNBinsPhi        = o.fNBinsPhi;
982     fPhiMin          = o.fPhiMin;
983     fPhiMax          = o.fPhiMax;
984     fHighPtThreshold = o.fHighPtThreshold;
985     fh2EtaPhi        = o.fh2EtaPhi;
986     fh1Pt            = o.fh1Pt;
987     fh2HighPtEtaPhi  = o.fh2HighPtEtaPhi;
988     fName            = o.fName;
989   }
990   
991   return *this;
992 }
993
994 //___________________________________________________________________
995 AliAnalysisTaskFragmentationFunction::AliFragFuncQATrackHistos::~AliFragFuncQATrackHistos()
996 {
997   // destructor 
998   
999   if(fh2EtaPhi)       delete fh2EtaPhi;
1000   if(fh2HighPtEtaPhi) delete fh2HighPtEtaPhi;
1001   if(fh1Pt)           delete fh1Pt;
1002 }
1003
1004 //______________________________________________________________________
1005 void AliAnalysisTaskFragmentationFunction::AliFragFuncQATrackHistos::DefineHistos()
1006 {
1007   // book track QA histos
1008
1009   fh2EtaPhi       = new TH2F(Form("fh2TrackQAEtaPhi%s", fName.Data()), Form("%s: #eta - #phi distribution", fName.Data()), fNBinsEta, fEtaMin, fEtaMax, fNBinsPhi, fPhiMin, fPhiMax);
1010   fh2HighPtEtaPhi = new TH2F(Form("fh2TrackQAHighPtEtaPhi%s", fName.Data()), Form("%s: #eta - #phi distribution for high-p_{T}", fName.Data()), fNBinsEta, fEtaMin, fEtaMax, fNBinsPhi, fPhiMin, fPhiMax);
1011   fh1Pt           = new TH1F(Form("fh1TrackQAPt%s", fName.Data()), Form("%s: p_{T} distribution", fName.Data()), fNBinsPt, fPtMin, fPtMax);
1012   
1013   AliAnalysisTaskFragmentationFunction::SetProperties(fh2EtaPhi, "#eta", "#phi"); 
1014   AliAnalysisTaskFragmentationFunction::SetProperties(fh2HighPtEtaPhi, "#eta", "#phi");
1015   AliAnalysisTaskFragmentationFunction::SetProperties(fh1Pt, "p_{T} [GeV/c]", "entries");
1016 }
1017
1018 //________________________________________________________________________________________________________
1019 void AliAnalysisTaskFragmentationFunction::AliFragFuncQATrackHistos::FillTrackQA(Float_t eta, Float_t phi, Float_t pt)
1020 {
1021   // fill track QA histos
1022     
1023   fh2EtaPhi->Fill( eta, phi);
1024   if(pt > fHighPtThreshold) fh2HighPtEtaPhi->Fill( eta, phi);
1025   fh1Pt->Fill( pt );    
1026 }
1027
1028 //______________________________________________________________________________________
1029 void AliAnalysisTaskFragmentationFunction::AliFragFuncQATrackHistos::AddToOutput(TList* list) const
1030 {
1031   // add histos to list
1032
1033   list->Add(fh2EtaPhi);
1034   list->Add(fh2HighPtEtaPhi);
1035   list->Add(fh1Pt);
1036 }
1037
1038 //______________________________________________________________________________________________________
1039 AliAnalysisTaskFragmentationFunction::AliFragFuncIntraJetHistos::AliFragFuncIntraJetHistos(const char* name, 
1040                                                          Int_t nJetPt, Float_t jetPtMin, Float_t jetPtMax,  
1041                                                          Int_t nPt, Float_t ptMin, Float_t ptMax,
1042                                                          Int_t nZ , Float_t zMin , Float_t zMax,
1043                                                          Int_t nCosTheta , Float_t costhetaMin , Float_t costhetaMax,
1044                                                          Int_t nTheta , Float_t thetaMin , Float_t thetaMax,
1045                                                          Int_t nJt , Float_t jtMin , Float_t jtMax)
1046   : TObject()
1047   ,fNBinsJetPt(nJetPt)
1048   ,fJetPtMin(jetPtMin)
1049   ,fJetPtMax(jetPtMax)
1050   ,fNBinsPt(nPt) 
1051   ,fPtMin(ptMin)   
1052   ,fPtMax(ptMax)   
1053   ,fNBinsZ(nZ) 
1054   ,fZMin(zMin)   
1055   ,fZMax(zMax)   
1056   ,fNBinsJt(nJt)
1057   ,fJtMin(jtMin)
1058   ,fJtMax(jtMax)
1059   ,fNBinsTheta(nTheta)
1060   ,fThetaMin(thetaMin)
1061   ,fThetaMax(thetaMax)
1062   ,fNBinsCosTheta(nCosTheta)
1063   ,fCosThetaMin(costhetaMin)
1064   ,fCosThetaMax(costhetaMax)
1065   ,fh2Theta(0)
1066   ,fh2CosTheta(0)
1067   ,fh2Jt(0)
1068   ,fh2PtvsZ(0)
1069   ,fhnIntraJet(0)
1070   ,fnDim(6)
1071   ,fName(name)
1072 {
1073   // default constructor
1074
1075 }
1076
1077 //___________________________________________________________________________
1078 AliAnalysisTaskFragmentationFunction::AliFragFuncIntraJetHistos::AliFragFuncIntraJetHistos(const AliFragFuncIntraJetHistos& copy)
1079   : TObject()
1080   ,fNBinsJetPt(copy.fNBinsJetPt)
1081   ,fJetPtMin(copy.fJetPtMin)
1082   ,fJetPtMax(copy.fJetPtMax)
1083   ,fNBinsPt(copy.fNBinsPt) 
1084   ,fPtMin(copy.fPtMin)   
1085   ,fPtMax(copy.fPtMax)   
1086   ,fNBinsZ(copy.fNBinsZ) 
1087   ,fZMin(copy.fZMin)   
1088   ,fZMax(copy.fZMax)   
1089   ,fNBinsJt(copy.fNBinsJt)
1090   ,fJtMin(copy.fJtMin)
1091   ,fJtMax(copy.fJtMax)
1092   ,fNBinsTheta(copy.fNBinsTheta)
1093   ,fThetaMin(copy.fThetaMin)
1094   ,fThetaMax(copy.fThetaMax)
1095   ,fNBinsCosTheta(copy.fNBinsCosTheta)
1096   ,fCosThetaMin(copy.fCosThetaMin)
1097   ,fCosThetaMax(copy.fCosThetaMax)
1098   ,fh2Theta(copy.fh2Theta)
1099   ,fh2CosTheta(copy.fh2CosTheta)
1100   ,fh2Jt(copy.fh2Jt)
1101   ,fh2PtvsZ(copy.fh2PtvsZ)
1102   ,fhnIntraJet(copy.fhnIntraJet)
1103   ,fnDim(copy.fnDim)
1104   ,fName(copy.fName)
1105 {
1106   // copy constructor
1107 }
1108
1109 //_______________________________________________________________________________________________________________________________________________________________
1110 AliAnalysisTaskFragmentationFunction::AliFragFuncIntraJetHistos& AliAnalysisTaskFragmentationFunction::AliFragFuncIntraJetHistos::operator=(const AliAnalysisTaskFragmentationFunction::AliFragFuncIntraJetHistos& o)
1111 {
1112   // assignment
1113   
1114   if(this!=&o){
1115     TObject::operator=(o);
1116     fNBinsJetPt       = o.fNBinsJetPt;
1117     fJetPtMin         = o.fJetPtMin;
1118     fJetPtMax         = o.fJetPtMax;
1119     fNBinsPt          = o.fNBinsPt; 
1120     fPtMin            = o.fPtMin;   
1121     fPtMax            = o.fPtMax;   
1122     fNBinsZ           = o.fNBinsZ; 
1123     fZMin             = o.fZMin;   
1124     fZMax             = o.fZMax;   
1125     fNBinsJt          = o.fNBinsJt;
1126     fJtMin            = o.fJtMin;
1127     fJtMax            = o.fJtMax;
1128     fNBinsTheta       = o.fNBinsTheta;
1129     fThetaMin         = o.fThetaMin;
1130     fThetaMax         = o.fThetaMax;
1131     fNBinsCosTheta    = o.fNBinsCosTheta;
1132     fCosThetaMin      = o.fCosThetaMin;
1133     fCosThetaMax      = o.fCosThetaMax;
1134     fh2Theta          = o.fh2Theta;
1135     fh2CosTheta       = o.fh2CosTheta;
1136     fh2Jt             = o.fh2Jt;
1137     fh2PtvsZ          = o.fh2PtvsZ;
1138     fhnIntraJet       = o.fhnIntraJet;
1139     fnDim             = o.fnDim;
1140     fName             = o.fName;
1141   }
1142     
1143   return *this;
1144 }
1145
1146 //_________________________________________________________
1147 AliAnalysisTaskFragmentationFunction::AliFragFuncIntraJetHistos::~AliFragFuncIntraJetHistos()
1148 {
1149   // destructor 
1150
1151
1152   if(fh2Theta)          delete fh2Theta;
1153   if(fh2CosTheta)       delete fh2CosTheta;
1154   if(fh2Jt)             delete fh2Jt;
1155   if(fh2PtvsZ)          delete fh2PtvsZ;
1156   if(fhnIntraJet)       delete fhnIntraJet;
1157
1158 }
1159
1160 //_________________________________________________________________
1161 void AliAnalysisTaskFragmentationFunction::AliFragFuncIntraJetHistos::DefineHistos()
1162 {
1163   // book FF histos
1164
1165   fh2Theta    = new TH2F(Form("fh2IJTheta%s",fName.Data()),"",fNBinsJetPt, fJetPtMin, fJetPtMax,fNBinsTheta, fThetaMin, fThetaMax);
1166   fh2CosTheta = new TH2F(Form("fh2IJcosTheta%s",fName.Data()),"",fNBinsJetPt, fJetPtMin, fJetPtMax,fNBinsCosTheta, fCosThetaMin, fCosThetaMax);
1167   fh2Jt       = new TH2F(Form("fh2IJJt%s",fName.Data()),"",fNBinsJetPt, fJetPtMin, fJetPtMax, fNBinsJt, fJtMin, fJtMax);
1168
1169   // Create 3D histograms
1170   Int_t    *iBin = new Int_t[fnDim];
1171   Double_t *min  = new Double_t[fnDim];
1172   Double_t *max  = new Double_t[fnDim];
1173
1174   iBin[0] = fNBinsJetPt; iBin[1] = fNBinsTheta; iBin[2] = fNBinsCosTheta; iBin[3] = fNBinsJt; iBin[4] = fNBinsZ; iBin[5] = fNBinsPt;
1175   min[0]  = fJetPtMin; min[1] = fThetaMin; min[2] = fCosThetaMin; min[3] = fJtMin; min[4] = fZMin; min[5] = fPtMin; 
1176   max[0]  = fJetPtMax; max[1] = fThetaMax; max[2] = fCosThetaMax; max[3] = fJtMax; max[4] = fZMax; max[5] = fPtMax;
1177
1178   const char* title = Form("fhnIntraJetPart%s",fName.Data());
1179   const char* comment = "THnSparseF p_{T} jet [GeV/c] : #Theta : cos(#Theta) : j_{T} : Z : p_{T} part [GeV/c]";
1180   fhnIntraJet = new THnSparseF(title,comment,fnDim,iBin,min,max);
1181
1182   const char** axisTitle = new const char*[fnDim];
1183   axisTitle[0] = "p_{T}^{jet} [GeV/c]";
1184   axisTitle[1] = "#Theta";
1185   axisTitle[2] = "Cos(#Theta)";
1186   axisTitle[3] = "j_{T} [GeV]";
1187   axisTitle[4] = "z = p_{T}^{had}/p_{T}^{jet}";
1188   axisTitle[5] = "p_{T}^{had} [GeV/c]";
1189   
1190   AliAnalysisTaskFragmentationFunction::SetProperties(fh2Theta,"jet p_{T} [GeV/c]","#Theta","entries");
1191   AliAnalysisTaskFragmentationFunction::SetProperties(fh2CosTheta,"jet p_{T} [GeV/c]","cos(#Theta)", "entries");
1192   AliAnalysisTaskFragmentationFunction::SetProperties(fh2Jt,"jet p_{T} [GeV/c]","j_{T}","entries");
1193   AliAnalysisTaskFragmentationFunction::SetProperties(fhnIntraJet,fnDim,axisTitle);
1194   delete[] iBin;
1195   delete[] min;
1196   delete[] max;
1197   delete[] axisTitle;
1198 }
1199
1200 //_______________________________________________________________________________________________________________
1201 void AliAnalysisTaskFragmentationFunction::AliFragFuncIntraJetHistos::FillIntraJet(TLorentzVector* trackV, TLorentzVector* jetV)
1202 {
1203   // fill IntraJet histos
1204  
1205   Float_t cosTheta = 0.; Float_t theta = 0.; 
1206   Float_t jt = 0.; Float_t z = 0.; 
1207   // For Theta distribution
1208   Float_t pxT  = trackV->Px();
1209   Float_t pyT  = trackV->Py();
1210   Float_t pzT  = trackV->Pz();
1211   Float_t ptT  = trackV->Pt();
1212   Float_t pT   = trackV->P();
1213   Float_t etaT = trackV->Eta();
1214   Float_t phiT = trackV->Phi(); // Check the value returned
1215   Float_t pxJ = jetV->Px();
1216   Float_t pyJ = jetV->Py();
1217   Float_t pzJ = jetV->Pz();
1218   Float_t ptJ = jetV->Pt();
1219   Float_t pJ  = jetV->P();
1220
1221   // Compute z
1222   z = (Float_t)(ptT/ptJ);
1223
1224   // Compute theta
1225   cosTheta = (pxT*pxJ+pyT*pyJ+pzT*pzJ)/(pT*pJ);
1226   theta = TMath::ACos(cosTheta);
1227
1228   // Compute jt
1229   TVector3 trackP; TVector3 jetP;
1230   jetP[0] = pxJ;
1231   jetP[1] = pyJ;
1232   jetP[2] = pzJ;
1233   trackP.SetPtEtaPhi(ptT,etaT,phiT);
1234   jt = TMath::Sin(trackP.Angle(jetP))*trackP.Mag();
1235
1236   // Fill histos and THnSparse
1237   fh2CosTheta->Fill(ptJ,cosTheta);
1238   fh2Theta->Fill(ptJ,theta);
1239   fh2Jt->Fill(ptJ,jt);
1240
1241   // Fill THnSparse
1242   Double_t *content = new Double_t[fnDim];
1243   content[0]= ptJ; content[1] = theta; content[2] = cosTheta; content[3] = jt; content[4] = z; content[5] = ptT; 
1244
1245   fhnIntraJet->Fill(content);
1246
1247   delete content;
1248
1249 }
1250
1251 //______________________________________________________________________________________________________
1252 AliAnalysisTaskFragmentationFunction::AliFragFuncDiJetHistos::AliFragFuncDiJetHistos(const char* name, Int_t kindSlices,
1253                                                          Int_t nJetInvMass, Float_t jetInvMassMin, Float_t jetInvMassMax,  
1254                                                          Int_t nJetPt, Float_t jetPtMin, Float_t jetPtMax,  
1255                                                          Int_t nPt, Float_t ptMin, Float_t ptMax,
1256                                                          Int_t nXi, Float_t xiMin, Float_t xiMax,
1257                                                          Int_t nZ , Float_t zMin , Float_t zMax)
1258   : TObject()
1259   ,fKindSlices(kindSlices)
1260   ,fNBinsJetInvMass(nJetInvMass)
1261   ,fJetInvMassMin(jetInvMassMin)
1262   ,fJetInvMassMax(jetInvMassMax)
1263   ,fNBinsJetPt(nJetPt)
1264   ,fJetPtMin(jetPtMin)
1265   ,fJetPtMax(jetPtMax)
1266   ,fNBinsPt(nPt) 
1267   ,fPtMin(ptMin)   
1268   ,fPtMax(ptMax)   
1269   ,fNBinsXi(nXi) 
1270   ,fXiMin(xiMin)   
1271   ,fXiMax(xiMax)   
1272   ,fNBinsZ(nZ)  
1273   ,fZMin(zMin)    
1274   ,fZMax(zMax)
1275   ,fh2TrackPtJet1(0)
1276   ,fh2TrackPtJet2(0)
1277   ,fh2TrackPtJet(0)
1278   ,fh1Jet1Pt(0)
1279   ,fh1Jet2Pt(0)
1280   ,fh1JetPt(0)
1281   ,fh2Xi1(0)
1282   ,fh2Xi2(0)
1283   ,fh2Xi(0)
1284   ,fh2Z1(0)
1285   ,fh2Z2(0)
1286   ,fh2Z(0)
1287   ,fh2Pt1(0)
1288   ,fh2Pt2(0)
1289   ,fh2Pt(0)
1290   ,fName(name)
1291 {
1292   // default constructor
1293
1294 }
1295
1296 //______________________________________________________________________________________________________
1297 AliAnalysisTaskFragmentationFunction::AliFragFuncDiJetHistos::AliFragFuncDiJetHistos(const AliFragFuncDiJetHistos& copy)
1298   : TObject()
1299   ,fKindSlices(copy.fKindSlices)
1300   ,fNBinsJetInvMass(copy.fNBinsJetInvMass)
1301   ,fJetInvMassMin(copy.fJetInvMassMin)
1302   ,fJetInvMassMax(copy.fJetInvMassMax)
1303   ,fNBinsJetPt(copy.fNBinsJetPt)
1304   ,fJetPtMin(copy.fJetPtMin)
1305   ,fJetPtMax(copy.fJetPtMax)
1306   ,fNBinsPt(copy.fNBinsPt) 
1307   ,fPtMin(copy.fPtMin)   
1308   ,fPtMax(copy.fPtMax)   
1309   ,fNBinsXi(copy.fNBinsXi) 
1310   ,fXiMin(copy.fXiMin)   
1311   ,fXiMax(copy.fXiMax)   
1312   ,fNBinsZ(copy.fNBinsZ)  
1313   ,fZMin(copy.fZMin)    
1314   ,fZMax(copy.fZMax)
1315   ,fh2TrackPtJet1(copy.fh2TrackPtJet1)
1316   ,fh2TrackPtJet2(copy.fh2TrackPtJet2)
1317   ,fh2TrackPtJet(copy.fh2TrackPtJet)
1318   ,fh1Jet1Pt(copy.fh1Jet1Pt)
1319   ,fh1Jet2Pt(copy.fh1Jet2Pt)
1320   ,fh1JetPt(copy.fh1JetPt)
1321   ,fh2Xi1(copy.fh2Xi1)
1322   ,fh2Xi2(copy.fh2Xi2)
1323   ,fh2Xi(copy.fh2Xi2)
1324   ,fh2Z1(copy.fh2Z1)
1325   ,fh2Z2(copy.fh2Z2)
1326   ,fh2Z(copy.fh2Z)
1327   ,fh2Pt1(copy.fh2Pt1)
1328   ,fh2Pt2(copy.fh2Pt2)
1329   ,fh2Pt(copy.fh2Pt)
1330   ,fName(copy.fName)
1331 {
1332   // default constructor
1333
1334 }
1335
1336 //_______________________________________________________________________________________________________________________________________________________________
1337 AliAnalysisTaskFragmentationFunction::AliFragFuncDiJetHistos& AliAnalysisTaskFragmentationFunction::AliFragFuncDiJetHistos::operator=(const AliAnalysisTaskFragmentationFunction::AliFragFuncDiJetHistos& o)
1338 {
1339   // assignment
1340   
1341   if(this!=&o){
1342     TObject::operator=(o);
1343     fKindSlices      = o.fKindSlices;
1344     fNBinsJetInvMass = o.fNBinsJetInvMass;
1345     fJetInvMassMin   = o.fJetInvMassMin;
1346     fJetInvMassMax   = o.fJetInvMassMax;
1347     fNBinsJetPt      = o.fNBinsJetPt;
1348     fJetPtMin        = o.fJetPtMin;
1349     fJetPtMax        = o.fJetPtMax;
1350     fNBinsPt         = o.fNBinsPt; 
1351     fPtMin           = o.fPtMin;   
1352     fPtMax           = o.fPtMax;   
1353     fNBinsXi         = o.fNBinsXi; 
1354     fXiMin           = o.fXiMin;   
1355     fXiMax           = o.fXiMax;   
1356     fNBinsZ          = o.fNBinsZ;  
1357     fZMin            = o.fZMin;    
1358     fZMax            = o.fZMax;   
1359     fh2TrackPtJet1   = o.fh2TrackPtJet1;
1360     fh2TrackPtJet2   = o.fh2TrackPtJet2;
1361     fh2TrackPtJet    = o.fh2TrackPtJet;
1362     fh1Jet1Pt        = o.fh1Jet1Pt;
1363     fh1Jet2Pt        = o.fh1Jet2Pt;
1364     fh1JetPt         = o.fh1JetPt;
1365     fh2Xi1           = o.fh2Xi1;
1366     fh2Xi2           = o.fh2Xi2;
1367     fh2Xi            = o.fh2Xi;
1368     fh2Z1            = o.fh2Z1;
1369     fh2Z2            = o.fh2Z2;
1370     fh2Z             = o.fh2Z;
1371     fh2Pt1           = o.fh2Pt1;
1372     fh2Pt2           = o.fh2Pt2;
1373     fh2Pt            = o.fh2Pt;
1374     fName            = o.fName;
1375   }
1376     
1377   return *this;
1378 }
1379
1380 //_________________________________________________________
1381 AliAnalysisTaskFragmentationFunction::AliFragFuncDiJetHistos::~AliFragFuncDiJetHistos()
1382 {
1383   // destructor 
1384
1385   if(fh2TrackPtJet1) delete fh2TrackPtJet1;
1386   if(fh2TrackPtJet2) delete fh2TrackPtJet2;
1387   if(fh2TrackPtJet ) delete fh2TrackPtJet;
1388   if(fh1Jet1Pt)      delete fh1Jet1Pt;
1389   if(fh1Jet2Pt)      delete fh1Jet2Pt;
1390   if(fh1JetPt)       delete fh1JetPt;
1391   if(fh2Xi1)         delete fh2Xi1;
1392   if(fh2Xi2)         delete fh2Xi2;
1393   if(fh2Xi)          delete fh2Xi;
1394   if(fh2Z1)          delete fh2Z1;
1395   if(fh2Z2)          delete fh2Z2;
1396   if(fh2Z)           delete fh2Z;
1397   if(fh2Pt1)         delete fh2Pt1;
1398   if(fh2Pt2)         delete fh2Pt2;
1399   if(fh2Pt)          delete fh2Pt;
1400 }
1401
1402 //________________________________________________________________________
1403 void AliAnalysisTaskFragmentationFunction::AliFragFuncDiJetHistos::DefineDiJetHistos()
1404 {
1405
1406   Int_t nBins = 0;
1407   Double_t min = 0.;
1408   Double_t max = 0.;
1409   const char *xaxis = "";
1410   if(fKindSlices == 1)
1411     {
1412       nBins = fNBinsJetInvMass;
1413       min   = fJetInvMassMin;
1414       max   = fJetInvMassMax;
1415       xaxis = "M_{JJ} [GeV]";
1416     }
1417   if(fKindSlices == 2 || fKindSlices == 3)
1418     {
1419       nBins = fNBinsJetPt;
1420       min   = fJetPtMin;
1421       max   = fJetPtMax;
1422       if(fKindSlices == 2) xaxis = "E_{Tmean} [GeV]";
1423       if(fKindSlices == 3) xaxis ="leading jet p_{T} [GeV/c]";
1424     }
1425   
1426   fh1Jet1Pt      = new TH1F(Form("fh1DJJet1Pt%s", fName.Data()), "", fNBinsJetPt, fJetPtMin, fJetPtMax);
1427   fh1Jet2Pt      = new TH1F(Form("fh1DJJet2Pt%s", fName.Data()), "", fNBinsJetPt, fJetPtMin, fJetPtMax);
1428   fh1JetPt       = new TH1F(Form("fh1DJJetPt%s",  fName.Data()), "", fNBinsJetPt, fJetPtMin, fJetPtMax);
1429   
1430   fh2TrackPtJet1 = new TH2F(Form("fh2DJTrackPtJet1%s", fName.Data()), "",nBins, min, max, fNBinsPt, fPtMin, fPtMax);
1431   fh2TrackPtJet2 = new TH2F(Form("fh2DJTrackPtJet2%s", fName.Data()), "",nBins, min, max, fNBinsPt, fPtMin, fPtMax);
1432   fh2TrackPtJet  = new TH2F(Form("fh2DJTrackPtJet%s", fName.Data()),  "",nBins, min, max, fNBinsPt, fPtMin, fPtMax);
1433   
1434   fh2Xi1         = new TH2F(Form("fh2DJXi1%s", fName.Data()), "",nBins, min, max, fNBinsXi, fXiMin, fXiMax);
1435   fh2Xi2         = new TH2F(Form("fh2DJXi2%s", fName.Data()), "",nBins, min, max, fNBinsXi, fXiMin, fXiMax);
1436   fh2Xi          = new TH2F(Form("fh2DJXi%s", fName.Data()),  "",nBins, min, max, fNBinsXi, fXiMin, fXiMax);
1437   
1438   fh2Z1          = new TH2F(Form("fh2DJZ1%s", fName.Data()), "",nBins, min, max, fNBinsZ, fZMin, fZMax);
1439   fh2Z2          = new TH2F(Form("fh2DJZ2%s", fName.Data()), "",nBins, min, max, fNBinsZ, fZMin, fZMax);
1440   fh2Z           = new TH2F(Form("fh2DJZ%s", fName.Data()),  "",nBins, min, max, fNBinsZ, fZMin, fZMax);
1441   
1442   fh2Pt1         = new TH2F(Form("fh2DJPt1%s", fName.Data()), "",nBins, min, max, fNBinsPt, fPtMin, fPtMax);
1443   fh2Pt2         = new TH2F(Form("fh2DJPt2%s", fName.Data()), "",nBins, min, max, fNBinsPt, fPtMin, fPtMax);
1444   fh2Pt          = new TH2F(Form("fh2DJPtZ%s", fName.Data()),  "",nBins, min, max, fNBinsPt, fPtMin, fPtMax);
1445       
1446   AliAnalysisTaskFragmentationFunction::SetProperties(fh1Jet1Pt, "p_{T} [GeV/c]", "entries");
1447   AliAnalysisTaskFragmentationFunction::SetProperties(fh1Jet2Pt, "p_{T} [GeV/c]", "entries");
1448   AliAnalysisTaskFragmentationFunction::SetProperties(fh1JetPt, "p_{T} [GeV/c]", "entries");
1449
1450   AliAnalysisTaskFragmentationFunction::SetProperties(fh2TrackPtJet1, xaxis, "p_{T} [GeV/c]", "Entries");
1451   AliAnalysisTaskFragmentationFunction::SetProperties(fh2TrackPtJet2, xaxis, "p_{T} [GeV/c]", "Entries");
1452   AliAnalysisTaskFragmentationFunction::SetProperties(fh2TrackPtJet, xaxis, "p_{T} [GeV/c]", "Entries");
1453   AliAnalysisTaskFragmentationFunction::SetProperties(fh2Xi1, xaxis, "#xi", "Entries");
1454   AliAnalysisTaskFragmentationFunction::SetProperties(fh2Xi2, xaxis, "#xi", "Entries");
1455   AliAnalysisTaskFragmentationFunction::SetProperties(fh2Xi, xaxis, "#xi", "Entries");
1456   AliAnalysisTaskFragmentationFunction::SetProperties(fh2Z1, xaxis, "z", "Entries");
1457   AliAnalysisTaskFragmentationFunction::SetProperties(fh2Z2, xaxis, "z", "Entries");
1458   AliAnalysisTaskFragmentationFunction::SetProperties(fh2Z, xaxis, "z", "Entries");
1459   AliAnalysisTaskFragmentationFunction::SetProperties(fh2Pt1, xaxis, "p_{T} [GeV/c]", "Entries");
1460   AliAnalysisTaskFragmentationFunction::SetProperties(fh2Pt2, xaxis, "p_{T} [GeV/c]", "Entries");
1461   AliAnalysisTaskFragmentationFunction::SetProperties(fh2Pt, xaxis, "p_{T} [GeV/c]", "Entries");
1462   
1463 }
1464
1465 //________________________________________________________________________
1466 void AliAnalysisTaskFragmentationFunction::AliFragFuncDiJetHistos::FillDiJetFF(Int_t jetType, Float_t trackPt, Float_t jetPt, Double_t jetBin, Bool_t incrementJetPt)
1467 {
1468   if(jetType == 0)
1469     {
1470       if(incrementJetPt) fh1JetPt->Fill(jetPt);  
1471       
1472       fh2TrackPtJet->Fill(jetBin, trackPt);
1473       
1474       Double_t z = trackPt / jetPt;
1475       Double_t xi = 0;
1476       if(z!=0) xi = TMath::Log(1/z);
1477       
1478       fh2Xi->Fill(jetBin, xi);
1479       fh2Z->Fill(jetBin, z);
1480     }
1481   if(jetType == 1)
1482     {
1483       if(incrementJetPt) fh1Jet1Pt->Fill(jetPt);
1484       
1485       fh2TrackPtJet1->Fill(jetBin, trackPt);
1486       
1487       Double_t z = trackPt / jetPt;
1488       Double_t xi = 0;
1489       if(z!=0) xi = TMath::Log(1/z);
1490       
1491       fh2Xi1->Fill(jetBin, xi);
1492       fh2Z1->Fill(jetBin, z);
1493     }
1494   if(jetType == 2)
1495     {
1496       if(incrementJetPt) fh1Jet2Pt->Fill(jetPt);
1497       
1498       fh2TrackPtJet2->Fill(jetBin, trackPt);
1499       
1500       Double_t z = trackPt / jetPt;
1501       Double_t xi = 0;
1502       if(z!=0) xi = TMath::Log(1/z);
1503       
1504       fh2Xi2->Fill(jetBin, xi);
1505       fh2Z2->Fill(jetBin, z);
1506     }
1507
1508
1509 }
1510
1511 //________________________________________________________________________
1512 void AliAnalysisTaskFragmentationFunction::AliFragFuncDiJetHistos::AddToOutput(TList* list)const
1513 {
1514   list->Add(fh1Jet1Pt);
1515   list->Add(fh1Jet2Pt);
1516   list->Add(fh1JetPt);
1517   list->Add(fh2TrackPtJet1);
1518   list->Add(fh2TrackPtJet2);
1519   list->Add(fh2TrackPtJet);
1520   list->Add(fh2Xi1);
1521   list->Add(fh2Xi2);
1522   list->Add(fh2Xi);
1523   list->Add(fh2Z1);
1524   list->Add(fh2Z2);
1525   list->Add(fh2Z);
1526 }
1527
1528 //______________________________________________________________________________________________________
1529 AliAnalysisTaskFragmentationFunction::AliFragFuncQADiJetHistos::AliFragFuncQADiJetHistos(const char* name, Int_t kindSlices,
1530                                             Int_t nInvMass, Float_t invMassMin, Float_t invMassMax, 
1531                                              Int_t nJetPt, Float_t jetPtMin, Float_t jetPtMax,  
1532                                             Int_t nDeltaPhi, Float_t deltaPhiMin, Float_t deltaPhiMax, 
1533                                            Int_t nDeltaEta, Float_t deltaEtaMin, Float_t deltaEtaMax, 
1534                                           Int_t nDeltaPt, Float_t deltaPtMin, Float_t deltaPtMax)
1535   : TObject()
1536   ,fKindSlices(kindSlices)
1537   ,fNBinsJetInvMass(nInvMass)
1538   ,fJetInvMassMin(invMassMin)
1539   ,fJetInvMassMax(invMassMax)
1540   ,fNBinsJetPt(nJetPt)
1541   ,fJetPtMin(jetPtMin)
1542   ,fJetPtMax(jetPtMax)
1543   ,fNBinsDeltaPhi(nDeltaPhi)
1544   ,fDeltaPhiMin(deltaPhiMin)
1545   ,fDeltaPhiMax(deltaPhiMax)
1546   ,fNBinsDeltaEta(nDeltaEta)
1547   ,fDeltaEtaMin(deltaEtaMin)
1548   ,fDeltaEtaMax(deltaEtaMax)
1549   ,fNBinsDeltaPt(nDeltaPt)
1550   ,fDeltaPtMin(deltaPtMin)
1551   ,fDeltaPtMax(deltaPtMax)
1552   ,fh2InvMass(0)
1553   ,fh2DeltaPhi(0)
1554   ,fh2DeltaEta(0)
1555   ,fh2DeltaPt(0)
1556   ,fName(name)
1557 {
1558   // default constructor
1559
1560 }
1561
1562 //______________________________________________________________________________________________________
1563 AliAnalysisTaskFragmentationFunction::AliFragFuncQADiJetHistos::AliFragFuncQADiJetHistos(const AliFragFuncQADiJetHistos& copy)
1564   : TObject()
1565   ,fKindSlices(copy.fKindSlices)
1566   ,fNBinsJetInvMass(copy.fNBinsJetInvMass)
1567   ,fJetInvMassMin(copy.fJetInvMassMin)
1568   ,fJetInvMassMax(copy.fJetInvMassMax)
1569   ,fNBinsJetPt(copy.fNBinsJetPt)
1570   ,fJetPtMin(copy.fJetPtMin)
1571   ,fJetPtMax(copy.fJetPtMax)
1572   ,fNBinsDeltaPhi(copy.fNBinsDeltaPhi)
1573   ,fDeltaPhiMin(copy.fDeltaPhiMin)
1574   ,fDeltaPhiMax(copy.fDeltaPhiMax)
1575   ,fNBinsDeltaEta(copy.fNBinsDeltaEta)
1576   ,fDeltaEtaMin(copy.fDeltaEtaMin)
1577   ,fDeltaEtaMax(copy.fDeltaEtaMax)
1578   ,fNBinsDeltaPt(copy.fNBinsDeltaPt)
1579   ,fDeltaPtMin(copy.fDeltaPtMin)
1580   ,fDeltaPtMax(copy.fDeltaPtMax)
1581   ,fh2InvMass(copy.fh2InvMass)
1582   ,fh2DeltaPhi(copy.fh2DeltaPhi)
1583   ,fh2DeltaEta(copy.fh2DeltaEta)
1584   ,fh2DeltaPt(copy.fh2DeltaPt)
1585   ,fName(copy.fName)
1586 {
1587   // default constructor
1588
1589 }
1590
1591 //_______________________________________________________________________________________________________________________________________________________________
1592 AliAnalysisTaskFragmentationFunction::AliFragFuncQADiJetHistos& AliAnalysisTaskFragmentationFunction::AliFragFuncQADiJetHistos::operator=(const AliAnalysisTaskFragmentationFunction::AliFragFuncQADiJetHistos& o)
1593 {
1594   // assignment
1595   
1596   if(this!=&o){
1597     TObject::operator=(o);
1598     fKindSlices       = o.fKindSlices;
1599     fNBinsJetInvMass  = o.fNBinsJetInvMass;
1600     fJetInvMassMin    = o.fJetInvMassMin;
1601     fJetInvMassMax    = o.fJetInvMassMax;
1602     fNBinsJetPt       = o.fNBinsJetPt;
1603     fJetPtMin         = o.fJetPtMin;
1604     fJetPtMax         = o.fJetPtMax;
1605     fNBinsDeltaPhi    = o.fNBinsDeltaPhi;
1606     fDeltaPhiMin      = o.fDeltaPhiMin;
1607     fDeltaPhiMax      = o.fDeltaPhiMax;
1608     fNBinsDeltaEta    = o.fNBinsDeltaEta;
1609     fDeltaEtaMin      = o.fDeltaEtaMin;
1610     fDeltaEtaMax      = o.fDeltaEtaMax;
1611     fNBinsDeltaPt     = o.fNBinsDeltaPt;
1612     fDeltaPtMin       = o.fDeltaPtMin;
1613     fDeltaPtMax       = o.fDeltaPtMax;
1614     fh2InvMass        = o.fh2InvMass;
1615     fh2DeltaPhi       = o.fh2DeltaPhi;
1616     fh2DeltaEta       = o.fh2DeltaEta;
1617     fh2DeltaPt        = o.fh2DeltaPt;
1618     fName             = o.fName;
1619   }
1620     
1621   return *this;
1622 }
1623
1624 //_________________________________________________________
1625 AliAnalysisTaskFragmentationFunction::AliFragFuncQADiJetHistos::~AliFragFuncQADiJetHistos()
1626 {
1627   // destructor 
1628
1629   if(fh2InvMass)  delete fh2InvMass;
1630   if(fh2DeltaPhi) delete fh2DeltaPhi;
1631   if(fh2DeltaEta) delete fh2DeltaEta;
1632   if(fh2DeltaPt)  delete fh2DeltaPt;
1633 }
1634
1635 //________________________________________________________________________
1636 void AliAnalysisTaskFragmentationFunction::AliFragFuncQADiJetHistos::DefineQADiJetHistos()
1637 {
1638
1639   Int_t nBins = 0;
1640   Double_t min = 0.;
1641   Double_t max = 0.;
1642   const char *xaxis = "";
1643   if(fKindSlices == 1)
1644     {
1645       nBins = fNBinsJetInvMass;
1646       min   = fJetInvMassMin;
1647       max   = fJetInvMassMax;
1648       xaxis = "M_{JJ} [GeV]";
1649     }
1650   if(fKindSlices == 2 || fKindSlices == 3)
1651     {
1652       nBins = fNBinsJetPt;
1653       min   = fJetPtMin;
1654       max   = fJetPtMax;
1655       if(fKindSlices == 2) xaxis = "E_{Tmean} [GeV]";
1656       if(fKindSlices == 3) xaxis ="leading jet p_{T} [GeV/c]";
1657     }
1658   
1659   
1660   fh2InvMass  = new TH2F(Form("fh2DJInvMassPositionCut%s",  fName.Data()), "",nBins, min, max, fNBinsJetInvMass, fJetInvMassMin, fJetInvMassMax);
1661   fh2DeltaPhi = new TH2F(Form("fh2DJDeltaPhiPositionCut%s", fName.Data()), "",nBins, min, max, fNBinsDeltaPhi, fDeltaPhiMin, fDeltaPhiMax);
1662   fh2DeltaEta = new TH2F(Form("fh2DJDeltaEtaPositionCut%s", fName.Data()), "",nBins, min, max, fNBinsDeltaEta, fDeltaEtaMin, fDeltaEtaMax);
1663   fh2DeltaPt  = new TH2F(Form("fh2DJDeltaPtPositionCut%s",  fName.Data()), "",nBins, min, max, fNBinsDeltaPt, fDeltaPtMin, fDeltaPtMax);
1664   
1665   AliAnalysisTaskFragmentationFunction::SetProperties(fh2InvMass, xaxis, "Invariant Mass", "Entries");
1666   AliAnalysisTaskFragmentationFunction::SetProperties(fh2DeltaPhi, xaxis, "#Delta #phi", "Entries");
1667   AliAnalysisTaskFragmentationFunction::SetProperties(fh2DeltaEta, xaxis, "#Delta #eta", "Entries");
1668   AliAnalysisTaskFragmentationFunction::SetProperties(fh2DeltaPt, xaxis, "#Delta p_{T}", "Entries");
1669
1670 }
1671
1672 //________________________________________________________________________
1673 void AliAnalysisTaskFragmentationFunction::AliFragFuncQADiJetHistos::FillDiJetQA(Double_t invMass, Double_t deltaPhi, Double_t deltaEta,Double_t deltaPt, Double_t jetBin)
1674 {
1675   fh2InvMass->Fill(jetBin, invMass);
1676   fh2DeltaPhi->Fill(jetBin, deltaPhi);
1677   fh2DeltaEta->Fill(jetBin, deltaEta);
1678   fh2DeltaPt->Fill(jetBin, deltaPt);
1679 }
1680
1681 //________________________________________________________________________
1682 void AliAnalysisTaskFragmentationFunction::AliFragFuncQADiJetHistos::AddToOutput(TList* list)const
1683 {
1684   list->Add(fh2InvMass);
1685   list->Add(fh2DeltaPhi);
1686   list->Add(fh2DeltaEta);
1687   list->Add(fh2DeltaPt);
1688 }
1689
1690 //_________________________________________________________________________________
1691 void AliAnalysisTaskFragmentationFunction::AliFragFuncIntraJetHistos::AddToOutput(TList* list) const
1692 {
1693   // add histos to list
1694
1695   list->Add(fh2CosTheta);
1696   list->Add(fh2Theta);
1697   list->Add(fh2Jt);
1698
1699   list->Add(fhnIntraJet);
1700
1701 }
1702
1703 //__________________________________________________________________
1704 void AliAnalysisTaskFragmentationFunction::UserCreateOutputObjects()
1705 {
1706   // create output objects
1707
1708   if(fDebug > 1) Printf("AliAnalysisTaskFragmentationFunction::UserCreateOutputObjects()");
1709  
1710   // create list of tracks and jets 
1711   
1712   fTracksRec = new TList();
1713   fTracksRec->SetOwner(kFALSE);  
1714
1715   fTracksRecCuts = new TList();
1716   fTracksRecCuts->SetOwner(kFALSE);  
1717
1718   fTracksGen = new TList();
1719   fTracksGen->SetOwner(kFALSE);
1720
1721   fJetsRec = new TList();
1722   fJetsRec->SetOwner(kFALSE);
1723
1724   fJetsRecCuts = new TList();
1725   fJetsRecCuts->SetOwner(kFALSE);
1726
1727   fJetsGen = new TList();
1728   fJetsGen->SetOwner(kFALSE);
1729
1730   // fJetsKine = new TList();
1731   // fJetsKine->SetOwner(kTRUE); // delete AOD jets using mom from Kine Tree via TList::Clear()
1732
1733
1734   //
1735   // Create histograms / output container
1736   //
1737
1738   OpenFile(1);
1739   fCommonHistList = new TList();
1740   
1741   Bool_t oldStatus = TH1::AddDirectoryStatus();
1742   TH1::AddDirectory(kFALSE);
1743   
1744   
1745   // Histograms 
1746   fh1EvtSelection            = new TH1F("fh1EvtSelection", "Event Selection", 6, -0.5, 5.5);
1747   fh1VertexNContributors     = new TH1F("fh1VertexNContributors", "Vertex N contributors", 11,-.5, 10.5);
1748   fh1VertexZ                 = new TH1F("fh1VertexZ", "Vertex z distribution", 30, -15., 15.);
1749   fh1EvtMult                 = new TH1F("fh1EvtMult","Event multiplicity, track pT cut > 150 MeV/c, |#eta| < 0.9",100,0.,100.);
1750   fh1nRecJetsCuts            = new TH1F("fh1nRecJetsCuts","reconstructed jets per event",10,-0.5,9.5);
1751   fh1nGenJets                = new TH1F("fh1nGenJets","generated jets per event",10,-0.5,9.5);
1752
1753
1754   fQATrackHistosRec          = new AliFragFuncQATrackHistos("Rec", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax, 
1755                                                             fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1756                                                             fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax, 
1757                                                             fQATrackHighPtThreshold);
1758   fQATrackHistosRecCuts      = new AliFragFuncQATrackHistos("RecCuts", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax, 
1759                                                             fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1760                                                             fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax, 
1761                                                             fQATrackHighPtThreshold);
1762   fQATrackHistosGen          = new AliFragFuncQATrackHistos("Gen", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax, 
1763                                                             fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1764                                                             fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax, 
1765                                                             fQATrackHighPtThreshold);
1766   
1767
1768   fQAJetHistosRec            = new AliFragFuncQAJetHistos("Rec", fQAJetNBinsPt, fQAJetPtMin, fQAJetPtMax, 
1769                                                           fQAJetNBinsEta, fQAJetEtaMin, fQAJetEtaMax,
1770                                                           fQAJetNBinsPhi, fQAJetPhiMin, fQAJetPhiMax);
1771   fQAJetHistosRecCuts        = new AliFragFuncQAJetHistos("RecCuts", fQAJetNBinsPt, fQAJetPtMin, fQAJetPtMax, 
1772                                                           fQAJetNBinsEta, fQAJetEtaMin, fQAJetEtaMax,
1773                                                           fQAJetNBinsPhi, fQAJetPhiMin, fQAJetPhiMax);
1774   fQAJetHistosRecCutsLeading = new AliFragFuncQAJetHistos("RecCutsLeading", fQAJetNBinsPt, fQAJetPtMin, fQAJetPtMax, 
1775                                                           fQAJetNBinsEta, fQAJetEtaMin, fQAJetEtaMax,
1776                                                           fQAJetNBinsPhi, fQAJetPhiMin, fQAJetPhiMax);
1777   fQAJetHistosGen            = new AliFragFuncQAJetHistos("Gen", fQAJetNBinsPt, fQAJetPtMin, fQAJetPtMax, 
1778                                                           fQAJetNBinsEta, fQAJetEtaMin, fQAJetEtaMax,
1779                                                           fQAJetNBinsPhi, fQAJetPhiMin, fQAJetPhiMax);
1780   fQAJetHistosGenLeading     = new AliFragFuncQAJetHistos("GenLeading", fQAJetNBinsPt, fQAJetPtMin, fQAJetPtMax, 
1781                                                           fQAJetNBinsEta, fQAJetEtaMin, fQAJetEtaMax,
1782                                                           fQAJetNBinsPhi, fQAJetPhiMin, fQAJetPhiMax);
1783   
1784
1785   fFFHistosRecCuts           = new AliFragFuncHistos("RecCuts", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax, 
1786                                                      fFFNBinsPt, fFFPtMin, fFFPtMax, 
1787                                                      fFFNBinsXi, fFFXiMin, fFFXiMax,  
1788                                                      fFFNBinsZ , fFFZMin , fFFZMax);
1789   fFFHistosRecLeading        = new AliFragFuncHistos("RecLeading", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax, 
1790                                                      fFFNBinsPt, fFFPtMin, fFFPtMax, 
1791                                                      fFFNBinsXi, fFFXiMin, fFFXiMax,  
1792                                                      fFFNBinsZ , fFFZMin , fFFZMax);
1793   fFFHistosRecLeadingTrack   = new AliFragFuncHistos("RecLeadingTrack", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax, 
1794                                                      fFFNBinsPt, fFFPtMin, fFFPtMax, 
1795                                                      fFFNBinsXi, fFFXiMin, fFFXiMax,  
1796                                                      fFFNBinsZ , fFFZMin , fFFZMax);
1797   fFFHistosGen               = new AliFragFuncHistos("Gen", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax, 
1798                                                      fFFNBinsPt, fFFPtMin, fFFPtMax, 
1799                                                      fFFNBinsXi, fFFXiMin, fFFXiMax,  
1800                                                      fFFNBinsZ , fFFZMin , fFFZMax);
1801   fFFHistosGenLeading        = new AliFragFuncHistos("GenLeading", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax, 
1802                                                      fFFNBinsPt, fFFPtMin, fFFPtMax, 
1803                                                      fFFNBinsXi, fFFXiMin, fFFXiMax,  
1804                                                      fFFNBinsZ , fFFZMin , fFFZMax);
1805   fFFHistosGenLeadingTrack   = new AliFragFuncHistos("GenLeadingTrack", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax, 
1806                                                      fFFNBinsPt, fFFPtMin, fFFPtMax, 
1807                                                      fFFNBinsXi, fFFXiMin, fFFXiMax,  
1808                                                      fFFNBinsZ , fFFZMin , fFFZMax);
1809
1810
1811   fIJHistosRecCuts           = new AliFragFuncIntraJetHistos("RecCuts", fIJNBinsJetPt, fIJJetPtMin, fIJJetPtMax, 
1812                                                              fIJNBinsPt, fIJPtMin, fIJPtMax, 
1813                                                              fIJNBinsZ, fIJZMin, fIJZMax,  
1814                                                              fIJNBinsCosTheta , fIJCosThetaMin , fIJCosThetaMax, 
1815                                                              fIJNBinsTheta , fIJThetaMin , fIJThetaMax,
1816                                                              fIJNBinsJt , fIJJtMin , fIJJtMax);
1817   fIJHistosRecLeading        = new AliFragFuncIntraJetHistos("RecLeading", fIJNBinsJetPt, fIJJetPtMin, fIJJetPtMax, 
1818                                                              fIJNBinsPt, fIJPtMin, fIJPtMax, 
1819                                                              fIJNBinsZ, fIJZMin, fIJZMax,  
1820                                                              fIJNBinsCosTheta , fIJCosThetaMin , fIJCosThetaMax, 
1821                                                              fIJNBinsTheta , fIJThetaMin , fIJThetaMax, 
1822                                                              fIJNBinsJt , fIJJtMin , fIJJtMax);
1823   fIJHistosRecLeadingTrack   = new AliFragFuncIntraJetHistos("RecLeadingTrack", fIJNBinsJetPt, fIJJetPtMin, fIJJetPtMax, 
1824                                                              fIJNBinsPt, fIJPtMin, fIJPtMax, 
1825                                                              fIJNBinsZ, fIJZMin, fIJZMax,  
1826                                                              fIJNBinsCosTheta , fIJCosThetaMin , fIJCosThetaMax, 
1827                                                              fIJNBinsTheta , fIJThetaMin , fIJThetaMax, 
1828                                                              fIJNBinsJt , fIJJtMin , fIJJtMax);
1829   fIJHistosGen               = new AliFragFuncIntraJetHistos("Gen", fIJNBinsJetPt, fIJJetPtMin, fIJJetPtMax, 
1830                                                              fIJNBinsPt, fIJPtMin, fIJPtMax, 
1831                                                              fIJNBinsZ, fIJZMin, fIJZMax,  
1832                                                              fIJNBinsCosTheta , fIJCosThetaMin , fIJCosThetaMax, 
1833                                                              fIJNBinsTheta , fIJThetaMin , fIJThetaMax, 
1834                                                              fIJNBinsJt , fIJJtMin , fIJJtMax);
1835   fIJHistosGenLeading        = new AliFragFuncIntraJetHistos("GenLeading", fIJNBinsJetPt, fIJJetPtMin, fIJJetPtMax, 
1836                                                              fIJNBinsPt, fIJPtMin, fIJPtMax, 
1837                                                              fIJNBinsZ, fIJZMin, fIJZMax,  
1838                                                              fIJNBinsCosTheta , fIJCosThetaMin , fIJCosThetaMax, 
1839                                                              fIJNBinsTheta , fIJThetaMin , fIJThetaMax, 
1840                                                              fIJNBinsJt , fIJJtMin , fIJJtMax);
1841   fIJHistosGenLeadingTrack   = new AliFragFuncIntraJetHistos("GenLeadingTrack", fIJNBinsJetPt, fIJJetPtMin, fIJJetPtMax, 
1842                                                              fIJNBinsPt, fIJPtMin, fIJPtMax, 
1843                                                              fIJNBinsZ, fIJZMin, fIJZMax,  
1844                                                              fIJNBinsCosTheta , fIJCosThetaMin , fIJCosThetaMax, 
1845                                                              fIJNBinsTheta , fIJThetaMin , fIJThetaMax, 
1846                                                              fIJNBinsJt , fIJJtMin , fIJJtMax);
1847
1848
1849   fFFDiJetHistosRecCuts         = new AliFragFuncDiJetHistos("RecCuts", fDiJetKindBins, 
1850                                                              fDiJetNBinsJetInvMass, fDiJetJetInvMassMin, fDiJetJetInvMassMax,
1851                                                              fDiJetNBinsJetPt, fDiJetJetPtMin, fDiJetJetPtMax,
1852                                                              fDiJetNBinsPt, fDiJetPtMin, fDiJetPtMax, 
1853                                                              fDiJetNBinsXi, fDiJetXiMin, fDiJetXiMax, 
1854                                                              fDiJetNBinsZ, fDiJetZMin, fDiJetZMax);
1855   fFFDiJetHistosRecLeading      = new AliFragFuncDiJetHistos("RecLeading", fDiJetKindBins, 
1856                                                              fDiJetNBinsJetInvMass, fDiJetJetInvMassMin, fDiJetJetInvMassMax,
1857                                                              fDiJetNBinsJetPt, fDiJetJetPtMin, fDiJetJetPtMax, 
1858                                                              fDiJetNBinsPt, fDiJetPtMin, fDiJetPtMax, 
1859                                                              fDiJetNBinsXi, fDiJetXiMin, fDiJetXiMax, 
1860                                                              fDiJetNBinsZ, fDiJetZMin, fDiJetZMax); 
1861   fFFDiJetHistosRecLeadingTrack = new AliFragFuncDiJetHistos("RecLeadingTrack", fDiJetKindBins, 
1862                                                              fDiJetNBinsJetInvMass, fDiJetJetInvMassMin, fDiJetJetInvMassMax,
1863                                                              fDiJetNBinsJetPt, fDiJetJetPtMin, fDiJetJetPtMax, 
1864                                                              fDiJetNBinsPt, fDiJetPtMin, fDiJetPtMax, 
1865                                                              fDiJetNBinsXi, fDiJetXiMin, fDiJetXiMax, 
1866                                                              fDiJetNBinsZ, fDiJetZMin, fDiJetZMax);
1867
1868   fFFDiJetHistosGen             = new AliFragFuncDiJetHistos("Gen", fDiJetKindBins, 
1869                                                              fDiJetNBinsJetInvMass, fDiJetJetInvMassMin, fDiJetJetInvMassMax,
1870                                                              fDiJetNBinsJetPt, fDiJetJetPtMin, fDiJetJetPtMax,
1871                                                              fDiJetNBinsPt, fDiJetPtMin, fDiJetPtMax, 
1872                                                              fDiJetNBinsXi, fDiJetXiMin, fDiJetXiMax,
1873                                                              fDiJetNBinsZ, fDiJetZMin, fDiJetZMax);
1874   fFFDiJetHistosGenLeading      = new AliFragFuncDiJetHistos("GenLeading", fDiJetKindBins, 
1875                                                              fDiJetNBinsJetInvMass, fDiJetJetInvMassMin, fDiJetJetInvMassMax,
1876                                                              fDiJetNBinsJetPt, fDiJetJetPtMin, fDiJetJetPtMax,
1877                                                              fDiJetNBinsPt, fDiJetPtMin, fDiJetPtMax,
1878                                                              fDiJetNBinsXi, fDiJetXiMin, fDiJetXiMax,
1879                                                              fDiJetNBinsZ, fDiJetZMin, fDiJetZMax);
1880   fFFDiJetHistosGenLeadingTrack = new AliFragFuncDiJetHistos("GenLeadingTrack", fDiJetKindBins, 
1881                                                              fDiJetNBinsJetInvMass, fDiJetJetInvMassMin, fDiJetJetInvMassMax,
1882                                                              fDiJetNBinsJetPt, fDiJetJetPtMin, fDiJetJetPtMax, 
1883                                                              fDiJetNBinsPt, fDiJetPtMin, fDiJetPtMax, 
1884                                                              fDiJetNBinsXi, fDiJetXiMin, fDiJetXiMax, 
1885                                                              fDiJetNBinsZ, fDiJetZMin, fDiJetZMax);
1886
1887   fQADiJetHistosRecCuts = new AliFragFuncQADiJetHistos("RecCuts", fDiJetKindBins, 
1888                                                        fQADiJetNBinsInvMass, fQADiJetInvMassMin, fQADiJetInvMassMax,
1889                                                        fQADiJetNBinsJetPt, fQADiJetJetPtMin, fQADiJetJetPtMax,
1890                                                        fQADiJetNBinsDeltaPhi, fQADiJetDeltaPhiMin, fQADiJetDeltaPhiMax , 
1891                                                        fQADiJetNBinsDeltaEta, fQADiJetDeltaEtaMin, fQADiJetDeltaEtaMax , 
1892                                                        fQADiJetNBinsDeltaPt, fQADiJetDeltaPtMin, fQADiJetDeltaPtMax);
1893   fQADiJetHistosGen     = new AliFragFuncQADiJetHistos("Gen", fDiJetKindBins, 
1894                                                        fQADiJetNBinsInvMass, fQADiJetInvMassMin,  fQADiJetInvMassMax, 
1895                                                        fDiJetNBinsJetPt, fDiJetJetPtMin, fDiJetJetPtMax,
1896                                                        fQADiJetNBinsDeltaPhi, fQADiJetDeltaPhiMin, fQADiJetDeltaPhiMax,
1897                                                        fQADiJetNBinsDeltaEta, fQADiJetDeltaEtaMin, fQADiJetDeltaEtaMax,
1898                                                        fQADiJetNBinsDeltaPt, fQADiJetDeltaPtMin, fQADiJetDeltaPtMax);
1899
1900
1901   fQATrackHistosRec->DefineHistos();
1902   fQATrackHistosRecCuts->DefineHistos();
1903   fQATrackHistosGen->DefineHistos();
1904
1905   fQAJetHistosRec->DefineHistos();
1906   fQAJetHistosRecCuts->DefineHistos();
1907   fQAJetHistosRecCutsLeading->DefineHistos();
1908   fQAJetHistosGen->DefineHistos();
1909   fQAJetHistosGenLeading->DefineHistos();
1910
1911   fFFHistosRecCuts->DefineHistos();
1912   fFFHistosRecLeading->DefineHistos();
1913   fFFHistosRecLeadingTrack->DefineHistos();
1914   fFFHistosGen->DefineHistos();
1915   fFFHistosGenLeading->DefineHistos();
1916   fFFHistosGenLeadingTrack->DefineHistos();
1917
1918   fIJHistosRecCuts->DefineHistos();
1919   fIJHistosRecLeading->DefineHistos();
1920   fIJHistosRecLeadingTrack->DefineHistos();
1921   fIJHistosGen->DefineHistos();
1922   fIJHistosGenLeading->DefineHistos();
1923   fIJHistosGenLeadingTrack->DefineHistos();
1924   
1925   fFFDiJetHistosRecCuts->DefineDiJetHistos();
1926   fFFDiJetHistosRecLeading->DefineDiJetHistos();
1927   fFFDiJetHistosRecLeadingTrack->DefineDiJetHistos();
1928   fFFDiJetHistosGen->DefineDiJetHistos();
1929   fFFDiJetHistosGenLeading->DefineDiJetHistos();
1930   fFFDiJetHistosGenLeadingTrack->DefineDiJetHistos();
1931   fQADiJetHistosRecCuts->DefineQADiJetHistos();
1932   fQADiJetHistosGen->DefineQADiJetHistos();
1933
1934   const Int_t saveLevel = 5;
1935   if(saveLevel>0){
1936     fCommonHistList->Add(fh1EvtSelection);
1937     fFFHistosRecCuts->AddToOutput(fCommonHistList);
1938     fFFHistosRecLeading->AddToOutput(fCommonHistList);
1939     fFFHistosRecLeadingTrack->AddToOutput(fCommonHistList);
1940     fFFHistosGen->AddToOutput(fCommonHistList);
1941     fFFHistosGenLeading->AddToOutput(fCommonHistList);
1942     fFFHistosGenLeadingTrack->AddToOutput(fCommonHistList);
1943   }
1944   if(saveLevel>1){
1945     fQATrackHistosRec->AddToOutput(fCommonHistList);
1946     fQATrackHistosRecCuts->AddToOutput(fCommonHistList);
1947     fQATrackHistosGen->AddToOutput(fCommonHistList);
1948     
1949     fQAJetHistosRec->AddToOutput(fCommonHistList);
1950     fQAJetHistosRecCuts->AddToOutput(fCommonHistList);
1951     fQAJetHistosRecCutsLeading->AddToOutput(fCommonHistList);
1952     fQAJetHistosGen->AddToOutput(fCommonHistList);
1953     fQAJetHistosGenLeading->AddToOutput(fCommonHistList);
1954     
1955     fCommonHistList->Add(fh1EvtMult);
1956     fCommonHistList->Add(fh1nRecJetsCuts);
1957     fCommonHistList->Add(fh1nGenJets);
1958   }
1959   if(saveLevel>2){
1960     fCommonHistList->Add(fh1VertexNContributors);
1961     fCommonHistList->Add(fh1VertexZ);    
1962   }
1963   if(saveLevel>3){
1964     fIJHistosRecCuts->AddToOutput(fCommonHistList);
1965     fIJHistosRecLeading->AddToOutput(fCommonHistList);
1966     fIJHistosRecLeadingTrack->AddToOutput(fCommonHistList);
1967     fIJHistosGen->AddToOutput(fCommonHistList);
1968     fIJHistosGenLeading->AddToOutput(fCommonHistList);
1969     fIJHistosGenLeadingTrack->AddToOutput(fCommonHistList);
1970   }
1971   if(saveLevel>4){
1972     fFFDiJetHistosRecCuts->AddToOutput(fCommonHistList);
1973     fFFDiJetHistosRecLeading->AddToOutput(fCommonHistList);
1974     fFFDiJetHistosRecLeadingTrack->AddToOutput(fCommonHistList);
1975     fFFDiJetHistosGen->AddToOutput(fCommonHistList);
1976     fFFDiJetHistosGenLeading->AddToOutput(fCommonHistList);
1977     fFFDiJetHistosGenLeadingTrack->AddToOutput(fCommonHistList);
1978     fQADiJetHistosRecCuts->AddToOutput(fCommonHistList);
1979     fQADiJetHistosGen->AddToOutput(fCommonHistList);    
1980   }
1981
1982   // =========== Switch on Sumw2 for all histos ===========
1983   for (Int_t i=0; i<fCommonHistList->GetEntries(); ++i){
1984     TH1 *h1 = dynamic_cast<TH1*>(fCommonHistList->At(i));
1985     if (h1) h1->Sumw2();
1986     else{
1987       THnSparse *hnSparse = dynamic_cast<THnSparse*>(fCommonHistList->At(i));
1988       if(hnSparse) hnSparse->Sumw2();
1989     }
1990   }
1991   
1992   TH1::AddDirectory(oldStatus);
1993 }
1994
1995 //_______________________________________________
1996 void AliAnalysisTaskFragmentationFunction::Init()
1997 {
1998   // Initialization
1999   if(fDebug > 1) Printf("AliAnalysisTaskFragmentationFunction::Init()");
2000
2001 }
2002
2003 //_____________________________________________________________
2004 void AliAnalysisTaskFragmentationFunction::UserExec(Option_t *) 
2005 {
2006   // Main loop
2007   // Called for each event
2008   if(fDebug > 1) Printf("AliAnalysisTaskFragmentationFunction::UserExec()");
2009         
2010
2011   if(fDebug > 1) Printf("Analysis event #%5d", (Int_t) fEntry);
2012   // Trigger selection
2013   
2014   AliInputEventHandler* inputHandler = (AliInputEventHandler*)
2015     ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
2016   if(inputHandler->IsEventSelected()&AliVEvent::kMB){
2017     if(fDebug > 1)  Printf(" Trigger Selection: event ACCEPTED ... ");
2018     fh1EvtSelection->Fill(1.);
2019   } else {
2020     fh1EvtSelection->Fill(0.);
2021     if(inputHandler->InheritsFrom("AliESDInputHandler")){ // PhysicsSelection only with ESD input
2022       if (fDebug > 1 ) Printf(" Trigger Selection: event REJECTED ... ");
2023       PostData(1, fCommonHistList);
2024       return;
2025     }
2026   }
2027   
2028   fESD = dynamic_cast<AliESDEvent*>(InputEvent());
2029   if(!fESD){
2030     if(fDebug>3) Printf("%s:%d ESDEvent not found in the input", (char*)__FILE__,__LINE__);
2031   }
2032   
2033   fMCEvent = MCEvent();
2034   if(!fMCEvent){
2035     if(fDebug>3) Printf("%s:%d MCEvent not found in the input", (char*)__FILE__,__LINE__);
2036   }
2037   
2038   // get AOD event from input/ouput
2039   TObject* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
2040   if( handler && handler->InheritsFrom("AliAODInputHandler") ) {
2041     fAOD  =  ((AliAODInputHandler*)handler)->GetEvent();
2042     if (fDebug > 1)  Printf("%s:%d AOD event from input", (char*)__FILE__,__LINE__);
2043   }
2044   else {
2045     handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
2046     if( handler && handler->InheritsFrom("AliAODHandler") ) {
2047       fAOD  = ((AliAODHandler*)handler)->GetAOD();
2048       if (fDebug > 1)  Printf("%s:%d AOD event from output", (char*)__FILE__,__LINE__);
2049     }
2050   }
2051   
2052   if(!fAOD){
2053     Printf("%s:%d AODEvent not found", (char*)__FILE__,__LINE__);
2054     return;
2055   }
2056   
2057
2058   // event selection (vertex) *****************************************
2059   
2060   AliAODVertex* primVtx = fAOD->GetPrimaryVertex();
2061   Int_t nTracksPrim = primVtx->GetNContributors();
2062   fh1VertexNContributors->Fill(nTracksPrim);
2063   
2064   if (fDebug > 1) Printf("%s:%d primary vertex selection: %d", (char*)__FILE__,__LINE__,nTracksPrim);
2065   if(!nTracksPrim){
2066     if (fDebug > 1) Printf("%s:%d primary vertex selection: event REJECTED...",(char*)__FILE__,__LINE__); 
2067     fh1EvtSelection->Fill(2.);
2068     PostData(1, fCommonHistList);
2069     return;
2070   }
2071
2072   fh1VertexZ->Fill(primVtx->GetZ());
2073   
2074   if(TMath::Abs(primVtx->GetZ())>10){
2075     if (fDebug > 1) Printf("%s:%d primary vertex z = %f: event REJECTED...",(char*)__FILE__,__LINE__,primVtx->GetZ()); 
2076     fh1EvtSelection->Fill(3.);
2077     PostData(1, fCommonHistList);
2078     return; 
2079   }
2080
2081   TString primVtxName(primVtx->GetName());
2082
2083   if(primVtxName.CompareTo("TPCVertex",TString::kIgnoreCase) == 1){
2084     if (fDebug > 1) Printf("%s:%d primary vertex selection: TPC vertex, event REJECTED...",(char*)__FILE__,__LINE__);
2085     fh1EvtSelection->Fill(4.);
2086     PostData(1, fCommonHistList);
2087     return;
2088   }
2089   if (fDebug > 1) Printf("%s:%d primary vertex selection: event ACCEPTED ...",(char*)__FILE__,__LINE__); 
2090   fh1EvtSelection->Fill(5.);
2091   
2092   
2093   //___ fetch jets __________________________________________________________________________
2094   
2095   Int_t nJ = GetListOfJets(fJetsRec, kJetsRec);
2096   Int_t nRecJets = 0;
2097   if(nJ>=0) nRecJets = fJetsRec->GetEntries();
2098   if(fDebug>2)Printf("%s:%d Selected Rec jets: %d %d",(char*)__FILE__,__LINE__,nJ,nRecJets);
2099   if(nJ != nRecJets) Printf("%s:%d Mismatch Selected Rec Jets: %d %d",(char*)__FILE__,__LINE__,nJ,nRecJets);
2100
2101   Int_t nJCuts = GetListOfJets(fJetsRecCuts, kJetsRecAcceptance);
2102   Int_t nRecJetsCuts = 0;
2103   if(nJCuts>=0) nRecJetsCuts = fJetsRecCuts->GetEntries();
2104   if(fDebug>2)Printf("%s:%d Selected Rec jets after cuts: %d %d",(char*)__FILE__,__LINE__,nJCuts,nRecJetsCuts);
2105   if(nRecJetsCuts != nJCuts) Printf("%s:%d Mismatch selected Rec jets after cuts: %d %d",(char*)__FILE__,__LINE__,nJCuts,nRecJetsCuts);
2106   fh1nRecJetsCuts->Fill(nRecJetsCuts);
2107
2108   
2109   if(fJetTypeGen==kJetsKine || fJetTypeGen == kJetsKineAcceptance) fJetsGen->SetOwner(kTRUE); // kine aod jets allocated on heap, delete them with TList::Clear() 
2110   Int_t nJGen  = GetListOfJets(fJetsGen, fJetTypeGen);
2111   Int_t nGenJets = 0;
2112   if(nJGen>=0) nGenJets = fJetsGen->GetEntries();
2113   if(fDebug>2)Printf("%s:%d Selected Gen jets: %d %d",(char*)__FILE__,__LINE__,nJGen,nGenJets);
2114   if(nJGen != nGenJets) Printf("%s:%d Mismatch selected Gen jets: %d %d",(char*)__FILE__,__LINE__,nJGen,nGenJets);
2115   fh1nGenJets->Fill(nGenJets);
2116
2117
2118   //____ fetch particles __________________________________________________________
2119   
2120   Int_t nT = GetListOfTracks(fTracksRec, kTrackAOD);
2121   Int_t nRecPart = 0;
2122   if(nT>=0) nRecPart = fTracksRec->GetEntries();
2123   if(fDebug>2)Printf("%s:%d Selected Rec tracks: %d %d",(char*)__FILE__,__LINE__,nT,nRecPart);
2124   if(nRecPart != nT) Printf("%s:%d Mismatch selected Rec tracks: %d %d",(char*)__FILE__,__LINE__,nT,nRecPart);
2125   
2126
2127   Int_t nTCuts = GetListOfTracks(fTracksRecCuts, kTrackAODCuts);
2128   Int_t nRecPartCuts = 0;
2129   if(nTCuts>=0) nRecPartCuts = fTracksRecCuts->GetEntries();
2130   if(fDebug>2)Printf("%s:%d Selected Rec tracks after cuts: %d %d",(char*)__FILE__,__LINE__,nTCuts,nRecPartCuts);
2131   if(nRecPartCuts != nTCuts) Printf("%s:%d Mismatch selected Rec tracks after cuts: %d %d",(char*)__FILE__,__LINE__,nTCuts,nRecPartCuts);
2132   fh1EvtMult->Fill(nRecPartCuts);
2133
2134
2135   Int_t nTGen = GetListOfTracks(fTracksGen,fTrackTypeGen);
2136   Int_t nGenPart = 0;
2137   if(nTGen>=0) nGenPart = fTracksGen->GetEntries();
2138   if(fDebug>2)Printf("%s:%d Selected Gen tracks: %d %d",(char*)__FILE__,__LINE__,nTGen,nGenPart);
2139   if(nGenPart != nTGen) Printf("%s:%d Mismatch selected Gen tracks: %d %d",(char*)__FILE__,__LINE__,nTGen,nGenPart);
2140   
2141   
2142   //____ analysis, fill histos ___________________________________________________
2143   
2144   // loop over tracks
2145
2146   for(Int_t it=0; it<nRecPart; ++it){
2147     AliVParticle *part = dynamic_cast<AliVParticle*>(fTracksRec->At(it));
2148     fQATrackHistosRec->FillTrackQA( part->Eta(), TVector2::Phi_0_2pi(part->Phi()), part->Pt());
2149   }
2150   for(Int_t it=0; it<nRecPartCuts; ++it){
2151     AliVParticle *part = dynamic_cast<AliVParticle*>(fTracksRecCuts->At(it));
2152     fQATrackHistosRecCuts->FillTrackQA( part->Eta(), TVector2::Phi_0_2pi(part->Phi()), part->Pt());
2153   }
2154   for(Int_t it=0; it<nGenPart; ++it){
2155     AliVParticle *part = dynamic_cast<AliVParticle*>(fTracksGen->At(it));
2156     fQATrackHistosGen->FillTrackQA( part->Eta(), TVector2::Phi_0_2pi(part->Phi()), part->Pt());
2157   }
2158   
2159   // loop over jets
2160
2161   for(Int_t ij=0; ij<nRecJets; ++ij){
2162
2163     AliAODJet* jet = dynamic_cast<AliAODJet*>(fJetsRec->At(ij));
2164     fQAJetHistosRec->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt());
2165   }
2166   
2167
2168   for(Int_t ij=0; ij<nRecJetsCuts; ++ij){
2169
2170     AliAODJet* jet = dynamic_cast<AliAODJet*>(fJetsRecCuts->At(ij));
2171     fQAJetHistosRecCuts->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt());
2172
2173     if(ij==0){ // leading jet
2174       
2175       fQAJetHistosRecCutsLeading->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt() );
2176       
2177       TList* jettracklist = new TList();
2178       Double_t sumPt = 0.;
2179       Float_t leadTrackPx = 0.;
2180       Float_t leadTrackPy = 0.;
2181       Float_t leadTrackPz = 0.;
2182       Float_t leadTrackP  = 0.;
2183       Float_t leadTrackPt = 0.;
2184       TLorentzVector* leadTrackV = new TLorentzVector();
2185       
2186       if(GetFFRadius()<=0){
2187         GetJetTracksTrackrefs(jettracklist, jet);
2188        } else {
2189         GetJetTracksPointing(fTracksRecCuts, jettracklist, jet, GetFFRadius(), sumPt);
2190       }
2191       
2192       for(Int_t it=0; it<jettracklist->GetSize(); ++it){
2193         Float_t trackPx = (dynamic_cast<AliVParticle*> (jettracklist->At(it)))->Px();
2194         Float_t trackPy = (dynamic_cast<AliVParticle*> (jettracklist->At(it)))->Py();
2195         Float_t trackPz = (dynamic_cast<AliVParticle*> (jettracklist->At(it)))->Pz();
2196         Float_t trackP = (dynamic_cast<AliVParticle*> (jettracklist->At(it)))->P();
2197         Float_t trackPt = (dynamic_cast<AliVParticle*> (jettracklist->At(it)))->Pt();
2198         Float_t jetPx = jet->Px();
2199         Float_t jetPy = jet->Py();
2200         Float_t jetPz = jet->Pz();
2201         Float_t jetP  = jet->P();
2202         Float_t jetPt = jet->Pt();
2203         TLorentzVector* trackV = new TLorentzVector();
2204         TLorentzVector *jetV = new TLorentzVector();
2205         trackV->SetPxPyPzE(trackPx,trackPy,trackPz,trackP);
2206         jetV->SetPxPyPzE(jetPx,jetPy,jetPz,jetP);
2207
2208         Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
2209         
2210         fFFHistosRecCuts->FillFF( trackPt, jetPt, incrementJetPt);
2211         fIJHistosRecCuts->FillIntraJet( trackV, jetV );
2212         
2213         if(it==0){ 
2214           leadTrackPx = trackPx;
2215           leadTrackPy = trackPy;
2216           leadTrackPz = trackPz;
2217           leadTrackP  = trackP;
2218           leadTrackPt = trackPt;
2219           fFFHistosRecLeadingTrack->FillFF( leadTrackPt, jetPt, kTRUE);
2220
2221           leadTrackV->SetPxPyPzE(leadTrackPx,leadTrackPy,leadTrackPz,leadTrackP);
2222           fIJHistosRecLeadingTrack->FillIntraJet( leadTrackV, jetV );
2223         }
2224         fFFHistosRecLeading->FillFF( trackPt, leadTrackPt , incrementJetPt);
2225         fIJHistosRecLeading->FillIntraJet( trackV, leadTrackV );
2226
2227         delete trackV;
2228         delete jetV;
2229       }
2230       
2231       delete leadTrackV;
2232       delete jettracklist;
2233     }
2234   }
2235         
2236
2237   for(Int_t ij=0; ij<nGenJets; ++ij){
2238
2239     AliAODJet* jet = dynamic_cast<AliAODJet*>(fJetsGen->At(ij));
2240     fQAJetHistosGen->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt());
2241     
2242     if(ij==0){ // leading jet
2243     
2244       fQAJetHistosGenLeading->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt());
2245       
2246       TList* jettracklist = new TList();
2247       Double_t sumPt = 0.;
2248       Float_t leadTrackPx = 0.;
2249       Float_t leadTrackPy = 0.;
2250       Float_t leadTrackPz = 0.;
2251       Float_t leadTrackP  = 0.;
2252       Float_t leadTrackPt = 0.;
2253       TLorentzVector* leadTrackV = new TLorentzVector();
2254
2255       if(GetFFRadius()<=0){
2256         GetJetTracksTrackrefs(jettracklist, jet);
2257       } else {
2258         GetJetTracksPointing(fTracksGen, jettracklist, jet, GetFFRadius(), sumPt);
2259       }
2260       
2261       for(Int_t it=0; it<jettracklist->GetSize(); ++it){
2262         Float_t trackPx = (dynamic_cast<AliVParticle*>(jettracklist->At(it)))->Px();
2263         Float_t trackPy = (dynamic_cast<AliVParticle*>(jettracklist->At(it)))->Py();
2264         Float_t trackPz = (dynamic_cast<AliVParticle*>(jettracklist->At(it)))->Pz();
2265         Float_t trackP  = (dynamic_cast<AliVParticle*>(jettracklist->At(it)))->P();
2266         Float_t trackPt = (dynamic_cast<AliVParticle*>(jettracklist->At(it)))->Pt();
2267         Float_t jetPx = jet->Px();
2268         Float_t jetPy = jet->Py();
2269         Float_t jetPz = jet->Pz();
2270         Float_t jetP  = jet->P();
2271         Float_t jetPt = jet->Pt();
2272         TLorentzVector* trackV = new TLorentzVector();
2273         TLorentzVector  *jetV = new TLorentzVector();
2274         trackV->SetPxPyPzE(trackPx,trackPy,trackPz,trackP);
2275         jetV->SetPxPyPzE(jetPx,jetPy,jetPz,jetP);
2276
2277         Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
2278
2279         fFFHistosGen->FillFF( trackPt, jetPt, incrementJetPt);
2280         fIJHistosGen->FillIntraJet( trackV, jetV );
2281         
2282         if(it==0){ 
2283           leadTrackPx = trackPx;
2284           leadTrackPy = trackPy;
2285           leadTrackPz = trackPz;
2286           leadTrackP  = trackP;
2287           leadTrackPt = trackPt;
2288           fFFHistosGenLeadingTrack->FillFF( leadTrackPt, jetPt, kTRUE);
2289
2290           leadTrackV->SetPxPyPzE(leadTrackPx,leadTrackPy,leadTrackPz,leadTrackP);
2291           fIJHistosGenLeadingTrack->FillIntraJet( leadTrackV, jetV );
2292         }
2293         fFFHistosGenLeading->FillFF( trackPt, leadTrackPt, incrementJetPt);
2294         fIJHistosGenLeading->FillIntraJet( trackV, leadTrackV );
2295
2296         delete trackV;
2297         delete jetV;
2298       }
2299       
2300       delete leadTrackV;
2301       delete jettracklist;
2302     }
2303   }
2304
2305   //_______ DiJet part _____________________________________________________
2306
2307   if (nRecJets > 1)
2308   {
2309
2310     AliAODJet* jet1 = dynamic_cast<AliAODJet*>(fJetsRecCuts->At(0));
2311     AliAODJet* jet2 = dynamic_cast<AliAODJet*>(fJetsRecCuts->At(1));
2312     
2313     // DiJet deltaphi calculation
2314     Double_t phi1 = TVector2::Phi_0_2pi(jet1->Phi());
2315     Double_t phi2 = TVector2::Phi_0_2pi(jet2->Phi());
2316     Double_t deltaPhi = TMath::Abs(phi1-phi2); 
2317     if (deltaPhi > TMath::Pi() && deltaPhi < 2*TMath::Pi()) deltaPhi = 2*TMath::Pi() - deltaPhi;
2318     
2319     // DiJet CDF cut calculation
2320     Double_t Et1     = TMath::Abs(jet1->E()*TMath::Sin(jet1->Theta()));
2321     Double_t Et2     = TMath::Abs(jet2->E()*TMath::Sin(jet2->Theta()));
2322     Double_t sumEt   = Et1 + Et2;
2323     Double_t normEt1PlusEt2   = TMath::Sqrt(Et1*Et1+Et2*Et2+2*Et1*Et2*TMath::Cos(deltaPhi));
2324     Double_t ratio = (Double_t)(normEt1PlusEt2/sumEt);
2325     
2326     // DiJet events selection
2327     Bool_t positionCut       = 0;
2328     Bool_t positionEnergyCut = 0;
2329     Bool_t cdfCut            = 0; 
2330
2331     // Position cut :
2332     if (deltaPhi > fDiJetDeltaPhiCut) positionCut = 1;
2333     // Position-Energy cut :
2334     if ((deltaPhi > fDiJetDeltaPhiCut) && ((jet2->Pt()) >= fDiJetPtFractionCut*(jet1->Pt()))) positionEnergyCut = 1;
2335     // CDF cut :
2336     if (ratio < fDiJetCDFCut) cdfCut = 1;
2337     
2338     Int_t go = 0;
2339     
2340     if (fDiJetCut == 1 && positionCut == 1) go = 1;
2341     if (fDiJetCut == 2 && positionEnergyCut == 1) go = 1;
2342     if (fDiJetCut == 3 && cdfCut == 1) go = 1;
2343
2344     if (go)
2345       {
2346         Double_t deltaEta      = TMath::Abs(jet1->Eta()-jet2->Eta());
2347         Double_t deltaPt       = TMath::Abs(jet1->Pt()-jet2->Pt());
2348         Double_t meanEt        = (Double_t)((Et1+Et2)/2.);
2349         Double_t invariantMass = (Double_t)InvMass(jet1,jet2);
2350         
2351         Double_t  jetBin = GetDiJetBin(invariantMass, jet1->Pt(), meanEt, fDiJetKindBins);
2352
2353         if (jetBin > 0)
2354           {
2355             fQADiJetHistosRecCuts->FillDiJetQA(invariantMass, deltaPhi, deltaEta, deltaPt, jetBin);
2356             
2357             TList* jettracklist1 = new TList();
2358             Double_t sumPt1      = 0.;
2359             Float_t leadTrackPt1 = 0;
2360             
2361             TList* jettracklist2 = new TList();
2362             Double_t sumPt2      = 0.;
2363             Float_t leadTrackPt2 = 0;
2364             
2365             if(GetFFRadius()<=0)
2366               {
2367                 GetJetTracksTrackrefs(jettracklist1, jet1);
2368                 GetJetTracksTrackrefs(jettracklist2, jet2);
2369               }
2370             else
2371               {
2372                 GetJetTracksPointing(fTracksRecCuts, jettracklist1, jet1, GetFFRadius(), sumPt1);
2373                 GetJetTracksPointing(fTracksRecCuts, jettracklist2, jet2, GetFFRadius(), sumPt2);
2374               }
2375             
2376             Int_t nTracks = jettracklist1->GetSize(); 
2377             if (jettracklist1->GetSize() < jettracklist2->GetSize()) nTracks = jettracklist2->GetSize();
2378             
2379             for(Int_t it=0; it<nTracks; ++it)
2380               {
2381                 if (it < jettracklist1->GetSize())
2382                   { 
2383                     Float_t trackPt1 = (dynamic_cast<AliVParticle*> (jettracklist1->At(it)))->Pt();
2384                     Float_t jetPt1   = jet1->Pt();
2385                     
2386                     Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
2387                     
2388                     fFFDiJetHistosRecCuts->FillDiJetFF(1, trackPt1, jetPt1, jetBin, incrementJetPt);
2389                     fFFDiJetHistosRecCuts->FillDiJetFF(0, trackPt1, jetPt1, jetBin, incrementJetPt);
2390                     
2391                     if (it == 0)
2392                       {
2393                         leadTrackPt1 = trackPt1;
2394                         
2395                         fFFDiJetHistosRecLeadingTrack->FillDiJetFF(1, leadTrackPt1, jetPt1, jetBin, kTRUE); 
2396                         fFFDiJetHistosRecLeadingTrack->FillDiJetFF(0, leadTrackPt1, jetPt1, jetBin, kTRUE); 
2397                       }
2398                     
2399                     fFFDiJetHistosRecLeading->FillDiJetFF(1, trackPt1, leadTrackPt1, jetBin, incrementJetPt); 
2400                     fFFDiJetHistosRecLeading->FillDiJetFF(0, trackPt1, leadTrackPt1, jetBin, incrementJetPt); 
2401                   }
2402                 
2403                 if (it < jettracklist2->GetSize())
2404                   { 
2405                     Float_t trackPt2   = (dynamic_cast<AliVParticle*>(jettracklist2->At(it)))->Pt();
2406                     Float_t jetPt2     = jet2->Pt();
2407                     
2408                     Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
2409                     
2410                     fFFDiJetHistosRecCuts->FillDiJetFF(2, trackPt2, jetPt2, jetBin, incrementJetPt);
2411                     fFFDiJetHistosRecCuts->FillDiJetFF(0, trackPt2, jetPt2, jetBin, incrementJetPt);
2412                     
2413                     if (it == 0)
2414                       {
2415                         leadTrackPt2 = trackPt2;
2416                         
2417                         fFFDiJetHistosRecLeadingTrack->FillDiJetFF(2, leadTrackPt2, jetPt2, jetBin, kTRUE);
2418                         fFFDiJetHistosRecLeadingTrack->FillDiJetFF(0, leadTrackPt2, jetPt2, jetBin, kTRUE);
2419                       }
2420                     
2421                     fFFDiJetHistosRecLeading->FillDiJetFF(2, trackPt2, leadTrackPt2, jetBin, incrementJetPt);
2422                     fFFDiJetHistosRecLeading->FillDiJetFF(0, trackPt2, leadTrackPt2, jetBin, incrementJetPt);
2423                   }
2424               } // End loop on tracks
2425
2426             delete jettracklist1;
2427             delete jettracklist2;
2428
2429           } // End if(jetBin > 0)
2430         else { Printf("Jet bins for di-jet studies not set !");}
2431       } // End if(go)
2432   } // End if(nRecJets > 1)
2433
2434   if (nGenJets > 1)
2435   {
2436     AliAODJet* jet1 = dynamic_cast<AliAODJet*>(fJetsGen->At(0));
2437     AliAODJet* jet2 = dynamic_cast<AliAODJet*>(fJetsGen->At(1));
2438
2439     Double_t deltaPhi = 0;
2440     Double_t phi1 = TVector2::Phi_0_2pi(jet1->Phi());
2441     Double_t phi2 = TVector2::Phi_0_2pi(jet2->Phi());
2442     deltaPhi      = TMath::Abs(phi1-phi2); 
2443     if (deltaPhi > TMath::Pi() && deltaPhi < 2*TMath::Pi()) deltaPhi = 2*TMath::Pi() - deltaPhi;
2444
2445     Double_t Et1            = TMath::Abs(jet1->E()*TMath::Sin(jet1->Theta()));
2446     Double_t Et2            = TMath::Abs(jet2->E()*TMath::Sin(jet2->Theta()));
2447     Double_t sumEt          = Et1 + Et2;
2448     Double_t normEt1PlusEt2 = TMath::Sqrt(Et1*Et1+Et2*Et2+2*Et1*Et2*TMath::Cos(deltaPhi));
2449     Double_t ratio          = (Double_t)(normEt1PlusEt2/sumEt);
2450
2451     // DiJet events selection
2452     Bool_t positionCut       = 0;
2453     Bool_t positionEnergyCut = 0;
2454     Bool_t cdfCut            = 0; 
2455
2456     // Position cut :
2457     if (deltaPhi > fDiJetDeltaPhiCut) positionCut = 1;
2458     // Position-Energy cut :
2459     if ((deltaPhi > fDiJetDeltaPhiCut) && ((jet2->Pt()) >= fDiJetPtFractionCut*(jet1->Pt()))) positionEnergyCut = 1;
2460     // CDF cut :
2461     if (ratio < fDiJetCDFCut) cdfCut = 1;    
2462
2463     Int_t go = 0;
2464
2465     if (fDiJetCut == 1 && positionCut == 1) go = 1;
2466     if (fDiJetCut == 2 && positionEnergyCut == 1) go = 1;
2467     if (fDiJetCut == 3 && cdfCut == 1) go = 1;
2468
2469     if (go)
2470     {
2471       Double_t deltaEta      = TMath::Abs(jet1->Eta()-jet2->Eta());
2472       Double_t deltaPt       = TMath::Abs(jet1->Pt()-jet2->Pt());
2473       Double_t meanEt        = (Double_t)((Et1+Et2)/2.);
2474       Double_t invariantMass = (Double_t)InvMass(jet1,jet2);
2475
2476       Double_t jetBin = GetDiJetBin(invariantMass, jet1->Pt(), meanEt, fDiJetKindBins);
2477
2478       if(jetBin > 0)
2479       {
2480         fQADiJetHistosGen->FillDiJetQA(invariantMass, deltaPhi, deltaEta, deltaPt, jetBin);
2481
2482         TList* jettracklist1 = new TList();
2483         Double_t sumPt1 = 0.;
2484         Float_t leadTrackPt1 = 0.;
2485
2486         TList* jettracklist2 = new TList();
2487         Double_t sumPt2 = 0.;
2488         Float_t leadTrackPt2 = 0.;
2489       
2490         if(GetFFRadius()<=0)
2491         {
2492           GetJetTracksTrackrefs(jettracklist1, jet1);
2493           GetJetTracksTrackrefs(jettracklist2, jet2);
2494         }
2495         else
2496         {
2497           GetJetTracksPointing(fTracksGen, jettracklist1, jet1, GetFFRadius(), sumPt1);
2498           GetJetTracksPointing(fTracksGen, jettracklist2, jet2, GetFFRadius(), sumPt2);
2499         }
2500       
2501         Int_t nTracks = jettracklist1->GetSize(); 
2502         if (jettracklist1->GetSize() < jettracklist2->GetSize()) nTracks = jettracklist2->GetSize();
2503
2504         for(Int_t it=0; it<nTracks; ++it)
2505         {
2506           if (it < jettracklist1->GetSize())
2507           { 
2508             Float_t trackPt1 = (dynamic_cast<AliAODTrack*>(jettracklist1->At(it)))->Pt();
2509             Float_t jetPt1 = jet1->Pt();
2510
2511             Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
2512
2513             fFFDiJetHistosGen->FillDiJetFF( 1, trackPt1, jetPt1, jetBin, incrementJetPt);
2514             fFFDiJetHistosGen->FillDiJetFF( 0, trackPt1, jetPt1, jetBin, incrementJetPt);
2515
2516             if(it==0)
2517             { 
2518               leadTrackPt1 = trackPt1;
2519
2520               fFFDiJetHistosGenLeadingTrack->FillDiJetFF( 1, leadTrackPt1, jetPt1, jetBin, kTRUE);
2521               fFFDiJetHistosGenLeadingTrack->FillDiJetFF( 0, leadTrackPt1, jetPt1, jetBin, kTRUE);
2522             }
2523
2524             fFFDiJetHistosGenLeading->FillDiJetFF( 1, trackPt1, leadTrackPt1, jetBin, incrementJetPt);
2525             fFFDiJetHistosGenLeading->FillDiJetFF( 0, trackPt1, leadTrackPt1, jetBin, incrementJetPt);
2526           }
2527       
2528           if (it < jettracklist2->GetSize())
2529           { 
2530             Float_t trackPt2 = (dynamic_cast<AliAODTrack*>(jettracklist2->At(it)))->Pt();
2531             Float_t jetPt2 = jet2->Pt();
2532
2533             Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
2534
2535             fFFDiJetHistosGen->FillDiJetFF( 2, trackPt2, jetPt2, jetBin, incrementJetPt);
2536             fFFDiJetHistosGen->FillDiJetFF( 0, trackPt2, jetPt2, jetBin, incrementJetPt);
2537         
2538             if (it==0)
2539             { 
2540               leadTrackPt2 = trackPt2;
2541
2542               fFFDiJetHistosGenLeadingTrack->FillDiJetFF( 2, leadTrackPt2, jetPt2, jetBin, kTRUE);
2543               fFFDiJetHistosGenLeadingTrack->FillDiJetFF( 0, leadTrackPt2, jetPt2, jetBin, kTRUE);
2544             }
2545
2546             fFFDiJetHistosGenLeading->FillDiJetFF( 2, trackPt2, leadTrackPt2, jetBin, incrementJetPt);
2547             fFFDiJetHistosGenLeading->FillDiJetFF( 0, trackPt2, leadTrackPt2, jetBin, incrementJetPt);
2548           }
2549         } // End loop on tracks
2550
2551         delete jettracklist1;
2552         delete jettracklist2;
2553
2554       } // End if(jetBin > 0)
2555       else { Printf("Jet bins for di-jet studies not set !");}
2556     } // End if (go)
2557   } // End if(nGenJets > 1)
2558
2559   
2560   fTracksRec->Clear();
2561   fTracksRecCuts->Clear();
2562   fTracksGen->Clear();
2563   fJetsRec->Clear();
2564   fJetsRecCuts->Clear();
2565   fJetsGen->Clear();
2566
2567   //Post output data.
2568   PostData(1, fCommonHistList);
2569   
2570 }
2571
2572 //________________________________________________________________________________________
2573 Double_t AliAnalysisTaskFragmentationFunction::InvMass(AliAODJet* jet1, AliAODJet* jet2)
2574 {
2575
2576   Double_t invMass = 0.;
2577   invMass = TMath::Sqrt(pow(jet1->E()+jet2->E(),2) - pow(jet1->Px()+jet2->Px(),2) - 
2578                         pow(jet1->Py()+jet2->Py(),2) - pow(jet1->Pz()+jet2->Pz(),2));
2579
2580   return invMass;
2581
2582 }
2583
2584 //________________________________________________________________________________________
2585 Double_t AliAnalysisTaskFragmentationFunction::GetDiJetBin(Double_t invMass, Double_t leadingJetPt, Double_t EtMean, Int_t kindBins)
2586 {
2587   Double_t jetBinOk = 0.;
2588   Double_t jetBin = 0.;
2589
2590   Float_t stepInvMass = (fDiJetJetInvMassMax - fDiJetJetInvMassMin)/fDiJetNBinsJetInvMass;
2591   Float_t stepPt = (fDiJetJetPtMax - fDiJetJetPtMin)/fDiJetNBinsJetPt;
2592
2593   if (kindBins == 1)
2594     {
2595       for(Int_t i=0; i<fDiJetNBinsJetInvMass; ++i)
2596         {
2597           jetBin = fDiJetJetInvMassMin + i*stepInvMass/2.;
2598           if(((fDiJetJetInvMassMin+i*stepInvMass) <= invMass) &&
2599              (fDiJetJetInvMassMin + (i+1)*stepInvMass) > invMass) jetBinOk = jetBin;
2600         }
2601       jetBinOk = -1.;
2602     }
2603   else if (kindBins == 3)
2604     {
2605       for(Int_t i=0; i<fDiJetNBinsJetPt; ++i)
2606         {
2607           jetBin = fDiJetJetPtMin + i*stepPt/2.;
2608           if(((fDiJetJetPtMin+i*stepPt) <= EtMean) &&
2609              (fDiJetJetPtMin + (i+1)*stepPt) > EtMean) jetBinOk = jetBin;
2610         }
2611       jetBinOk = -1.;
2612     }
2613   else if (kindBins == 2)
2614     {
2615       for(Int_t i=0; i<fDiJetNBinsJetPt; ++i)
2616         {
2617           jetBin = fDiJetJetPtMin + i*stepPt/2.;
2618           if(((fDiJetJetPtMin+i*stepPt) <= leadingJetPt) &&
2619              (fDiJetJetPtMin + (i+1)*stepPt) > leadingJetPt) jetBinOk = jetBin;
2620         }
2621       jetBinOk = -1.;
2622     }
2623   else {Printf("WARNING: kindBins wrongly set ! Please make sure to call SetKindSlices() and set the kind parameter to 1, 2 or 3.\n");}
2624
2625   return jetBinOk;
2626
2627 }
2628
2629
2630 //______________________________________________________________
2631 void AliAnalysisTaskFragmentationFunction::Terminate(Option_t *) 
2632 {
2633   // terminated
2634
2635   if(fDebug > 1) printf("AliAnalysisTaskFragmentationFunction::Terminate() \n");
2636 }  
2637
2638 //_________________________________________________________________________________
2639 Int_t AliAnalysisTaskFragmentationFunction::GetListOfTracks(TList *list, Int_t type)
2640 {
2641   // fill list of tracks selected according to type
2642
2643   if(fDebug > 2) Printf("%s:%d Selecting tracks with %d", (char*)__FILE__,__LINE__,type);
2644   
2645   if(!list){
2646     if(fDebug>1) Printf("%s:%d no input list", (char*)__FILE__,__LINE__);
2647     return -1;
2648   }
2649
2650   if(type==kTrackUndef) return 0;
2651   
2652   Int_t iCount = 0;
2653   if(type==kTrackAODCuts || type==kTrackAOD){
2654
2655     // all rec. tracks, esd filter mask, eta range
2656     if(!fAOD) return -1;
2657     
2658     for(Int_t it=0; it<fAOD->GetNumberOfTracks(); ++it){
2659       AliAODTrack *tr = fAOD->GetTrack(it);
2660       
2661       if(type == kTrackAODCuts){
2662         if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask)))   continue;
2663         if(tr->Eta() < fTrackEtaMin || tr->Eta() > fTrackEtaMax) continue;
2664         if(tr->Phi() < fTrackPhiMin || tr->Phi() > fTrackPhiMax) continue;
2665         if(tr->Pt()  < fTrackPtCut) continue;
2666       }
2667       list->Add(tr);
2668       iCount++;
2669     }
2670   }
2671   else if (type==kTrackKineAll || type==kTrackKineCharged || type==kTrackKineChargedAcceptance){
2672     // kine particles, all or rather charged
2673     if(!fMCEvent) return iCount;
2674     
2675     for(Int_t it=0; it<fMCEvent->GetNumberOfTracks(); ++it){
2676       AliMCParticle* part = (AliMCParticle*) fMCEvent->GetTrack(it);
2677       
2678       if(type == kTrackKineCharged || type == kTrackKineChargedAcceptance){
2679         if(part->Charge()==0) continue;
2680         
2681         if(type == kTrackKineChargedAcceptance && 
2682            (       part->Eta() < fTrackEtaMin
2683                 || part->Eta() > fTrackEtaMax
2684                 || part->Phi() < fTrackPhiMin
2685                 || part->Phi() > fTrackPhiMax 
2686                 || part->Pt()  < fTrackPtCut)) continue;
2687       }
2688       
2689       list->Add(part);
2690       iCount++;
2691     }
2692   }
2693   else if (type==kTrackAODMCCharged || type==kTrackAODMCAll || type==kTrackAODMCChargedAcceptance) {
2694     // MC particles (from AOD), physical primaries, all or rather charged or rather charged within acceptance
2695     if(!fAOD) return -1;
2696     
2697     TClonesArray *tca = dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName()));
2698     if(!tca)return iCount;
2699     
2700     for(int it=0; it<tca->GetEntriesFast(); ++it){
2701       AliAODMCParticle *part = dynamic_cast<AliAODMCParticle*>(tca->At(it));
2702       if(!part->IsPhysicalPrimary())continue;
2703       
2704       if (type==kTrackAODMCCharged || type==kTrackAODMCChargedAcceptance){
2705         if(part->Charge()==0) continue;
2706         if(type==kTrackAODMCChargedAcceptance && 
2707            (     part->Eta() > fTrackEtaMax
2708               || part->Eta() < fTrackEtaMin
2709               || part->Phi() > fTrackPhiMax
2710               || part->Phi() < fTrackPhiMin
2711               || part->Pt()  < fTrackPtCut)) continue;
2712       }
2713       
2714       list->Add(part);
2715       iCount++;
2716     }
2717   }
2718   
2719   list->Sort();
2720   return iCount;
2721   
2722 }
2723 // _______________________________________________________________________________
2724 Int_t AliAnalysisTaskFragmentationFunction::GetListOfJets(TList *list, Int_t type)
2725 {
2726   // fill list of jets selected according to type
2727   
2728   if(!list){
2729     if(fDebug>1) Printf("%s:%d no input list", (char*)__FILE__,__LINE__);
2730     return -1;
2731   }
2732
2733   if(type == kJetsRec || type == kJetsRecAcceptance){ // reconstructed jets
2734
2735     if(fBranchRecJets.Length()==0){
2736       Printf("%s:%d no rec jet branch specified", (char*)__FILE__,__LINE__);
2737       if(fDebug>1)fAOD->Print();
2738       return 0;
2739     }
2740
2741     TClonesArray *aodRecJets = new TClonesArray();
2742     if(fBranchRecJets.Length()) aodRecJets = dynamic_cast<TClonesArray*>(fAOD->FindListObject(fBranchRecJets.Data()));\r
2743     if(!aodRecJets)             aodRecJets = dynamic_cast<TClonesArray*>(fAOD->GetList()->FindObject(fBranchRecJets.Data()));\r
2744
2745     if(!aodRecJets){
2746       if(fBranchRecJets.Length()) Printf("%s:%d no reconstructed jet array with name %s in AOD", (char*)__FILE__,__LINE__,fBranchRecJets.Data());
2747
2748       if(fDebug>1)fAOD->Print();
2749       return 0;
2750     }
2751
2752     Int_t nRecJets = 0;
2753     
2754     for(Int_t ij=0; ij<aodRecJets->GetEntries(); ++ij){
2755
2756       AliAODJet *tmp = dynamic_cast<AliAODJet*>(aodRecJets->At(ij));
2757       if(!tmp) continue;
2758         
2759       if( tmp->Pt() < fJetPtCut ) continue;
2760       if( type == kJetsRecAcceptance &&
2761           (    tmp->Eta() < fJetEtaMin
2762             || tmp->Eta() > fJetEtaMax
2763             || tmp->Phi() < fJetPhiMin
2764             || tmp->Phi() > fJetPhiMax )) continue;
2765       
2766       list->Add(tmp);
2767           
2768       nRecJets++;
2769     }
2770
2771     list->Sort();
2772     return nRecJets;
2773     delete aodRecJets;
2774   }
2775   else if(type == kJetsKine || type == kJetsKineAcceptance){
2776     
2777     // generated jets
2778     Int_t nGenJets = 0;
2779     
2780     if(!fMCEvent){
2781       if(fDebug>1) Printf("%s:%d no mcEvent",(char*)__FILE__,__LINE__);
2782       return 0;
2783     }
2784     
2785     AliGenPythiaEventHeader*  pythiaGenHeader = AliAnalysisHelperJetTasks::GetPythiaEventHeader(fMCEvent);
2786     if(!pythiaGenHeader){
2787       Printf("%s:%d no pythiaGenHeader found", (char*)__FILE__,__LINE__);
2788       return 0;
2789     }
2790     
2791     // fetch the pythia generated jets
2792     for(int ip=0; ip<pythiaGenHeader->NTriggerJets(); ++ip){
2793       
2794       Float_t p[4];
2795       AliAODJet *jet = new AliAODJet();
2796       pythiaGenHeader->TriggerJet(ip, p);
2797       jet->SetPxPyPzE(p[0], p[1], p[2], p[3]);
2798
2799       if( type == kJetsKineAcceptance &&
2800           (    jet->Eta() < fJetEtaMin
2801             || jet->Eta() > fJetEtaMax
2802             || jet->Phi() < fJetPhiMin
2803             || jet->Phi() > fJetPhiMax )) continue;
2804       
2805       list->Add(jet);
2806       nGenJets++;
2807     }
2808     list->Sort();
2809     return nGenJets;
2810   }
2811   else if(type == kJetsGen || type == kJetsGenAcceptance ){
2812
2813     if(fBranchGenJets.Length()==0){
2814       if(fDebug>1) Printf("%s:%d no gen jet branch specified", (char*)__FILE__,__LINE__);
2815       return 0;
2816     }
2817     
2818     TClonesArray *aodGenJets = new TClonesArray();
2819     if(fBranchGenJets.Length()) aodGenJets = dynamic_cast<TClonesArray*>(fAOD->FindListObject(fBranchGenJets.Data()));\r
2820     if(!aodGenJets)             aodGenJets = dynamic_cast<TClonesArray*>(fAOD->GetList()->FindObject(fBranchGenJets.Data()));\r
2821
2822     if(!aodGenJets){
2823       if(fDebug>0){
2824         if(fBranchGenJets.Length())         Printf("%s:%d Generated jet branch %s not found",(char*)__FILE__,__LINE__,fBranchGenJets.Data());
2825       }
2826       if(fDebug>1)fAOD->Print();
2827       return 0;
2828     }
2829
2830     Int_t nGenJets = 0;
2831     
2832     for(Int_t ig=0; ig<aodGenJets->GetEntries(); ++ig){
2833           
2834       AliAODJet *tmp = dynamic_cast<AliAODJet*>(aodGenJets->At(ig));
2835       if(!tmp) continue;
2836           
2837       if( tmp->Pt() < fJetPtCut ) continue;
2838       if( type == kJetsGenAcceptance &&
2839           (    tmp->Eta() < fJetEtaMin
2840             || tmp->Eta() > fJetEtaMax
2841             || tmp->Phi() < fJetPhiMin
2842             || tmp->Phi() > fJetPhiMax )) continue;
2843       
2844       list->Add(tmp);
2845       
2846       nGenJets++;
2847     }
2848     list->Sort();
2849     return nGenJets;
2850     delete aodGenJets;
2851   } 
2852   else{
2853     if(fDebug>0)Printf("%s:%d no such type %d",(char*)__FILE__,__LINE__,type);
2854     return 0;
2855   }
2856 }
2857
2858 // _________________________________________________________________________________________________________
2859 void AliAnalysisTaskFragmentationFunction::SetProperties(THnSparse* h,const Int_t dim, const char** labels)
2860 {
2861   //Set properties of THnSparse 
2862
2863   for(Int_t i=0; i<dim; i++){
2864
2865     h->GetAxis(i)->SetTitle(labels[i]);
2866     h->GetAxis(i)->SetTitleColor(1);
2867   }
2868 }
2869
2870 // __________________________________________________________________________________________
2871 void AliAnalysisTaskFragmentationFunction::SetProperties(TH1* h,const char* x, const char* y)
2872 {
2873   //Set properties of histos (x and y title)
2874
2875   h->SetXTitle(x);
2876   h->SetYTitle(y);
2877   h->GetXaxis()->SetTitleColor(1);
2878   h->GetYaxis()->SetTitleColor(1);
2879 }
2880
2881 // _________________________________________________________________________________________________________
2882 void AliAnalysisTaskFragmentationFunction::SetProperties(TH2* h,const char* x, const char* y, const char* z)
2883 {
2884   //Set properties of histos (x,y and z title)
2885
2886   h->SetXTitle(x);
2887   h->SetYTitle(y);
2888   h->SetZTitle(z);
2889   h->GetXaxis()->SetTitleColor(1);
2890   h->GetYaxis()->SetTitleColor(1);
2891   h->GetZaxis()->SetTitleColor(1);
2892 }
2893
2894 // ________________________________________________________________________________________________________________________________________________________
2895 void AliAnalysisTaskFragmentationFunction::GetJetTracksPointing(TList* inputlist, TList* outputlist, AliAODJet* jet, const Double_t radius,Double_t& sumPt)
2896 {
2897   // fill list of tracks in cone around jet axis  
2898
2899   sumPt = 0;
2900
2901   Double_t jetMom[3];
2902   jet->PxPyPz(jetMom);
2903   TVector3 jet3mom(jetMom);
2904
2905   for (Int_t itrack=0; itrack<inputlist->GetSize(); itrack++){
2906
2907     AliVParticle* track = dynamic_cast<AliVParticle*>(inputlist->At(itrack));
2908
2909     Double_t trackMom[3];
2910     track->PxPyPz(trackMom);
2911     TVector3 track3mom(trackMom);
2912
2913     Double_t dR = jet3mom.DeltaR(track3mom);
2914
2915     if(dR<radius){
2916
2917       outputlist->Add(track);
2918       
2919       sumPt += track->Pt();
2920     }
2921   }
2922   
2923   outputlist->Sort();
2924 }
2925
2926 // ___________________________________________________________________________________________
2927 void AliAnalysisTaskFragmentationFunction::GetJetTracksTrackrefs(TList* list, AliAODJet* jet)
2928 {
2929   // list of jet tracks from trackrefs
2930   
2931   Int_t nTracks = jet->GetRefTracks()->GetEntriesFast();
2932
2933   for (Int_t itrack=0; itrack<nTracks; itrack++) {
2934     
2935     AliVParticle* track = dynamic_cast<AliVParticle*>(jet->GetRefTracks()->At(itrack));
2936     if(!track){
2937       AliError("expected ref track not found ");
2938       continue;
2939     }
2940         
2941     list->Add(track);
2942   }
2943   
2944   list->Sort();
2945 }