AliPhysics  9b6b435 (9b6b435)
AliAnalysisTaskEmcalJetCDF.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-2016, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Author: The ALICE Off-line Project. *
5  * Contributors are mentioned in the code where appropriate. *
6  * *
7  * Permission to use, copy, modify and distribute this software and its *
8  * documentation strictly for non-commercial purposes is hereby granted *
9  * without fee, provided that the above copyright notice appears in all *
10  * copies and that both the copyright notice and this permission notice *
11  * appear in the supporting documentation. The authors make no claims *
12  * about the suitability of this software for any purpose. It is *
13  * provided "as is" without express or implied warranty. *
14  **************************************************************************/
15 
16 #include <cmath>
17 
18 #include <TSystem.h>
19 #include <TObject.h>
20 #include <TObjArray.h>
21 #include <TObjString.h>
22 #include <TString.h>
23 #include <TFile.h>
24 #include <TKey.h>
25 #include <TChain.h>
26 #include <TFileCollection.h>
27 #include <TCollection.h>
28 #include <THashList.h>
29 #include <TRegexp.h>
30 #include <TFileInfo.h>
31 
32 #include <TClonesArray.h>
33 #include <TH1D.h>
34 #include <TH2D.h>
35 #include <TArrayD.h>
36 #include <TString.h>
37 
38 
39 #include <AliVCluster.h>
40 #include <AliVParticle.h>
41 #include <AliLog.h>
42 
43 #include "AliTLorentzVector.h"
44 #include "AliEmcalJet.h"
45 #include "AliRhoParameter.h"
46 #include "AliJetContainer.h"
47 #include "AliParticleContainer.h"
48 #include "AliClusterContainer.h"
49 #include "AliAnalysisManager.h"
50 #include "AliVEventHandler.h"
51 #include "AliAnalysisDataContainer.h"
52 
53 #include "AliEmcalList.h"
54 
56 
57 using std::cout;
58 using std::endl;
59 
61 ClassImp ( AliAnalysisTaskEmcalJetCDF );
63 
69  fHistManager()
70 {}
71 
78  AliAnalysisTaskEmcalJet ( name, kTRUE ),
79  fHistManager(name)
80  {
81  // Standard constructor.
82  SetMakeGeneralHistograms ( kTRUE );
83  }
84 
87 
96  {
97  return kTRUE;
98  }
99 
100 //________________________________________________________________________
102  {
103  TH1::SetDefaultSumw2(kTRUE);
104  TH2::SetDefaultSumw2(kTRUE);
105 
106  namespace CDF = NS_AliAnalysisTaskEmcalJetCDF;
107  TString histname = "", groupname = "", fullgroupname = "";
108 
109  AliJetContainer* jetCont = NULL;
110  TIter next(&fJetCollArray);
111  while ( (jetCont = static_cast<AliJetContainer*>(next())) ) {
112  //##### EARLY VALIDITY CHECKS - BAIL OUT FAST
113 
114  // get particles connected to jets
115  AliParticleContainer* fTracksCont = jetCont->GetParticleContainer();
116  if (!fTracksCont) { std::cout << "********* JET CONTAINER WITHOUT TRACKS CONTAINER *********" << std::endl; continue; }
117  TClonesArray* fTracksContArray = fTracksCont->GetArray();
118 
119  // Number of Jets found in event - accepted cuts applied by JetContainer
120  Int_t fNJets_accepted = jetCont->GetNJets();
121 
122  // Multiplicity in event - accepted tracks in tracks container
123  Int_t fNaccPart = fTracksCont->GetNAcceptedParticles();
124 
125  // protection
126  if ( ( fNJets_accepted < 1 ) || ( fNaccPart < 1 ) ) {
127  if ( fDebug > 1 ) { std::cout << "accepted (fNJets || fNPart) == 0" << std::endl; }
128  continue;
129  }
130 
131  if ( fDebug > 1 ) { std::cout << "fNJets = " << fNJets_accepted << " ; fNPart = " << fNaccPart << std::endl; }
132 
133  // get jet1 - if there is no leading jet there is no point to continue
134  AliEmcalJet* jet1 = jetCont->GetLeadingJet(); // internaly checked for AcceptedJet
135  if (!jet1) {
136  if ( fDebug > 1 ) { Printf ( "Jet1 not found (did not survive cuts?)\n" ); }
137  continue;
138  }
139 
140 //######################################################################################################
141  groupname = jetCont->GetName();
142 
143  Double_t jet_pt_min = jetCont->GetMinPt();
144  Double_t jet_pt_max = jetCont->GetMaxPt();
145 
146  TString jetstrmin = TString::Itoa((Int_t)jet_pt_min,10);
147  TString jetstrmax = TString::Itoa((Int_t)jet_pt_max,10);
148 
149  // add to groupname the min,max pt cuts of jets in the container
150  groupname = groupname + "_" + "ptbin" + "_" + jetstrmin + "_" + jetstrmax;
151 
152 //######################################################################################################
153 // Get histo pointers from Hist Manager
154  histname = TString::Format("%s/histo1_%d", groupname.Data(), fCentBin);
155  TH1D* fH1 = (TH1D*)GetHistogram(histname.Data());
156 
157  histname = TString::Format("%s/histo2_%d", groupname.Data(), fCentBin);
158  TH1D* fH2 = (TH1D*)GetHistogram(histname.Data());
159 
160  histname = TString::Format("%s/histo3_%d", groupname.Data(), fCentBin);
161  TH1D* fH3 = (TH1D*)GetHistogram(histname.Data());
162 
163  histname = TString::Format("%s/histo4_%d", groupname.Data(), fCentBin);
164  TH1D* fH4 = (TH1D*)GetHistogram(histname.Data());
165 
166  histname = TString::Format("%s/histo4c_%d", groupname.Data(), fCentBin);
167  TH1D* fH4c = (TH1D*)GetHistogram(histname.Data());
168 
169  histname = TString::Format("%s/histo5_%d", groupname.Data(), fCentBin);
170  TH1D* fH5 = (TH1D*)GetHistogram(histname.Data());
171 
172  histname = TString::Format("%s/histo6_%d", groupname.Data(), fCentBin);
173  TH1D* fH6 = (TH1D*)GetHistogram(histname.Data());
174 
175  histname = TString::Format("%s/histo6c_%d", groupname.Data(), fCentBin);
176  TH1D* fH6c = (TH1D*)GetHistogram(histname.Data());
177 
178  histname = TString::Format("%s/histo7_%d", groupname.Data(), fCentBin);
179  TH2D* fH7 = (TH2D*)GetHistogram(histname.Data());
180 
181  histname = TString::Format("%s/histo7all_%d", groupname.Data(), fCentBin);
182  TH2D* fH7all = (TH2D*)GetHistogram(histname.Data());
183 //######################################################################################################
184  histname = TString::Format("%s/histo8_%d", groupname.Data(), fCentBin);
185  TH1D* fH8 = (TH1D*)GetHistogram(histname.Data());
186 
187  histname = TString::Format("%s/histo8_all_%d", groupname.Data(), fCentBin);
188  TH1D* fH8_all = (TH1D*)GetHistogram(histname.Data());
189 
190  histname = TString::Format("%s/histo8_p_%d", groupname.Data(), fCentBin);
191  TH1D* fH8_p = (TH1D*)GetHistogram(histname.Data());
192 
193  histname = TString::Format("%s/histo8_all_p_%d", groupname.Data(), fCentBin);
194  TH1D* fH8_all_p = (TH1D*)GetHistogram(histname.Data());
195 
196  histname = TString::Format("%s/histo8_pt_%d", groupname.Data(), fCentBin);
197  TH1D* fH8_pt = (TH1D*)GetHistogram(histname.Data());
198 
199  histname = TString::Format("%s/histo8_all_pt_%d", groupname.Data(), fCentBin);
200  TH1D* fH8_all_pt = (TH1D*)GetHistogram(histname.Data());
201 
202  histname = TString::Format("%s/histo8xi_%d", groupname.Data(), fCentBin);
203  TH1D* fH8xi = (TH1D*)GetHistogram(histname.Data());
204 
205  histname = TString::Format("%s/histo8xi_all_%d", groupname.Data(), fCentBin);
206  TH1D* fH8xi_all = (TH1D*)GetHistogram(histname.Data());
207 
208  histname = TString::Format("%s/histo8xi_p_%d", groupname.Data(), fCentBin);
209  TH1D* fH8xi_p = (TH1D*)GetHistogram(histname.Data());
210 
211  histname = TString::Format("%s/histo8xi_all_p_%d", groupname.Data(), fCentBin);
212  TH1D* fH8xi_all_p = (TH1D*)GetHistogram(histname.Data());
213 
214  histname = TString::Format("%s/histo8xi_pt_%d", groupname.Data(), fCentBin);
215  TH1D* fH8xi_pt = (TH1D*)GetHistogram(histname.Data());
216 
217  histname = TString::Format("%s/histo8xi_all_pt_%d", groupname.Data(), fCentBin);
218  TH1D* fH8xi_all_pt = (TH1D*)GetHistogram(histname.Data());
219 //######################################################################################################
220 
221  histname = TString::Format("%s/histo15_%d", groupname.Data(), fCentBin);
222  TH2D* fH15 = (TH2D*)GetHistogram(histname.Data());
223 
224  histname = TString::Format("%s/histo15_n70_%d", groupname.Data(), fCentBin);
225  TH2D* fH15_n70 = (TH2D*)GetHistogram(histname.Data());
226 
227  histname = TString::Format("%s/histo15_n75_%d", groupname.Data(), fCentBin);
228  TH2D* fH15_n75 = (TH2D*)GetHistogram(histname.Data());
229 
230  histname = TString::Format("%s/histo15_n80_%d", groupname.Data(), fCentBin);
231  TH2D* fH15_n80 = (TH2D*)GetHistogram(histname.Data());
232 
233  histname = TString::Format("%s/histo15_n85_%d", groupname.Data(), fCentBin);
234  TH2D* fH15_n85 = (TH2D*)GetHistogram(histname.Data());
235 
236  histname = TString::Format("%s/histo15_n90_%d", groupname.Data(), fCentBin);
237  TH2D* fH15_n90 = (TH2D*)GetHistogram(histname.Data());
238 
239  histname = TString::Format("%s/histo15_pt70_%d", groupname.Data(), fCentBin);
240  TH2D* fH15_pt70 = (TH2D*)GetHistogram(histname.Data());
241 
242  histname = TString::Format("%s/histo15_pt75_%d", groupname.Data(), fCentBin);
243  TH2D* fH15_pt75 = (TH2D*)GetHistogram(histname.Data());
244 
245  histname = TString::Format("%s/histo15_pt80_%d", groupname.Data(), fCentBin);
246  TH2D* fH15_pt80 = (TH2D*)GetHistogram(histname.Data());
247 
248  histname = TString::Format("%s/histo15_pt85_%d", groupname.Data(), fCentBin);
249  TH2D* fH15_pt85 = (TH2D*)GetHistogram(histname.Data());
250 
251  histname = TString::Format("%s/histo15_pt90_%d", groupname.Data(), fCentBin);
252  TH2D* fH15_pt90 = (TH2D*)GetHistogram(histname.Data());
253 
254  histname = TString::Format("%s/histo15all_%d", groupname.Data(), fCentBin);
255  TH2D* fH15all = (TH2D*)GetHistogram(histname.Data());
256 
257  histname = TString::Format("%s/histo15all_n70_%d", groupname.Data(), fCentBin);
258  TH2D* fH15all_n70 = (TH2D*)GetHistogram(histname.Data());
259 
260  histname = TString::Format("%s/histo15all_n75_%d", groupname.Data(), fCentBin);
261  TH2D* fH15all_n75 = (TH2D*)GetHistogram(histname.Data());
262 
263  histname = TString::Format("%s/histo15all_n80_%d", groupname.Data(), fCentBin);
264  TH2D* fH15all_n80 = (TH2D*)GetHistogram(histname.Data());
265 
266  histname = TString::Format("%s/histo15all_n85_%d", groupname.Data(), fCentBin);
267  TH2D* fH15all_n85 = (TH2D*)GetHistogram(histname.Data());
268 
269  histname = TString::Format("%s/histo15all_n90_%d", groupname.Data(), fCentBin);
270  TH2D* fH15all_n90 = (TH2D*)GetHistogram(histname.Data());
271 
272  histname = TString::Format("%s/histo15all_pt70_%d", groupname.Data(), fCentBin);
273  TH2D* fH15all_pt70 = (TH2D*)GetHistogram(histname.Data());
274 
275  histname = TString::Format("%s/histo15all_pt75_%d", groupname.Data(), fCentBin);
276  TH2D* fH15all_pt75 = (TH2D*)GetHistogram(histname.Data());
277 
278  histname = TString::Format("%s/histo15all_pt80_%d", groupname.Data(), fCentBin);
279  TH2D* fH15all_pt80 = (TH2D*)GetHistogram(histname.Data());
280 
281  histname = TString::Format("%s/histo15all_pt85_%d", groupname.Data(), fCentBin);
282  TH2D* fH15all_pt85 = (TH2D*)GetHistogram(histname.Data());
283 
284  histname = TString::Format("%s/histo15all_pt90_%d", groupname.Data(), fCentBin);
285  TH2D* fH15all_pt90 = (TH2D*)GetHistogram(histname.Data());
286 //######################################################################################################
287 
288  histname = TString::Format("%s/histo20_%d", groupname.Data(), fCentBin);
289  TH1D* fH20 = (TH1D*)GetHistogram(histname.Data());
290 
291  histname = TString::Format("%s/histo20_n70_%d", groupname.Data(), fCentBin);
292  TH1D* fH20_n70 = (TH1D*)GetHistogram(histname.Data());
293 
294  histname = TString::Format("%s/histo20_n75_%d", groupname.Data(), fCentBin);
295  TH1D* fH20_n75 = (TH1D*)GetHistogram(histname.Data());
296 
297  histname = TString::Format("%s/histo20_n80_%d", groupname.Data(), fCentBin);
298  TH1D* fH20_n80 = (TH1D*)GetHistogram(histname.Data());
299 
300  histname = TString::Format("%s/histo20_n85_%d", groupname.Data(), fCentBin);
301  TH1D* fH20_n85 = (TH1D*)GetHistogram(histname.Data());
302 
303  histname = TString::Format("%s/histo20_n90_%d", groupname.Data(), fCentBin);
304  TH1D* fH20_n90 = (TH1D*)GetHistogram(histname.Data());
305 
306  histname = TString::Format("%s/histo20_pt70_%d", groupname.Data(), fCentBin);
307  TH1D* fH20_pt70 = (TH1D*)GetHistogram(histname.Data());
308 
309  histname = TString::Format("%s/histo20_pt75_%d", groupname.Data(), fCentBin);
310  TH1D* fH20_pt75 = (TH1D*)GetHistogram(histname.Data());
311 
312  histname = TString::Format("%s/histo20_pt80_%d", groupname.Data(), fCentBin);
313  TH1D* fH20_pt80 = (TH1D*)GetHistogram(histname.Data());
314 
315  histname = TString::Format("%s/histo20_pt85_%d", groupname.Data(), fCentBin);
316  TH1D* fH20_pt85 = (TH1D*)GetHistogram(histname.Data());
317 
318  histname = TString::Format("%s/histo20_pt90_%d", groupname.Data(), fCentBin);
319  TH1D* fH20_pt90 = (TH1D*)GetHistogram(histname.Data());
320 
321  histname = TString::Format("%s/histo20all_%d", groupname.Data(), fCentBin);
322  TH1D* fH20all = (TH1D*)GetHistogram(histname.Data());
323 
324  histname = TString::Format("%s/histo20all_n70_%d", groupname.Data(), fCentBin);
325  TH1D* fH20all_n70 = (TH1D*)GetHistogram(histname.Data());
326 
327  histname = TString::Format("%s/histo20all_n75_%d", groupname.Data(), fCentBin);
328  TH1D* fH20all_n75 = (TH1D*)GetHistogram(histname.Data());
329 
330  histname = TString::Format("%s/histo20all_n80_%d", groupname.Data(), fCentBin);
331  TH1D* fH20all_n80 = (TH1D*)GetHistogram(histname.Data());
332 
333  histname = TString::Format("%s/histo20all_n85_%d", groupname.Data(), fCentBin);
334  TH1D* fH20all_n85 = (TH1D*)GetHistogram(histname.Data());
335 
336  histname = TString::Format("%s/histo20all_n90_%d", groupname.Data(), fCentBin);
337  TH1D* fH20all_n90 = (TH1D*)GetHistogram(histname.Data());
338 
339  histname = TString::Format("%s/histo20all_pt70_%d", groupname.Data(), fCentBin);
340  TH1D* fH20all_pt70 = (TH1D*)GetHistogram(histname.Data());
341 
342  histname = TString::Format("%s/histo20all_pt75_%d", groupname.Data(), fCentBin);
343  TH1D* fH20all_pt75 = (TH1D*)GetHistogram(histname.Data());
344 
345  histname = TString::Format("%s/histo20all_pt80_%d", groupname.Data(), fCentBin);
346  TH1D* fH20all_pt80 = (TH1D*)GetHistogram(histname.Data());
347 
348  histname = TString::Format("%s/histo20all_pt85_%d", groupname.Data(), fCentBin);
349  TH1D* fH20all_pt85 = (TH1D*)GetHistogram(histname.Data());
350 
351  histname = TString::Format("%s/histo20all_pt90_%d", groupname.Data(), fCentBin);
352  TH1D* fH20all_pt90 = (TH1D*)GetHistogram(histname.Data());
353 //######################################################################################################
354 
355  histname = TString::Format("%s/histo_g_%d", groupname.Data(), fCentBin);
356  TH1D* fHg = (TH1D*)GetHistogram(histname.Data());
357 
358  histname = TString::Format("%s/histo_g_n70_%d", groupname.Data(), fCentBin);
359  TH1D* fHg_n70 = (TH1D*)GetHistogram(histname.Data());
360 
361  histname = TString::Format("%s/histo_g_n75_%d", groupname.Data(), fCentBin);
362  TH1D* fHg_n75 = (TH1D*)GetHistogram(histname.Data());
363 
364  histname = TString::Format("%s/histo_g_n80_%d", groupname.Data(), fCentBin);
365  TH1D* fHg_n80 = (TH1D*)GetHistogram(histname.Data());
366 
367  histname = TString::Format("%s/histo_g_n85_%d", groupname.Data(), fCentBin);
368  TH1D* fHg_n85 = (TH1D*)GetHistogram(histname.Data());
369 
370  histname = TString::Format("%s/histo_g_n90_%d", groupname.Data(), fCentBin);
371  TH1D* fHg_n90 = (TH1D*)GetHistogram(histname.Data());
372 
373  histname = TString::Format("%s/histo_g_pt70_%d", groupname.Data(), fCentBin);
374  TH1D* fHg_pt70 = (TH1D*)GetHistogram(histname.Data());
375 
376  histname = TString::Format("%s/histo_g_pt75_%d", groupname.Data(), fCentBin);
377  TH1D* fHg_pt75 = (TH1D*)GetHistogram(histname.Data());
378 
379  histname = TString::Format("%s/histo_g_pt80_%d", groupname.Data(), fCentBin);
380  TH1D* fHg_pt80 = (TH1D*)GetHistogram(histname.Data());
381 
382  histname = TString::Format("%s/histo_g_pt85_%d", groupname.Data(), fCentBin);
383  TH1D* fHg_pt85 = (TH1D*)GetHistogram(histname.Data());
384 
385  histname = TString::Format("%s/histo_g_pt90_%d", groupname.Data(), fCentBin);
386  TH1D* fHg_pt90 = (TH1D*)GetHistogram(histname.Data());
387 //######################################################################################################
388 
389  histname = TString::Format("%s/histo_ptd_%d", groupname.Data(), fCentBin);
390  TH1D* fHptd = (TH1D*)GetHistogram(histname.Data());
391 
392  histname = TString::Format("%s/histo_ptd_n70_%d", groupname.Data(), fCentBin);
393  TH1D* fHptd_n70 = (TH1D*)GetHistogram(histname.Data());
394 
395  histname = TString::Format("%s/histo_ptd_n75_%d", groupname.Data(), fCentBin);
396  TH1D* fHptd_n75 = (TH1D*)GetHistogram(histname.Data());
397 
398  histname = TString::Format("%s/histo_ptd_n80_%d", groupname.Data(), fCentBin);
399  TH1D* fHptd_n80 = (TH1D*)GetHistogram(histname.Data());
400 
401  histname = TString::Format("%s/histo_ptd_n85_%d", groupname.Data(), fCentBin);
402  TH1D* fHptd_n85 = (TH1D*)GetHistogram(histname.Data());
403 
404  histname = TString::Format("%s/histo_ptd_n90_%d", groupname.Data(), fCentBin);
405  TH1D* fHptd_n90 = (TH1D*)GetHistogram(histname.Data());
406 
407  histname = TString::Format("%s/histo_ptd_pt70_%d", groupname.Data(), fCentBin);
408  TH1D* fHptd_pt70 = (TH1D*)GetHistogram(histname.Data());
409 
410  histname = TString::Format("%s/histo_ptd_pt75_%d", groupname.Data(), fCentBin);
411  TH1D* fHptd_pt75 = (TH1D*)GetHistogram(histname.Data());
412 
413  histname = TString::Format("%s/histo_ptd_pt80_%d", groupname.Data(), fCentBin);
414  TH1D* fHptd_pt80 = (TH1D*)GetHistogram(histname.Data());
415 
416  histname = TString::Format("%s/histo_ptd_pt85_%d", groupname.Data(), fCentBin);
417  TH1D* fHptd_pt85 = (TH1D*)GetHistogram(histname.Data());
418 
419  histname = TString::Format("%s/histo_ptd_pt90_%d", groupname.Data(), fCentBin);
420  TH1D* fHptd_pt90 = (TH1D*)GetHistogram(histname.Data());
421 //######################################################################################################
422 
423  histname = TString::Format("%s/histo_Rjt_%d", groupname.Data(), fCentBin);
424  TH2D* fH_Rjt = (TH2D*)GetHistogram(histname.Data());
425 
426  histname = TString::Format("%s/histo_Rjt_n70_%d", groupname.Data(), fCentBin);
427  TH2D* fH_Rjt_n70 = (TH2D*)GetHistogram(histname.Data());
428 
429  histname = TString::Format("%s/histo_Rjt_n75_%d", groupname.Data(), fCentBin);
430  TH2D* fH_Rjt_n75 = (TH2D*)GetHistogram(histname.Data());
431 
432  histname = TString::Format("%s/histo_Rjt_n80_%d", groupname.Data(), fCentBin);
433  TH2D* fH_Rjt_n80 = (TH2D*)GetHistogram(histname.Data());
434 
435  histname = TString::Format("%s/histo_Rjt_n85_%d", groupname.Data(), fCentBin);
436  TH2D* fH_Rjt_n85 = (TH2D*)GetHistogram(histname.Data());
437 
438  histname = TString::Format("%s/histo_Rjt_n90_%d", groupname.Data(), fCentBin);
439  TH2D* fH_Rjt_n90 = (TH2D*)GetHistogram(histname.Data());
440 
441  histname = TString::Format("%s/histo_Rjt_pt70_%d", groupname.Data(), fCentBin);
442  TH2D* fH_Rjt_pt70 = (TH2D*)GetHistogram(histname.Data());
443 
444  histname = TString::Format("%s/histo_Rjt_pt75_%d", groupname.Data(), fCentBin);
445  TH2D* fH_Rjt_pt75 = (TH2D*)GetHistogram(histname.Data());
446 
447  histname = TString::Format("%s/histo_Rjt_pt80_%d", groupname.Data(), fCentBin);
448  TH2D* fH_Rjt_pt80 = (TH2D*)GetHistogram(histname.Data());
449 
450  histname = TString::Format("%s/histo_Rjt_pt85_%d", groupname.Data(), fCentBin);
451  TH2D* fH_Rjt_pt85 = (TH2D*)GetHistogram(histname.Data());
452 
453  histname = TString::Format("%s/histo_Rjt_pt90_%d", groupname.Data(), fCentBin);
454  TH2D* fH_Rjt_pt90 = (TH2D*)GetHistogram(histname.Data());
455 //######################################################################################################
456 
457  histname = TString::Format("%s/histo_jt_%d", groupname.Data(), fCentBin);
458  TH1D* fH_jt = (TH1D*)GetHistogram(histname.Data());
459 
460  histname = TString::Format("%s/histo_jt_n70_%d", groupname.Data(), fCentBin);
461  TH1D* fH_jt_n70 = (TH1D*)GetHistogram(histname.Data());
462 
463  histname = TString::Format("%s/histo_jt_n75_%d", groupname.Data(), fCentBin);
464  TH1D* fH_jt_n75 = (TH1D*)GetHistogram(histname.Data());
465 
466  histname = TString::Format("%s/histo_jt_n80_%d", groupname.Data(), fCentBin);
467  TH1D* fH_jt_n80 = (TH1D*)GetHistogram(histname.Data());
468 
469  histname = TString::Format("%s/histo_jt_n85_%d", groupname.Data(), fCentBin);
470  TH1D* fH_jt_n85 = (TH1D*)GetHistogram(histname.Data());
471 
472  histname = TString::Format("%s/histo_jt_n90_%d", groupname.Data(), fCentBin);
473  TH1D* fH_jt_n90 = (TH1D*)GetHistogram(histname.Data());
474 
475  histname = TString::Format("%s/histo_jt_pt70_%d", groupname.Data(), fCentBin);
476  TH1D* fH_jt_pt70 = (TH1D*)GetHistogram(histname.Data());
477 
478  histname = TString::Format("%s/histo_jt_pt75_%d", groupname.Data(), fCentBin);
479  TH1D* fH_jt_pt75 = (TH1D*)GetHistogram(histname.Data());
480 
481  histname = TString::Format("%s/histo_jt_pt80_%d", groupname.Data(), fCentBin);
482  TH1D* fH_jt_pt80 = (TH1D*)GetHistogram(histname.Data());
483 
484  histname = TString::Format("%s/histo_jt_pt85_%d", groupname.Data(), fCentBin);
485  TH1D* fH_jt_pt85 = (TH1D*)GetHistogram(histname.Data());
486 
487  histname = TString::Format("%s/histo_jt_pt90_%d", groupname.Data(), fCentBin);
488  TH1D* fH_jt_pt90 = (TH1D*)GetHistogram(histname.Data());
489 
490 
491 //######################################################################################################
492  // get clusters connected to jets
493  AliClusterContainer* fCaloClustersCont = jetCont->GetClusterContainer();
494  // accepted clusters in cluster container
495  Int_t fNaccClus = -1;
496  if (fCaloClustersCont) { fNaccClus = fCaloClustersCont->GetNAcceptedClusters(); }
497 
498  fH5->Fill ( fNJets_accepted ); // Distribution of jets in events;
499 
500  UShort_t counter_part = 0; Double_t counter_pt = 0.; // counter for npart and pt recording
501  UShort_t jet_n90 = -99 ; Double_t jet_pt90 = -99.99 ;
502  UShort_t jet_n85 = -99 ; Double_t jet_pt85 = -99.99 ;
503  UShort_t jet_n80 = -99 ; Double_t jet_pt80 = -99.99 ;
504  UShort_t jet_n75 = -99 ; Double_t jet_pt75 = -99.99 ;
505  UShort_t jet_n70 = -99 ; Double_t jet_pt70 = -99.99 ;
506 
507  // variables used to compute g and ptD
508  Double_t g_tot = 0.; Double_t sum_part_pt_tot = 0.; Double_t sum_part_pt2_tot = 0.;
509  Double_t g_n90 = 0.; Double_t sum_part_pt_n90 = 0.; Double_t sum_part_pt2_n90 = 0.;
510  Double_t g_n85 = 0.; Double_t sum_part_pt_n85 = 0.; Double_t sum_part_pt2_n85 = 0.;
511  Double_t g_n80 = 0.; Double_t sum_part_pt_n80 = 0.; Double_t sum_part_pt2_n80 = 0.;
512  Double_t g_n75 = 0.; Double_t sum_part_pt_n75 = 0.; Double_t sum_part_pt2_n75 = 0.;
513  Double_t g_n70 = 0.; Double_t sum_part_pt_n70 = 0.; Double_t sum_part_pt2_n70 = 0.;
514  Double_t g_pt90 = 0.; Double_t sum_part_pt_pt90 = 0.; Double_t sum_part_pt2_pt90 = 0.;
515  Double_t g_pt85 = 0.; Double_t sum_part_pt_pt85 = 0.; Double_t sum_part_pt2_pt85 = 0.;
516  Double_t g_pt80 = 0.; Double_t sum_part_pt_pt80 = 0.; Double_t sum_part_pt2_pt80 = 0.;
517  Double_t g_pt75 = 0.; Double_t sum_part_pt_pt75 = 0.; Double_t sum_part_pt2_pt75 = 0.;
518  Double_t g_pt70 = 0.; Double_t sum_part_pt_pt70 = 0.; Double_t sum_part_pt2_pt70 = 0.;
519 
520  // **************************************************************
521  // LEADING JETS
522  // **************************************************************
523  if ( jet1 ) {
524  if ( fDebug > 1 ) { std::cout << "+++++++++++++++++>>>>>>>>> Leading jet found" << std::endl; jet1->Print(); }
525 
526  // vector of sorted indexes of particles in leading jet - jet1 : sorting by p_T jet constituents
527  std::vector< int > jet1_sorted_idxvec = jet1->GetPtSortedTrackConstituentIndexes ( fTracksContArray );
528 
529  Double_t jet1_pt = jet1->Pt();
530  UShort_t jet1_npart = jet1->GetNumberOfTracks();
531  UShort_t jet1_nconst = jet1->GetNumberOfConstituents();
532 
533  UShort_t jet1_n90 = ( UShort_t ) ( 0.90 * jet1_npart );
534  Double_t jet1_pt90 = 0.90 * jet1_pt;
535 
536  UShort_t jet1_n85 = ( UShort_t ) ( 0.85 * jet1_npart );
537  Double_t jet1_pt85 = 0.85 * jet1_pt;
538 
539  UShort_t jet1_n80 = ( UShort_t ) ( 0.80 * jet1_npart );
540  Double_t jet1_pt80 = 0.80 * jet1_pt;
541 
542  UShort_t jet1_n75 = ( UShort_t ) ( 0.75 * jet1_npart );
543  Double_t jet1_pt75 = 0.75 * jet1_pt;
544 
545  UShort_t jet1_n70 = ( UShort_t ) ( 0.70 * jet1_npart );
546  Double_t jet1_pt70 = 0.70 * jet1_pt;
547 
548  fH6->Fill ( jet1_npart ); // Multiplicity of jet1 - charged tracks
549  fH6c->Fill ( jet1_nconst ); // Multiplicity of jet1 - all constituents
550  fH7->Fill ( jet1_pt, jet1_nconst ); // N(jet) vs P_{T}(jet1)
551 
552  Int_t track_idx = -999; // index variable for tracks
553  counter_part = 0; counter_pt = 0.; // reset counters
554  //___________________________________________________________________________
555  // parsing tracks of jet1 (leading jet) in decreasing order of Pt
556  for (std::size_t i = 0; i < jet1_npart; i++ ) {
557  track_idx = jet1_sorted_idxvec.at (i);
558  //track = dynamic_cast<AliVParticle*>(fTracksContArray->At( track_idx ));
559  AliVParticle* track = jet1->TrackAt ( track_idx, fTracksContArray );
560  if (!track) { std::cout << "Parsing tracks of jet1 :: track not defined but it should!!!" << std::endl; continue; }
561 
562  Double_t dpart = jet1->DeltaR ( track );
563  Double_t track_pt = track->Pt();
564 
565  fH8->Fill ( jet1->GetZ ( track ) ); // Momentum distribution for leading jet (FF)
566  fH8xi->Fill ( jet1->GetXi ( track ) ); // Momentum distribution for leading jet (FF) xi
567 
568  Double_t z_p = CDF::Z_ptot(jet1, track);
569  fH8_p->Fill ( z_p ); // Momentum distribution for jets (FF)
570  fH8xi_p->Fill ( CDF::Xi (z_p) ); // Momentum distribution for jets (FF) xi
571 
572  Double_t z_pt = CDF::Z_pt(jet1, track);
573  fH8_pt->Fill ( z_pt ); // Momentum distribution for jets (FF)
574  fH8xi_pt->Fill ( CDF::Xi (z_pt) ); // Momentum distribution for jets (FF) xi
575 
576  fH15->Fill ( dpart, track_pt, track_pt ); // <p_{T}> track vs the Distance R from Jet1
577  fH20->Fill ( dpart ); // Distribution of R in leading jet
578 
579  // fill histograms for 70% of particles with highest pt
580  if ( counter_part <= jet1_n70 ) {
581  fH15_n70->Fill ( dpart, track_pt, track_pt ); // <p_{T}> track vs the Distance R from Jet1 - 80% of particles
582  fH20_n70->Fill ( dpart ); // Distribution of R in leading jet
583  }
584  // fill histograms for 75% of particles with highest pt
585  if ( counter_part <= jet1_n75 ) {
586  fH15_n75->Fill ( dpart, track_pt, track_pt ); // <p_{T}> track vs the Distance R from Jet1 - 80% of particles
587  fH20_n75->Fill ( dpart ); // Distribution of R in leading jet
588  }
589  // fill histograms for 80% of particles with highest pt
590  if ( counter_part <= jet1_n80 ) {
591  fH15_n80->Fill ( dpart, track_pt, track_pt ); // <p_{T}> track vs the Distance R from Jet1 - 80% of particles
592  fH20_n80->Fill ( dpart ); // Distribution of R in leading jet
593  }
594  // fill histograms for 85% of particles with highest pt
595  if ( counter_part <= jet1_n85 ) {
596  fH15_n85->Fill ( dpart, track_pt, track_pt ); // <p_{T}> track vs the Distance R from Jet1 - 80% of particles
597  fH20_n85->Fill ( dpart ); // Distribution of R in leading jet
598  }
599  // fill histograms for 90% of particles with highest pt
600  if ( counter_part <= jet1_n90 ) {
601  fH15_n90->Fill ( dpart, track_pt, track_pt ); // <p_{T}> track vs the Distance R from Jet1 - 80% of particles
602  fH20_n90->Fill ( dpart ); // Distribution of R in leading jet
603  }
604 
605  // fill histograms for particles that have first 70% of pt
606  if ( counter_pt <= jet1_pt70 ) {
607  fH15_pt70->Fill ( dpart, track_pt, track_pt ); // <p_{T}> track vs the Distance R from Jet1 - 80% of pt
608  fH20_pt70->Fill ( dpart ); // Distribution of R in leading jet
609  }
610  // fill histograms for particles that have first 75% of pt
611  if ( counter_pt <= jet1_pt75 ) {
612  fH15_pt75->Fill ( dpart, track_pt, track_pt ); // <p_{T}> track vs the Distance R from Jet1 - 80% of pt
613  fH20_pt75->Fill ( dpart ); // Distribution of R in leading jet
614  }
615  // fill histograms for particles that have first 80% of pt
616  if ( counter_pt <= jet1_pt80 ) {
617  fH15_pt80->Fill ( dpart, track_pt, track_pt ); // <p_{T}> track vs the Distance R from Jet1 - 80% of pt
618  fH20_pt80->Fill ( dpart ); // Distribution of R in leading jet
619  }
620  // fill histograms for particles that have first 80% of pt
621  if ( counter_pt <= jet1_pt85 ) {
622  fH15_pt85->Fill ( dpart, track_pt, track_pt ); // <p_{T}> track vs the Distance R from Jet1 - 80% of pt
623  fH20_pt85->Fill ( dpart ); // Distribution of R in leading jet
624  }
625  // fill histograms for particles that have first 80% of pt
626  if ( counter_pt <= jet1_pt90 ) {
627  fH15_pt90->Fill ( dpart, track_pt, track_pt ); // <p_{T}> track vs the Distance R from Jet1 - 80% of pt
628  fH20_pt90->Fill ( dpart ); // Distribution of R in leading jet
629  }
630 
631  ++counter_part; counter_pt += track_pt;
632  } // end of loop over jet1 tracks
633  } // end of jet1 (leading jet) processing
634 
635  // **************************************************************
636  // ALL JETS
637  // **************************************************************
638  Double_t jet_pt = 0. ; UShort_t jet_npart = 0; UShort_t jet_nconst = 0;
639 
640  // loop over all jets
641  for( auto jet : jetCont->accepted()) {
642  if (!jet) { continue; }
643 
644  // vector of sorted indexes of particles in jet
645  std::vector< int > jet_sorted_idxvec ;
646 
647  Int_t track_idx = -999; // index variable for tracks
648  jet_pt = 0. ; jet_npart = 0; jet_nconst = 0;
649  counter_part = 0; counter_pt = 0.; // reset counters
650 
651  // jet : Sorting by p_T jet constituents
652  jet_sorted_idxvec.clear();
653  jet_sorted_idxvec = jet->GetPtSortedTrackConstituentIndexes ( fTracksContArray );
654 
655  jet_pt = jet->Pt();
656  jet_npart = jet->GetNumberOfTracks();
657  jet_nconst = jet->GetNumberOfConstituents();
658 
659  // variables for g and pdt computations
660  g_tot = 0.; sum_part_pt_tot = 0.; sum_part_pt2_tot = 0.;
661  g_n90 = 0.; sum_part_pt_n90 = 0.; sum_part_pt2_n90 = 0.;
662  g_n85 = 0.; sum_part_pt_n85 = 0.; sum_part_pt2_n85 = 0.;
663  g_n80 = 0.; sum_part_pt_n80 = 0.; sum_part_pt2_n80 = 0.;
664  g_n75 = 0.; sum_part_pt_n75 = 0.; sum_part_pt2_n75 = 0.;
665  g_n70 = 0.; sum_part_pt_n70 = 0.; sum_part_pt2_n70 = 0.;
666  g_pt90 = 0.; sum_part_pt_pt90 = 0.; sum_part_pt2_pt90 = 0.;
667  g_pt85 = 0.; sum_part_pt_pt85 = 0.; sum_part_pt2_pt85 = 0.;
668  g_pt80 = 0.; sum_part_pt_pt80 = 0.; sum_part_pt2_pt80 = 0.;
669  g_pt75 = 0.; sum_part_pt_pt75 = 0.; sum_part_pt2_pt75 = 0.;
670  g_pt70 = 0.; sum_part_pt_pt70 = 0.; sum_part_pt2_pt70 = 0.;
671 
672  // sentinels for the pt tail cuts
673  jet_n90 = ( UShort_t ) ( 0.9 * jet_npart );
674  jet_pt90 = 0.9 * jet_pt;
675 
676  jet_n85 = ( UShort_t ) ( 0.85 * jet_npart );
677  jet_pt85 = 0.85 * jet_pt;
678 
679  jet_n80 = ( UShort_t ) ( 0.8 * jet_npart );
680  jet_pt80 = 0.8 * jet_pt;
681 
682  jet_n75 = ( UShort_t ) ( 0.75 * jet_npart );
683  jet_pt75 = 0.75 * jet_pt;
684 
685  jet_n70 = ( UShort_t ) ( 0.7 * jet_npart );
686  jet_pt70 = 0.7 * jet_pt;
687 
688  fH1->Fill ( jet_pt ); // Pt distribution of jets
689  fH2->Fill ( jet->Eta() ); // Eta distribution of jets
690  fH3->Fill ( jet->Phi() ); // Phi distribution of jets
691  fH4->Fill ( jet_npart ); // Multiplicity of jets
692  fH4c->Fill ( jet_nconst ); // Multiplicity of jets - all constituents
693  fH7all->Fill ( jet_pt, jet_nconst ); // N(jet) vs P_{T} - all jets
694 
695  // parsing all jet tracks
696  for (std::size_t i = 0; i < jet_npart; i++ ) {
697  track_idx = jet_sorted_idxvec.at (i);
698  //track = dynamic_cast<AliVParticle*>(fTracksContArray->At( track_idx ));
699  AliVParticle* track = jet->TrackAt ( track_idx, fTracksContArray );
700  if (!track) { std::cout << "Parsing tracks of jets :: track not defined but it should!!!" << std::endl; continue; }
701 
702  Double_t dpart = jet->DeltaR ( track );
703  Double_t track_pt = track->Pt();
704  Double_t jt = CDF::Perp (*track, *jet);
705 
706  fH8_all->Fill ( jet->GetZ ( track ) ); // Momentum distribution for jets (FF)
707  fH8xi_all->Fill ( jet->GetXi ( track ) ); // Momentum distribution for jets (FF) xi
708 
709  Double_t z_p = CDF::Z_ptot(jet,track);
710  fH8_all_p->Fill ( z_p ); // Momentum distribution for jets (FF)
711  fH8xi_all_p->Fill ( CDF::Xi (z_p) ); // Momentum distribution for jets (FF) xi
712 
713  Double_t z_pt = CDF::Z_pt(jet,track);
714  fH8_all_pt->Fill ( z_pt ); // Momentum distribution for jets (FF)
715  fH8xi_all_pt->Fill ( CDF::Xi (z_pt) ); // Momentum distribution for jets (FF) xi
716 
717  fH15all->Fill ( dpart, track_pt, track_pt ); // p_{T} track vs the Distance R from jet
718  fH20all->Fill ( dpart ); // Distribution of R in leading jet
719 
720  fH_Rjt->Fill ( dpart, jt, jt ); // jt track vs dR weighted with jt
721  fH_jt->Fill ( jt ); // jt track vs dR
722 
723  // computing components for g and ptD in the jet tracks loop
724  g_tot += (track_pt * dpart)/jet_pt;
725  sum_part_pt_tot += TMath::Abs(track_pt);
726  sum_part_pt2_tot += TMath::Power( track_pt, 2 );
727 
728  //#############################################################################################
729  if ( counter_part <= jet_n90 ) {
730  fH15all_n90->Fill ( dpart, track_pt ); // p_{T} track vs the Distance R from Jet - 80% of particles
731  fH20all_n90->Fill ( dpart );
732  fH_Rjt_n90->Fill ( dpart, jt, jt );
733  fH_jt_n90->Fill ( jt ); // jt track vs dR
734 
735  // computing components for g and ptD in the jet tracks loop
736  g_n90 += (track_pt * dpart)/jet_pt;
737  sum_part_pt_n90 += track_pt;
738  sum_part_pt2_n90 += TMath::Power( track_pt, 2 );
739  }
740 
741  if ( counter_pt <= jet_pt90 ) {
742  fH15all_pt90->Fill ( dpart, track_pt ); // p_{T} track vs the Distance R from Jet - 80% of pt
743  fH20all_pt90->Fill ( dpart );
744  fH_Rjt_pt90->Fill ( dpart, jt, jt );
745  fH_jt_pt90->Fill ( jt ); // jt track vs dR
746 
747  // computing components for g and ptD in the jet tracks loop
748  g_pt90 += (track_pt * dpart)/jet_pt;
749  sum_part_pt_pt90 += track_pt;
750  sum_part_pt2_pt90 += TMath::Power( track_pt, 2 );
751  }
752 
753  //#############################################################################################
754  if ( counter_part <= jet_n85 ) {
755  fH15all_n85->Fill ( dpart, track_pt ); // p_{T} track vs the Distance R from Jet - 80% of particles
756  fH20all_n85->Fill ( dpart );
757  fH_Rjt_n85->Fill ( dpart, jt, jt );
758  fH_jt_n85->Fill ( jt ); // jt track vs dR
759 
760  // computing components for g and ptD in the jet tracks loop
761  g_n85 += (track_pt * dpart)/jet_pt;
762  sum_part_pt_n85 += track_pt;
763  sum_part_pt2_n85 += TMath::Power( track_pt, 2 );
764  }
765 
766  if ( counter_pt <= jet_pt85 ) {
767  fH15all_pt85->Fill ( dpart, track_pt ); // p_{T} track vs the Distance R from Jet - 80% of pt
768  fH20all_pt85->Fill ( dpart );
769  fH_Rjt_pt85->Fill ( dpart, jt, jt );
770  fH_jt_pt85->Fill ( jt ); // jt track vs dR
771 
772  // computing components for g and ptD in the jet tracks loop
773  g_pt85 += (track_pt * dpart)/jet_pt;
774  sum_part_pt_pt85 += track_pt;
775  sum_part_pt2_pt85 += TMath::Power( track_pt, 2 );
776  }
777 
778  //#############################################################################################
779  if ( counter_part <= jet_n80 ) {
780  fH15all_n80->Fill ( dpart, track_pt ); // p_{T} track vs the Distance R from Jet - 80% of particles
781  fH20all_n80->Fill ( dpart );
782  fH_Rjt_n80->Fill ( dpart, jt, jt );
783  fH_jt_n80->Fill ( jt ); // jt track vs dR
784 
785  // computing components for g and ptD in the jet tracks loop
786  g_n80 += (track_pt * dpart)/jet_pt;
787  sum_part_pt_n80 += track_pt;
788  sum_part_pt2_n80 += TMath::Power( track_pt, 2 );
789  }
790 
791  if ( counter_pt <= jet_pt80 ) {
792  fH15all_pt80->Fill ( dpart, track_pt ); // p_{T} track vs the Distance R from Jet - 80% of pt
793  fH20all_pt80->Fill ( dpart );
794  fH_Rjt_pt80->Fill ( dpart, jt, jt );
795  fH_jt_pt80->Fill ( jt ); // jt track vs dR
796 
797  // computing components for g and ptD in the jet tracks loop
798  g_pt80 += (track_pt * dpart)/jet_pt;
799  sum_part_pt_pt80 += track_pt;
800  sum_part_pt2_pt80 += TMath::Power( track_pt, 2 );
801  }
802 
803  //#############################################################################################
804  if ( counter_part <= jet_n75 ) {
805  fH15all_n75->Fill ( dpart, track_pt ); // p_{T} track vs the Distance R from Jet - 80% of particles
806  fH20all_n75->Fill ( dpart );
807  fH_Rjt_n75->Fill ( dpart, jt, jt );
808  fH_jt_n75->Fill ( jt ); // jt track vs dR
809 
810  // computing components for g and ptD in the jet tracks loop
811  g_n75 += (track_pt * dpart)/jet_pt;
812  sum_part_pt_n75 += track_pt;
813  sum_part_pt2_n75 += TMath::Power( track_pt, 2 );
814  }
815 
816  if ( counter_pt <= jet_pt75 ) {
817  fH15all_pt75->Fill ( dpart, track_pt ); // p_{T} track vs the Distance R from Jet - 80% of pt
818  fH20all_pt75->Fill ( dpart );
819  fH_Rjt_pt75->Fill ( dpart, jt, jt );
820  fH_jt_pt75->Fill ( jt ); // jt track vs dR
821 
822  // computing components for g and ptD in the jet tracks loop
823  g_pt75 += (track_pt * dpart)/jet_pt;
824  sum_part_pt_pt75 += track_pt;
825  sum_part_pt2_pt75 += TMath::Power( track_pt, 2 );
826  }
827 
828  //#############################################################################################
829  if ( counter_part <= jet_n70 ) {
830  fH15all_n70->Fill ( dpart, track_pt ); // p_{T} track vs the Distance R from Jet - 80% of particles
831  fH20all_n70->Fill ( dpart );
832  fH_Rjt_n70->Fill ( dpart, jt, jt );
833  fH_jt_n70->Fill ( jt ); // jt track vs dR
834 
835  // computing components for g and ptD in the jet tracks loop
836  g_n70 += (track_pt * dpart)/jet_pt;
837  sum_part_pt_n70 += track_pt;
838  sum_part_pt2_n70 += TMath::Power( track_pt, 2 );
839  }
840 
841  if ( counter_pt <= jet_pt70 ) {
842  fH15all_pt70->Fill ( dpart, track_pt ); // p_{T} track vs the Distance R from Jet - 80% of pt
843  fH20all_pt70->Fill ( dpart );
844  fH_Rjt_pt70->Fill ( dpart, jt, jt );
845  fH_jt_pt70->Fill ( jt ); // jt track vs dR
846 
847  // computing components for g and ptD in the jet tracks loop
848  g_pt70 += (track_pt * dpart)/jet_pt;
849  sum_part_pt_pt70 += track_pt;
850  sum_part_pt2_pt70 += TMath::Power( track_pt, 2 );
851  }
852  ++counter_part; counter_pt += track_pt;
853  } // end of loop over jet tracks
854 
855  fHg->Fill ( g_tot );
856  fHg_n70->Fill ( g_n70 ); fHg_pt70->Fill ( g_pt70 );
857  fHg_n75->Fill ( g_n75 ); fHg_pt75->Fill ( g_pt75 );
858  fHg_n80->Fill ( g_n80 ); fHg_pt80->Fill ( g_pt80 );
859  fHg_n85->Fill ( g_n85 ); fHg_pt85->Fill ( g_pt85 );
860  fHg_n90->Fill ( g_n90 ); fHg_pt90->Fill ( g_pt90 );
861 
862  if ( sum_part_pt_tot > 1e-8 )
863  { fHptd->Fill( TMath::Sqrt(sum_part_pt2_tot)/sum_part_pt_tot ); }
864  else
865  { if ( fDebug > 2 ) cout << "sum_part_pt_tot aprox 0!!!!" << endl; }
866 
867  if ( sum_part_pt_n70 > 1e-8 )
868  { fHptd_n70->Fill( TMath::Sqrt(sum_part_pt2_n70)/sum_part_pt_n70 ); }
869  else
870  { if ( fDebug > 2 ) cout << "sum_part_pt_n70 aprox 0!!!!" << endl; }
871 
872  if ( sum_part_pt_n75 > 1e-8 )
873  { fHptd_n75->Fill( TMath::Sqrt(sum_part_pt2_n75)/sum_part_pt_n75 ); }
874  else
875  { if ( fDebug > 2 ) cout << "sum_part_pt_n75 aprox 0!!!!" << endl; }
876 
877  if ( sum_part_pt_n80 > 1e-8 )
878  { fHptd_n80->Fill( TMath::Sqrt(sum_part_pt2_n80)/sum_part_pt_n80 ); }
879  else
880  { if ( fDebug > 2 ) cout << "sum_part_pt_n80 aprox 0!!!!" << endl; }
881 
882  if ( sum_part_pt_n85 > 1e-8 )
883  { fHptd_n85->Fill( TMath::Sqrt(sum_part_pt2_n85)/sum_part_pt_n85 ); }
884  else
885  { if ( fDebug > 2 ) cout << "sum_part_pt_n85 aprox 0!!!!" << endl; }
886 
887  if ( sum_part_pt_n90 > 1e-8 )
888  { fHptd_n90->Fill( TMath::Sqrt(sum_part_pt2_n90)/sum_part_pt_n90 ); }
889  else
890  { if ( fDebug > 2 ) cout << "sum_part_pt_n90 aprox 0!!!!" << endl; }
891 
892  if ( sum_part_pt_pt70 > 1e-8 )
893  { fHptd_pt70->Fill( TMath::Sqrt(sum_part_pt2_pt70)/sum_part_pt_pt70 ); }
894  else
895  { if ( fDebug > 2 ) cout << "sum_part_pt_pt70 aprox 0!!!!" << endl; }
896 
897  if ( sum_part_pt_pt75 > 1e-8 )
898  { fHptd_pt75->Fill( TMath::Sqrt(sum_part_pt2_pt75)/sum_part_pt_pt75 ); }
899  else
900  { if ( fDebug > 2 ) cout << "sum_part_pt_pt75 aprox 0!!!!" << endl; }
901 
902  if ( sum_part_pt_pt80 > 1e-8 )
903  { fHptd_pt80->Fill( TMath::Sqrt(sum_part_pt2_pt80)/sum_part_pt_pt80 ); }
904  else
905  { if ( fDebug > 2 ) cout << "sum_part_pt_pt80 aprox 0!!!!" << endl; }
906 
907  if ( sum_part_pt_pt85 > 1e-8 )
908  { fHptd_pt85->Fill( TMath::Sqrt(sum_part_pt2_pt85)/sum_part_pt_pt85 ); }
909  else
910  { if ( fDebug > 2 ) cout << "sum_part_pt_pt85 aprox 0!!!!" << endl; }
911 
912  if ( sum_part_pt_pt90 > 1e-8 )
913  { fHptd_pt90->Fill( TMath::Sqrt(sum_part_pt2_pt90)/sum_part_pt_pt90 ); }
914  else
915  { if ( fDebug > 2 ) cout << "sum_part_pt_pt90 aprox 0!!!!" << endl; }
916 
917  } // end of loopt over all jets
918  } // end of loop over jet container collection
919 
920  // post data at every processing
921  PostData ( 1, fOutput ); // Post data for ALL output slots > 0 here.
922  return kTRUE;
923  }
924 
925 //________________________________________________________________________
927  {
928  // Create user output.
930 
931  TString histname = "", histtitle = "", groupname = "", fullgroupname = "";
932  AliJetContainer* jetCont = 0;
933  TIter next(&fJetCollArray);
934  while ((jetCont = static_cast<AliJetContainer*>(next())))
935  {
936  groupname = jetCont->GetName();
937 
938  Double_t jet_pt_min = jetCont->GetMinPt();
939  Double_t jet_pt_max = jetCont->GetMaxPt();
940 
941  TString jetstrmin = TString::Itoa((Int_t)jet_pt_min,10);
942  TString jetstrmax = TString::Itoa((Int_t)jet_pt_max,10);
943 
944  // add to groupname the min,max pt cuts of jets in the container
945  groupname = groupname + "_" + "ptbin" + "_" + jetstrmin + "_" + jetstrmax;
946 
947  fHistManager.CreateHistoGroup(groupname);
948  for (Int_t cent = 0; cent < fNcentBins; cent++)
949  {
950  //=====================================================================================
951  Int_t h1_nbin = 300; Double_t h1_binwidth = 1; Double_t h1_low = 0;
952  Double_t h1_high = h1_low + h1_binwidth * h1_nbin; // 1GeV/bin
953  histname = TString::Format("%s/histo1_%d", groupname.Data(), cent);
954  histtitle = TString::Format("%s;#it{p}_{T,jet} (GeV/#it{c}) (accepted);Jets", histname.Data()); // Pt distro of jets
955  fHistManager.CreateTH1(histname, histtitle, h1_nbin, h1_low, h1_high);
956 
957  //=====================================================================================
958  Int_t h2_nbin = 200; Double_t h2_binwidth = 0.01; Double_t h2_low = -1;
959  Double_t h2_high = h2_low + h2_binwidth * h2_nbin;
960  histname = TString::Format("%s/histo2_%d", groupname.Data(), cent);
961  histtitle = TString::Format("%s;#it{#eta}_{jet};Jets", histname.Data()); // Eta distro of jets
962  fHistManager.CreateTH1(histname, histtitle, h2_nbin, h2_low, h2_high); // 1 unit of rapidity / 100 bin
963 
964  //=====================================================================================
965  Int_t h3_nbin = 126; Double_t h3_binwidth = 0.05; Double_t h3_low = 0.;
966  Double_t h3_high = h3_low + h3_binwidth * h3_nbin;
967  histname = TString::Format("%s/histo3_%d", groupname.Data(), cent);
968  histtitle = TString::Format("%s;#it{#phi}_{jet};Jets", histname.Data()); // Phi distro of jets
969  fHistManager.CreateTH1(histname, histtitle, h3_nbin, h3_low, h3_high);
970 
971  //=====================================================================================
972  Int_t h4_nbin = 100; Double_t h4_binwidth = 1; Double_t h4_low = 0;
973  Double_t h4_high = h4_low + h4_binwidth * h4_nbin; // 1 unit of multiplicity /bin
974  histname = TString::Format("%s/histo4_%d", groupname.Data(), cent);
975  histtitle = TString::Format("%s;N_{tracks}(jet);Jets", histname.Data()); // Multiplicity of all jets; chg tracks
976  fHistManager.CreateTH1(histname, histtitle, h4_nbin, h4_low, h4_high);
977 
978  histname = TString::Format("%s/histo4c_%d", groupname.Data(), cent);
979  histtitle = TString::Format("%s;N_{tracks}(jet);Jets", histname.Data());
980  fHistManager.CreateTH1(histname, histtitle, h4_nbin, h4_low, h4_high); // Multiplicity of all jets; all tracks
981 
982  histname = TString::Format("%s/histo6_%d", groupname.Data(), cent);
983  histtitle = TString::Format("%s;N_{tracks}(jet1);Jets", histname.Data()); // Multiplicity of jet1; chg tracks
984  fHistManager.CreateTH1(histname, histtitle, h4_nbin, h4_low, h4_high);
985 
986  histname = TString::Format("%s/histo6c_%d", groupname.Data(), cent);
987  histtitle = TString::Format("%s;N_{tracks}(jet1);Jets", histname.Data()); // Multiplicity of jet1; all tracks
988  fHistManager.CreateTH1(histname, histtitle, h4_nbin, h4_low, h4_high);
989  //#####################################
990 
991  //=====================================================================================
992  Int_t h5_nbin = 100; Double_t h5_binwidth = 1; Double_t h5_low = 0;
993  Double_t h5_high = h5_low + h5_binwidth * h5_nbin;
994  histname = TString::Format("%s/histo5_%d", groupname.Data(), cent);
995  histtitle = TString::Format("%s;N_{jets};Events", histname.Data()); // Distribution of jets in events
996  fHistManager.CreateTH1(histname, histtitle, h5_nbin, h5_low, h5_high);
997 
998  //=====================================================================================
999  Int_t h7_xnbin = 300; Double_t h7_xbinwidth = 1; Double_t h7_xlow = 0;
1000  Double_t h7_xhigh = h7_xlow + h7_xbinwidth * h7_xnbin;
1001  Int_t h7_ynbin = 100; Double_t h7_ybinwidth = 1; Double_t h7_ylow = 0;
1002  Double_t h7_yhigh = h7_ylow + h7_ybinwidth * h7_ynbin;
1003 
1004  histname = TString::Format("%s/histo7_%d", groupname.Data(), cent);
1005  histtitle = TString::Format("%s;#it{p}_{T,jet1} (GeV/c);N_{tracks}(jet1)", histname.Data()); // N vs pt jet1
1006  fHistManager.CreateTH2(histname, histtitle, h7_xnbin, h7_xlow, h7_xhigh, h7_ynbin, h7_ylow, h7_yhigh);
1007 
1008  histname = TString::Format("%s/histo7all_%d", groupname.Data(), cent);
1009  histtitle = TString::Format("%s;#it{p}_{T,jet} (GeV/c);N_{tracks}(jets)", histname.Data()); // N vs pt all jets
1010  fHistManager.CreateTH2(histname, histtitle, h7_xnbin, h7_xlow, h7_xhigh, h7_ynbin, h7_ylow, h7_yhigh);
1011  //#####################################
1012 
1013  //=====================================================================================
1014  Int_t h8_nbin = 101; Double_t h8_binwidth = 0.01; Double_t h8_low = 0;
1015  Double_t h8_high = h8_low + h8_binwidth * h8_nbin;
1016 
1017  // Standard implementation of Z
1018  histname = TString::Format("%s/histo8_%d", groupname.Data(), cent);
1019  histtitle = TString::Format("%s - jet1;z;F(Z) = 1/N_{jet1} dN/dz", histname.Data()); // scalar z ; jet1
1020  fHistManager.CreateTH1(histname, histtitle, h8_nbin, h8_low, h8_high);
1021 
1022  histname = TString::Format("%s/histo8_all_%d", groupname.Data(), cent);
1023  histtitle = TString::Format("%s - all jets;z;F(Z) = 1/N_{jets} dN/dz", histname.Data()); // scalar z ; all jets
1024  fHistManager.CreateTH1(histname, histtitle, h8_nbin, h8_low, h8_high);
1025 
1026  //########################################################
1027  // P_tot implementation of Z
1028  histname = TString::Format("%s/histo8_p_%d", groupname.Data(), cent);
1029  histtitle = TString::Format("%s - jet1 P_tot;z = p_{track}/p_{jet1};F(Z) = 1/N_{jet1} dN/dz", histname.Data());
1030  fHistManager.CreateTH1(histname, histtitle, h8_nbin, h8_low, h8_high);
1031 
1032  histname = TString::Format("%s/histo8_all_p_%d", groupname.Data(), cent);
1033  histtitle = TString::Format("%s - all jets P_tot;z = p_{track}/p_{jet};F(Z) = 1/N_{jets} dN/dz", histname.Data());
1034  fHistManager.CreateTH1(histname, histtitle, h8_nbin, h8_low, h8_high);
1035 
1036  //########################################################
1037  // Pt implementation of Z
1038  histname = TString::Format("%s/histo8_pt_%d", groupname.Data(), cent);
1039  histtitle = TString::Format("%s - jet1 Pt;z = p_{T,track}/p_{T,jet1};F(Z) = 1/N_{jet1} dN/dz", histname.Data());
1040  fHistManager.CreateTH1(histname, histtitle, h8_nbin, h8_low, h8_high);
1041 
1042  histname = TString::Format("%s/histo8_all_pt_%d", groupname.Data(), cent);
1043  histtitle = TString::Format("%s - all jets Pt;z = p_{T,track}/p_{T,jet1};F(Z) = 1/N_{jets} dN/dz", histname.Data());
1044  fHistManager.CreateTH1(histname, histtitle, h8_nbin, h8_low, h8_high);
1045  //########################################################
1046 
1047  //=====================================================================================
1048  Int_t h8xi_nbin = 140; Double_t h8xi_binwidth = 0.05; Double_t h8xi_low = 0;
1049  Double_t h8xi_high = h8xi_low + h8xi_binwidth * h8xi_nbin;
1050  histname = TString::Format("%s/histo8xi_%d", groupname.Data(), cent);
1051  histtitle = TString::Format("%s - jet1;#xi = ln(1/z);D(#xi) = 1/N_{jet1} dN/d#xi", histname.Data());
1052  fHistManager.CreateTH1(histname, histtitle, h8xi_nbin, h8xi_low, h8xi_high);
1053 
1054  histname = TString::Format("%s/histo8xi_all_%d", groupname.Data(), cent);
1055  histtitle = TString::Format("%s - all jets;#xi = ln(1/z);D(#xi) = 1/N_{jets} dN/d#xi", histname.Data());
1056  fHistManager.CreateTH1(histname, histtitle, h8xi_nbin, h8xi_low, h8xi_high);
1057 
1058  //########################################################
1059  histname = TString::Format("%s/histo8xi_p_%d", groupname.Data(), cent);
1060  histtitle = TString::Format("%s - jet1 P_tot;#xi = ln(1/z);D(#xi) = 1/N_{jet1} dN/d#xi", histname.Data());
1061  fHistManager.CreateTH1(histname, histtitle, h8xi_nbin, h8xi_low, h8xi_high);
1062 
1063  histname = TString::Format("%s/histo8xi_all_p_%d", groupname.Data(), cent);
1064  histtitle = TString::Format("%s - all jets P_tot;#xi = ln(1/z);D(#xi) = 1/N_{jets} dN/d#xi", histname.Data());
1065  fHistManager.CreateTH1(histname, histtitle, h8xi_nbin, h8xi_low, h8xi_high);
1066 
1067  //########################################################
1068  histname = TString::Format("%s/histo8xi_pt_%d", groupname.Data(), cent);
1069  histtitle = TString::Format("%s - jet1 Pt;#xi = ln(1/z);D(#xi) = 1/N_{jet1} dN/d#xi", histname.Data());
1070  fHistManager.CreateTH1(histname, histtitle, h8xi_nbin, h8xi_low, h8xi_high);
1071 
1072  histname = TString::Format("%s/histo8xi_all_pt_%d", groupname.Data(), cent);
1073  histtitle = TString::Format("%s - all jets Pt;#xi = ln(1/z);D(#xi) = 1/N_{jets} dN/d#xi", histname.Data());
1074  fHistManager.CreateTH1(histname, histtitle, h8xi_nbin, h8xi_low, h8xi_high);
1075  //########################################################
1076 
1077  //=====================================================================================
1078  Int_t h15_xnbin = 60; Double_t h15_xbinwidth = 0.01; Double_t h15_xlow = 0.;
1079  Double_t h15_xhigh = h15_xlow + h15_xbinwidth * h15_xnbin;
1080  Int_t h15_ynbin = 400; Double_t h15_ybinwidth = 1.; Double_t h15_ylow = 0.;
1081  Double_t h15_yhigh = h15_ylow + h15_ybinwidth * h15_ynbin;
1082 
1083  //########################################################
1084  histname = TString::Format("%s/histo15_%d", groupname.Data(), cent);
1085  histtitle = TString::Format("%s - jet1;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
1086  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
1087 
1088  //########################################################
1089  histname = TString::Format("%s/histo15_n70_%d", groupname.Data(), cent);
1090  histtitle = TString::Format("%s - jet1;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
1091  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
1092 
1093  histname = TString::Format("%s/histo15_n75_%d", groupname.Data(), cent);
1094  histtitle = TString::Format("%s - jet1;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
1095  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
1096 
1097  histname = TString::Format("%s/histo15_n80_%d", groupname.Data(), cent);
1098  histtitle = TString::Format("%s - jet1;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
1099  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
1100 
1101  histname = TString::Format("%s/histo15_n85_%d", groupname.Data(), cent);
1102  histtitle = TString::Format("%s - jet1;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
1103  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
1104 
1105  histname = TString::Format("%s/histo15_n90_%d", groupname.Data(), cent);
1106  histtitle = TString::Format("%s - jet1;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
1107  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
1108 
1109  //########################################################
1110  histname = TString::Format("%s/histo15_pt70_%d", groupname.Data(), cent);
1111  histtitle = TString::Format("%s - jet1;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
1112  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
1113 
1114  histname = TString::Format("%s/histo15_pt75_%d", groupname.Data(), cent);
1115  histtitle = TString::Format("%s - jet1;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
1116  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
1117 
1118  histname = TString::Format("%s/histo15_pt80_%d", groupname.Data(), cent);
1119  histtitle = TString::Format("%s - jet1;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
1120  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
1121 
1122  histname = TString::Format("%s/histo15_pt85_%d", groupname.Data(), cent);
1123  histtitle = TString::Format("%s - jet1;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
1124  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
1125 
1126  histname = TString::Format("%s/histo15_pt90_%d", groupname.Data(), cent);
1127  histtitle = TString::Format("%s - jet1;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
1128  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
1129  //########################################################
1130 
1131  //########################################################
1132  histname = TString::Format("%s/histo15all_%d", groupname.Data(), cent);
1133  histtitle = TString::Format("%s - all jets;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
1134  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
1135 
1136  //########################################################
1137  histname = TString::Format("%s/histo15all_n70_%d", groupname.Data(), cent);
1138  histtitle = TString::Format("%s - all jets;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
1139  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
1140 
1141  histname = TString::Format("%s/histo15all_n75_%d", groupname.Data(), cent);
1142  histtitle = TString::Format("%s - all jets;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
1143  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
1144 
1145  histname = TString::Format("%s/histo15all_n80_%d", groupname.Data(), cent);
1146  histtitle = TString::Format("%s - all jets;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
1147  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
1148 
1149  histname = TString::Format("%s/histo15all_n85_%d", groupname.Data(), cent);
1150  histtitle = TString::Format("%s - all jets;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
1151  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
1152 
1153  histname = TString::Format("%s/histo15all_n90_%d", groupname.Data(), cent);
1154  histtitle = TString::Format("%s - all jets;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
1155  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
1156 
1157  //########################################################
1158  histname = TString::Format("%s/histo15all_pt70_%d", groupname.Data(), cent);
1159  histtitle = TString::Format("%s - all jets;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
1160  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
1161 
1162  histname = TString::Format("%s/histo15all_pt75_%d", groupname.Data(), cent);
1163  histtitle = TString::Format("%s - all jets;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
1164  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
1165 
1166  histname = TString::Format("%s/histo15all_pt80_%d", groupname.Data(), cent);
1167  histtitle = TString::Format("%s - all jets;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
1168  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
1169 
1170  histname = TString::Format("%s/histo15all_pt85_%d", groupname.Data(), cent);
1171  histtitle = TString::Format("%s - all jets;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
1172  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
1173 
1174  histname = TString::Format("%s/histo15all_pt90_%d", groupname.Data(), cent);
1175  histtitle = TString::Format("%s - all jets;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
1176  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
1177  //########################################################
1178 
1179  //=====================================================================================
1180  Int_t h20_nbin = 60; Double_t h20_binwidth = 0.01; Double_t h20_low = 0.;
1181  Double_t h20_high = h20_low + h20_binwidth * h20_nbin;
1182 
1183  //########################################################
1184  histname = TString::Format("%s/histo20_%d", groupname.Data(), cent);
1185  histtitle = TString::Format("%s - jet1;R_{tracks};dN/dR", histname.Data());
1186  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
1187 
1188  //########################################################
1189  histname = TString::Format("%s/histo20_n70_%d", groupname.Data(), cent);
1190  histtitle = TString::Format("%s - jet1;R_{tracks};dN/dR", histname.Data());
1191  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
1192 
1193  histname = TString::Format("%s/histo20_n75_%d", groupname.Data(), cent);
1194  histtitle = TString::Format("%s - jet1;R_{tracks};dN/dR", histname.Data());
1195  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
1196 
1197  histname = TString::Format("%s/histo20_n80_%d", groupname.Data(), cent);
1198  histtitle = TString::Format("%s - jet1;R_{tracks};dN/dR", histname.Data());
1199  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
1200 
1201  histname = TString::Format("%s/histo20_n85_%d", groupname.Data(), cent);
1202  histtitle = TString::Format("%s - jet1;R_{tracks};dN/dR", histname.Data());
1203  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
1204 
1205  histname = TString::Format("%s/histo20_n90_%d", groupname.Data(), cent);
1206  histtitle = TString::Format("%s - jet1;R_{tracks};dN/dR", histname.Data());
1207  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
1208  //########################################################
1209 
1210  //########################################################
1211  histname = TString::Format("%s/histo20_pt70_%d", groupname.Data(), cent);
1212  histtitle = TString::Format("%s - jet1;R_{tracks};dN/dR", histname.Data());
1213  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
1214 
1215  histname = TString::Format("%s/histo20_pt75_%d", groupname.Data(), cent);
1216  histtitle = TString::Format("%s - jet1;R_{tracks};dN/dR", histname.Data());
1217  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
1218 
1219  histname = TString::Format("%s/histo20_pt80_%d", groupname.Data(), cent);
1220  histtitle = TString::Format("%s - jet1;R_{tracks};dN/dR", histname.Data());
1221  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
1222 
1223  histname = TString::Format("%s/histo20_pt85_%d", groupname.Data(), cent);
1224  histtitle = TString::Format("%s - jet1;R_{tracks};dN/dR", histname.Data());
1225  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
1226 
1227  histname = TString::Format("%s/histo20_pt90_%d", groupname.Data(), cent);
1228  histtitle = TString::Format("%s - jet1;R_{tracks};dN/dR", histname.Data());
1229  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
1230  //########################################################
1231 
1232  //########################################################
1233  histname = TString::Format("%s/histo20all_%d", groupname.Data(), cent);
1234  histtitle = TString::Format("%s - all jets;R_{tracks};dN/dR", histname.Data());
1235  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
1236 
1237  //########################################################
1238  histname = TString::Format("%s/histo20all_n70_%d", groupname.Data(), cent);
1239  histtitle = TString::Format("%s - all jets;R_{tracks};dN/dR", histname.Data());
1240  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
1241 
1242  histname = TString::Format("%s/histo20all_n75_%d", groupname.Data(), cent);
1243  histtitle = TString::Format("%s - all jets;R_{tracks};dN/dR", histname.Data());
1244  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
1245 
1246  histname = TString::Format("%s/histo20all_n80_%d", groupname.Data(), cent);
1247  histtitle = TString::Format("%s - all jets;R_{tracks};dN/dR", histname.Data());
1248  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
1249 
1250  histname = TString::Format("%s/histo20all_n85_%d", groupname.Data(), cent);
1251  histtitle = TString::Format("%s - all jets;R_{tracks};dN/dR", histname.Data());
1252  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
1253 
1254  histname = TString::Format("%s/histo20all_n90_%d", groupname.Data(), cent);
1255  histtitle = TString::Format("%s - all jets;R_{tracks};dN/dR", histname.Data());
1256  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
1257  //########################################################
1258 
1259  //########################################################
1260  histname = TString::Format("%s/histo20all_pt70_%d", groupname.Data(), cent);
1261  histtitle = TString::Format("%s - all jets;R_{tracks};dN/dR", histname.Data());
1262  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
1263 
1264  histname = TString::Format("%s/histo20all_pt75_%d", groupname.Data(), cent);
1265  histtitle = TString::Format("%s - all jets;R_{tracks};dN/dR", histname.Data());
1266  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
1267 
1268  histname = TString::Format("%s/histo20all_pt80_%d", groupname.Data(), cent);
1269  histtitle = TString::Format("%s - all jets;R_{tracks};dN/dR", histname.Data());
1270  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
1271 
1272  histname = TString::Format("%s/histo20all_pt85_%d", groupname.Data(), cent);
1273  histtitle = TString::Format("%s - all jets;R_{tracks};dN/dR", histname.Data());
1274  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
1275 
1276  histname = TString::Format("%s/histo20all_pt90_%d", groupname.Data(), cent);
1277  histtitle = TString::Format("%s - all jets;R_{tracks};dN/dR", histname.Data());
1278  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
1279  //########################################################
1280 
1281  //=====================================================================================
1282  // Distribution of girth (radial girth) g = sum_jet_parts ( r_i * ( pt_i/pt_jet ) )
1283  Int_t hg_nbin = 50; Double_t hg_binwidth = 0.01; Double_t hg_low = 0.;
1284  Double_t hg_high = hg_low + hg_binwidth * hg_nbin;
1285 
1286  //########################################################
1287  histname = TString::Format("%s/histo_g_%d", groupname.Data(), cent);
1288  histtitle = TString::Format("%s - all jets;g;1/N_{jets} dN/dg", histname.Data());
1289  fHistManager.CreateTH1(histname, histtitle, hg_nbin, hg_low, hg_high);
1290  //########################################################
1291 
1292  //########################################################
1293  histname = TString::Format("%s/histo_g_n70_%d", groupname.Data(), cent);
1294  histtitle = TString::Format("%s - all jets;g;1/N_{jets} dN/dg", histname.Data());
1295  fHistManager.CreateTH1(histname, histtitle, hg_nbin, hg_low, hg_high);
1296 
1297  histname = TString::Format("%s/histo_g_n75_%d", groupname.Data(), cent);
1298  histtitle = TString::Format("%s - all jets;g;1/N_{jets} dN/dg", histname.Data());
1299  fHistManager.CreateTH1(histname, histtitle, hg_nbin, hg_low, hg_high);
1300 
1301  histname = TString::Format("%s/histo_g_n80_%d", groupname.Data(), cent);
1302  histtitle = TString::Format("%s - all jets;g;1/N_{jets} dN/dg", histname.Data());
1303  fHistManager.CreateTH1(histname, histtitle, hg_nbin, hg_low, hg_high);
1304 
1305  histname = TString::Format("%s/histo_g_n85_%d", groupname.Data(), cent);
1306  histtitle = TString::Format("%s - all jets;g;1/N_{jets} dN/dg", histname.Data());
1307  fHistManager.CreateTH1(histname, histtitle, hg_nbin, hg_low, hg_high);
1308 
1309  histname = TString::Format("%s/histo_g_n90_%d", groupname.Data(), cent);
1310  histtitle = TString::Format("%s - all jets;g;1/N_{jets} dN/dg", histname.Data());
1311  fHistManager.CreateTH1(histname, histtitle, hg_nbin, hg_low, hg_high);
1312  //########################################################
1313 
1314  //########################################################
1315  histname = TString::Format("%s/histo_g_pt70_%d", groupname.Data(), cent);
1316  histtitle = TString::Format("%s - all jets;g;1/N_{jets} dN/dg", histname.Data());
1317  fHistManager.CreateTH1(histname, histtitle, hg_nbin, hg_low, hg_high);
1318 
1319  histname = TString::Format("%s/histo_g_pt75_%d", groupname.Data(), cent);
1320  histtitle = TString::Format("%s - all jets;g;1/N_{jets} dN/dg", histname.Data());
1321  fHistManager.CreateTH1(histname, histtitle, hg_nbin, hg_low, hg_high);
1322 
1323  histname = TString::Format("%s/histo_g_pt80_%d", groupname.Data(), cent);
1324  histtitle = TString::Format("%s - all jets;g;1/N_{jets} dN/dg", histname.Data());
1325  fHistManager.CreateTH1(histname, histtitle, hg_nbin, hg_low, hg_high);
1326 
1327  histname = TString::Format("%s/histo_g_pt85_%d", groupname.Data(), cent);
1328  histtitle = TString::Format("%s - all jets;g;1/N_{jets} dN/dg", histname.Data());
1329  fHistManager.CreateTH1(histname, histtitle, hg_nbin, hg_low, hg_high);
1330 
1331  histname = TString::Format("%s/histo_g_pt90_%d", groupname.Data(), cent);
1332  histtitle = TString::Format("%s - all jets;g;1/N_{jets} dN/dg", histname.Data());
1333  fHistManager.CreateTH1(histname, histtitle, hg_nbin, hg_low, hg_high);
1334  //########################################################
1335 
1336  //=====================================================================================
1337  // Distribution of dispersion d pt_D = sqrt ( sum (pt_i^2) )/sum (pt_i)
1338  Int_t hptd_nbin = 40; Double_t hptd_binwidth = 0.05; Double_t hptd_low = 0.;
1339  Double_t hptd_high = hptd_low + hptd_binwidth * hptd_nbin;
1340 
1341  //########################################################
1342  histname = TString::Format("%s/histo_ptd_%d", groupname.Data(), cent);
1343  histtitle = TString::Format("%s - all jets;ptd;1/N_{jets} dN/dp_{T}D", histname.Data());
1344  fHistManager.CreateTH1(histname, histtitle, hptd_nbin, hptd_low, hptd_high);
1345  //########################################################
1346 
1347  //########################################################
1348  histname = TString::Format("%s/histo_ptd_n70_%d", groupname.Data(), cent);
1349  histtitle = TString::Format("%s - all jets;ptd;1/N_{jets} dN/dp_{T}D", histname.Data());
1350  fHistManager.CreateTH1(histname, histtitle, hptd_nbin, hptd_low, hptd_high);
1351 
1352  histname = TString::Format("%s/histo_ptd_n75_%d", groupname.Data(), cent);
1353  histtitle = TString::Format("%s - all jets;ptd;1/N_{jets} dN/dp_{T}D", histname.Data());
1354  fHistManager.CreateTH1(histname, histtitle, hptd_nbin, hptd_low, hptd_high);
1355 
1356  histname = TString::Format("%s/histo_ptd_n80_%d", groupname.Data(), cent);
1357  histtitle = TString::Format("%s - all jets;ptd;1/N_{jets} dN/dp_{T}D", histname.Data());
1358  fHistManager.CreateTH1(histname, histtitle, hptd_nbin, hptd_low, hptd_high);
1359 
1360  histname = TString::Format("%s/histo_ptd_n85_%d", groupname.Data(), cent);
1361  histtitle = TString::Format("%s - all jets;ptd;1/N_{jets} dN/dp_{T}D", histname.Data());
1362  fHistManager.CreateTH1(histname, histtitle, hptd_nbin, hptd_low, hptd_high);
1363 
1364  histname = TString::Format("%s/histo_ptd_n90_%d", groupname.Data(), cent);
1365  histtitle = TString::Format("%s - all jets;ptd;1/N_{jets} dN/dp_{T}D", histname.Data());
1366  fHistManager.CreateTH1(histname, histtitle, hptd_nbin, hptd_low, hptd_high);
1367  //########################################################
1368 
1369  //########################################################
1370  histname = TString::Format("%s/histo_ptd_pt70_%d", groupname.Data(), cent);
1371  histtitle = TString::Format("%s - all jets;ptd;1/N_{jets} dN/dp_{T}D", histname.Data());
1372  fHistManager.CreateTH1(histname, histtitle, hptd_nbin, hptd_low, hptd_high);
1373 
1374  histname = TString::Format("%s/histo_ptd_pt75_%d", groupname.Data(), cent);
1375  histtitle = TString::Format("%s - all jets;ptd;1/N_{jets} dN/dp_{T}D", histname.Data());
1376  fHistManager.CreateTH1(histname, histtitle, hptd_nbin, hptd_low, hptd_high);
1377 
1378  histname = TString::Format("%s/histo_ptd_pt80_%d", groupname.Data(), cent);
1379  histtitle = TString::Format("%s - all jets;ptd;1/N_{jets} dN/dp_{T}D", histname.Data());
1380  fHistManager.CreateTH1(histname, histtitle, hptd_nbin, hptd_low, hptd_high);
1381 
1382  histname = TString::Format("%s/histo_ptd_pt85_%d", groupname.Data(), cent);
1383  histtitle = TString::Format("%s - all jets;ptd;1/N_{jets} dN/dp_{T}D", histname.Data());
1384  fHistManager.CreateTH1(histname, histtitle, hptd_nbin, hptd_low, hptd_high);
1385 
1386  histname = TString::Format("%s/histo_ptd_pt90_%d", groupname.Data(), cent);
1387  histtitle = TString::Format("%s - all jets;ptd;1/N_{jets} dN/dp_{T}D", histname.Data());
1388  fHistManager.CreateTH1(histname, histtitle, hptd_nbin, hptd_low, hptd_high);
1389  //########################################################
1390 
1391  //=====================================================================================
1392  Int_t h_Rjt_xnbin = 60; Double_t h_Rjt_xbinwidth = 0.01; Double_t h_Rjt_xlow = 0.;
1393  Double_t h_Rjt_xhigh = h_Rjt_xlow + h_Rjt_xbinwidth * h_Rjt_xnbin;
1394  Int_t h_Rjt_ynbin = 500; Double_t h_Rjt_ybinwidth = 0.01; Double_t h_Rjt_ylow = 0.;
1395  Double_t h_Rjt_yhigh = h_Rjt_ylow + h_Rjt_ybinwidth * h_Rjt_ynbin;
1396 
1397  //########################################################
1398  histname = TString::Format("%s/histo_Rjt_%d", groupname.Data(), cent);
1399  histtitle = TString::Format("%s ;dR;j_{T} (GeV/c);", histname.Data()); // j_T track vs dR
1400  fHistManager.CreateTH2(histname, histtitle, h_Rjt_xnbin, h_Rjt_xlow, h_Rjt_xhigh, h_Rjt_ynbin, h_Rjt_ylow, h_Rjt_yhigh);
1401 
1402  //########################################################
1403  histname = TString::Format("%s/histo_Rjt_n70_%d", groupname.Data(), cent);
1404  histtitle = TString::Format("%s ;dR;j_{T} (GeV/c);", histname.Data()); // j_T track vs dR
1405  fHistManager.CreateTH2(histname, histtitle, h_Rjt_xnbin, h_Rjt_xlow, h_Rjt_xhigh, h_Rjt_ynbin, h_Rjt_ylow, h_Rjt_yhigh);
1406 
1407  histname = TString::Format("%s/histo_Rjt_n75_%d", groupname.Data(), cent);
1408  histtitle = TString::Format("%s ;dR;j_{T} (GeV/c);", histname.Data()); // j_T track vs dR
1409  fHistManager.CreateTH2(histname, histtitle, h_Rjt_xnbin, h_Rjt_xlow, h_Rjt_xhigh, h_Rjt_ynbin, h_Rjt_ylow, h_Rjt_yhigh);
1410 
1411  histname = TString::Format("%s/histo_Rjt_n80_%d", groupname.Data(), cent);
1412  histtitle = TString::Format("%s ;dR;j_{T} (GeV/c);", histname.Data()); // j_T track vs dR
1413  fHistManager.CreateTH2(histname, histtitle, h_Rjt_xnbin, h_Rjt_xlow, h_Rjt_xhigh, h_Rjt_ynbin, h_Rjt_ylow, h_Rjt_yhigh);
1414 
1415  histname = TString::Format("%s/histo_Rjt_n85_%d", groupname.Data(), cent);
1416  histtitle = TString::Format("%s ;dR;j_{T} (GeV/c);", histname.Data()); // j_T track vs dR
1417  fHistManager.CreateTH2(histname, histtitle, h_Rjt_xnbin, h_Rjt_xlow, h_Rjt_xhigh, h_Rjt_ynbin, h_Rjt_ylow, h_Rjt_yhigh);
1418 
1419  histname = TString::Format("%s/histo_Rjt_n90_%d", groupname.Data(), cent);
1420  histtitle = TString::Format("%s ;dR;j_{T} (GeV/c);", histname.Data()); // j_T track vs dR
1421  fHistManager.CreateTH2(histname, histtitle, h_Rjt_xnbin, h_Rjt_xlow, h_Rjt_xhigh, h_Rjt_ynbin, h_Rjt_ylow, h_Rjt_yhigh);
1422 
1423  //########################################################
1424  histname = TString::Format("%s/histo_Rjt_pt70_%d", groupname.Data(), cent);
1425  histtitle = TString::Format("%s ;dR;j_{T} (GeV/c);", histname.Data()); // j_T track vs dR
1426  fHistManager.CreateTH2(histname, histtitle, h_Rjt_xnbin, h_Rjt_xlow, h_Rjt_xhigh, h_Rjt_ynbin, h_Rjt_ylow, h_Rjt_yhigh);
1427 
1428  histname = TString::Format("%s/histo_Rjt_pt75_%d", groupname.Data(), cent);
1429  histtitle = TString::Format("%s ;dR;j_{T} (GeV/c);", histname.Data()); // j_T track vs dR
1430  fHistManager.CreateTH2(histname, histtitle, h_Rjt_xnbin, h_Rjt_xlow, h_Rjt_xhigh, h_Rjt_ynbin, h_Rjt_ylow, h_Rjt_yhigh);
1431 
1432  histname = TString::Format("%s/histo_Rjt_pt80_%d", groupname.Data(), cent);
1433  histtitle = TString::Format("%s ;dR;j_{T} (GeV/c);", histname.Data()); // j_T track vs dR
1434  fHistManager.CreateTH2(histname, histtitle, h_Rjt_xnbin, h_Rjt_xlow, h_Rjt_xhigh, h_Rjt_ynbin, h_Rjt_ylow, h_Rjt_yhigh);
1435 
1436  histname = TString::Format("%s/histo_Rjt_pt85_%d", groupname.Data(), cent);
1437  histtitle = TString::Format("%s ;dR;j_{T} (GeV/c);", histname.Data()); // j_T track vs dR
1438  fHistManager.CreateTH2(histname, histtitle, h_Rjt_xnbin, h_Rjt_xlow, h_Rjt_xhigh, h_Rjt_ynbin, h_Rjt_ylow, h_Rjt_yhigh);
1439 
1440  histname = TString::Format("%s/histo_Rjt_pt90_%d", groupname.Data(), cent);
1441  histtitle = TString::Format("%s ;dR;j_{T} (GeV/c);", histname.Data()); // j_T track vs dR
1442  fHistManager.CreateTH2(histname, histtitle, h_Rjt_xnbin, h_Rjt_xlow, h_Rjt_xhigh, h_Rjt_ynbin, h_Rjt_ylow, h_Rjt_yhigh);
1443  //########################################################
1444 
1445  //=====================================================================================
1446  Int_t h_jt_xnbin = 500; Double_t h_jt_xbinwidth = 0.01; Double_t h_jt_xlow = 0.;
1447  Double_t h_jt_xhigh = h_jt_xlow + h_jt_xbinwidth * h_jt_xnbin;
1448 
1449  //########################################################
1450  histname = TString::Format("%s/histo_jt_%d", groupname.Data(), cent);
1451  histtitle = TString::Format("%s ;j_{T} (GeV/c);1/N_{jets} dN/dj_{T};", histname.Data()); // j_T track vs dR
1452  fHistManager.CreateTH1(histname, histtitle, h_jt_xnbin, h_jt_xlow, h_jt_xhigh);
1453 
1454  //########################################################
1455  histname = TString::Format("%s/histo_jt_n70_%d", groupname.Data(), cent);
1456  histtitle = TString::Format("%s ;j_{T} (GeV/c);1/N_{jets} dN/dj_{T};", histname.Data()); // j_T track vs dR
1457  fHistManager.CreateTH1(histname, histtitle, h_jt_xnbin, h_jt_xlow, h_jt_xhigh);
1458 
1459  histname = TString::Format("%s/histo_jt_n75_%d", groupname.Data(), cent);
1460  histtitle = TString::Format("%s ;j_{T} (GeV/c);1/N_{jets} dN/dj_{T};", histname.Data()); // j_T track vs dR
1461  fHistManager.CreateTH1(histname, histtitle, h_jt_xnbin, h_jt_xlow, h_jt_xhigh);
1462 
1463  histname = TString::Format("%s/histo_jt_n80_%d", groupname.Data(), cent);
1464  histtitle = TString::Format("%s ;j_{T} (GeV/c);1/N_{jets} dN/dj_{T};", histname.Data()); // j_T track vs dR
1465  fHistManager.CreateTH1(histname, histtitle, h_jt_xnbin, h_jt_xlow, h_jt_xhigh);
1466 
1467  histname = TString::Format("%s/histo_jt_n85_%d", groupname.Data(), cent);
1468  histtitle = TString::Format("%s ;j_{T} (GeV/c);1/N_{jets} dN/dj_{T};", histname.Data()); // j_T track vs dR
1469  fHistManager.CreateTH1(histname, histtitle, h_jt_xnbin, h_jt_xlow, h_jt_xhigh);
1470 
1471  histname = TString::Format("%s/histo_jt_n90_%d", groupname.Data(), cent);
1472  histtitle = TString::Format("%s ;j_{T} (GeV/c);1/N_{jets} dN/dj_{T};", histname.Data()); // j_T track vs dR
1473  fHistManager.CreateTH1(histname, histtitle, h_jt_xnbin, h_jt_xlow, h_jt_xhigh);
1474 
1475  //########################################################
1476  histname = TString::Format("%s/histo_jt_pt70_%d", groupname.Data(), cent);
1477  histtitle = TString::Format("%s ;j_{T} (GeV/c);1/N_{jets} dN/dj_{T};", histname.Data()); // j_T track vs dR
1478  fHistManager.CreateTH1(histname, histtitle, h_jt_xnbin, h_jt_xlow, h_jt_xhigh);
1479 
1480  histname = TString::Format("%s/histo_jt_pt75_%d", groupname.Data(), cent);
1481  histtitle = TString::Format("%s ;j_{T} (GeV/c);1/N_{jets} dN/dj_{T};", histname.Data()); // j_T track vs dR
1482  fHistManager.CreateTH1(histname, histtitle, h_jt_xnbin, h_jt_xlow, h_jt_xhigh);
1483 
1484  histname = TString::Format("%s/histo_jt_pt80_%d", groupname.Data(), cent);
1485  histtitle = TString::Format("%s ;j_{T} (GeV/c);1/N_{jets} dN/dj_{T};", histname.Data()); // j_T track vs dR
1486  fHistManager.CreateTH1(histname, histtitle, h_jt_xnbin, h_jt_xlow, h_jt_xhigh);
1487 
1488  histname = TString::Format("%s/histo_jt_pt85_%d", groupname.Data(), cent);
1489  histtitle = TString::Format("%s ;j_{T} (GeV/c);1/N_{jets} dN/dj_{T};", histname.Data()); // j_T track vs dR
1490  fHistManager.CreateTH1(histname, histtitle, h_jt_xnbin, h_jt_xlow, h_jt_xhigh);
1491 
1492  histname = TString::Format("%s/histo_jt_pt90_%d", groupname.Data(), cent);
1493  histtitle = TString::Format("%s ;j_{T} (GeV/c);1/N_{jets} dN/dj_{T};", histname.Data()); // j_T track vs dR
1494  fHistManager.CreateTH1(histname, histtitle, h_jt_xnbin, h_jt_xlow, h_jt_xhigh);
1495  //########################################################
1496 
1497  }
1498  //end of loop over fNcentBins
1499  }
1500  // end of loop over jet containers
1501 
1502  // =========== Switch on Sumw2 for all histos ===========
1503  TH1::SetDefaultSumw2(kTRUE);
1504  TH2::SetDefaultSumw2(kTRUE);
1505 
1506  // add all fHistManager content to fOutput
1507  TIter nexthist(fHistManager.GetListOfHistograms());
1508  TObject* obj = NULL;
1509  while ((obj = nexthist())) { fOutput->Add(obj); }
1510 
1511  PostData ( 1, fOutput ); // Post data for ALL output slots > 0 here.
1512  }
1513 
1519  {
1521  }
1522 
1527  {
1528  }
1529 
1530 //________________________________________________________________________
1532 {
1533  return fHistManager.FindObject(histName);
1534 }
1535 
1536 //########################################################################
1537 // Namespace AliAnalysisTaskEmcalJetCDF
1538 //########################################################################
1539 
1540 //__________________________________________________________________________________________________
1541 std::vector<Int_t> NS_AliAnalysisTaskEmcalJetCDF::SortTracksPt ( AliVEvent* event )
1542  {
1543  // Sorting by p_T (decreasing) event tracks
1544  Int_t entries = event->GetNumberOfTracks();
1545 
1546  // Create vector for Pt sorting
1547  std::vector<ptidx_pair> pair_list;
1548  pair_list.reserve ( entries );
1549 
1550  for ( Int_t i_entry = 0; i_entry < entries; i_entry++ )
1551  {
1552  AliVParticle* track = event->GetTrack ( i_entry );
1553  if ( !track ) { std::cout << Form ("Unable to find track %d in collection %s", i_entry, event->GetName()) << std::endl ; continue; }
1554 
1555  pair_list.push_back ( std::make_pair ( track->Pt(), i_entry ) );
1556  }
1557 
1558  std::stable_sort ( pair_list.begin(), pair_list.end(), sort_descend() );
1559 
1560  // return an vector of indexes of constituents (sorted descending by pt)
1561  std::vector<Int_t> index_sorted_list;
1562  index_sorted_list.reserve ( entries );
1563 
1564  // populating the return object with indexes of sorted tracks
1565  for (auto it : pair_list) { index_sorted_list.push_back(it.second); }
1566 
1567  return index_sorted_list;
1568  }
1569 
1570 //__________________________________________________________________________________________________
1572  {
1573  // Sorting by p_T (decreasing) event tracks
1574  Int_t entries = trackscont->GetNEntries();
1575 
1576  // Create vector for Pt sorting
1577  std::vector<ptidx_pair> pair_list;
1578  pair_list.reserve ( entries );
1579 
1580  UInt_t i_entry = 0;
1581  AliParticleContainer* partCont = 0;
1582  for(auto part : partCont->all()) {
1583  if (!part) {continue;}
1584  pair_list.push_back ( std::make_pair ( part->Pt(), i_entry++ ) );
1585  }
1586 
1587  std::stable_sort ( pair_list.begin(), pair_list.end(), sort_descend() );
1588 
1589  // return an vector of indexes of constituents (sorted descending by pt)
1590  std::vector<Int_t> index_sorted_list;
1591  index_sorted_list.reserve ( i_entry );
1592 
1593  // populating the return object with indexes of sorted tracks
1594  for (auto it : pair_list) { index_sorted_list.push_back(it.second); }
1595 
1596  return index_sorted_list;
1597  }
1598 
1599 
1606 AliAnalysisTaskEmcalJetCDF* NS_AliAnalysisTaskEmcalJetCDF::AddTaskEmcalJetCDF ( const char* ntracks, const char* nclusters, const char* ncells, const char* tag)
1607  {
1608  // Get the pointer to the existing analysis manager via the static access method.
1609  //==============================================================================
1610  AliAnalysisManager* mgr = AliAnalysisManager::GetAnalysisManager();
1611  if ( !mgr ) { ::Error ( "AddTaskEmcalJetCDF", "No analysis manager to connect to." ); return NULL; }
1612 
1613  // Check the analysis type using the event handlers connected to the analysis manager.
1614  //==============================================================================
1615  AliVEventHandler* handler = mgr->GetInputEventHandler();
1616  if (!handler) { ::Error ( "AddTaskEmcalJetCDF", "This task requires an input event handler" ); return NULL; }
1617 
1618  enum EDataType_t { kUnknown, kESD, kAOD }; EDataType_t dataType = kUnknown;
1619 
1620  if (handler->InheritsFrom("AliESDInputHandler")) { dataType = kESD; }
1621  else
1622  if (handler->InheritsFrom("AliAODInputHandler")) { dataType = kAOD; }
1623 
1624  //-------------------------------------------------------
1625  // Init the task and do settings
1626  //-------------------------------------------------------
1627  TString suffix ( tag );
1628  TString tracks ( ntracks );
1629  TString clusters ( nclusters );
1630  TString cells ( ncells );
1631 
1632  if ( tracks.EqualTo("usedefault") )
1633  {
1634  if ( dataType == kESD ) { tracks = "Tracks"; }
1635  else
1636  if ( dataType == kAOD ) { tracks = "tracks"; }
1637  else
1638  { tracks = ""; }
1639  }
1640 
1641  if ( clusters.EqualTo("usedefault") )
1642  {
1643  if ( dataType == kESD ) { clusters = "CaloClusters"; }
1644  else
1645  if ( dataType == kAOD ) { clusters = "caloClusters"; }
1646  else
1647  { clusters = ""; }
1648  }
1649 
1650  if ( cells.EqualTo("usedefault") )
1651  {
1652  if (dataType == kESD) { cells = "EMCALCells"; }
1653  else
1654  if (dataType == kAOD) { cells = "emcalCells"; }
1655  else
1656  { cells = ""; }
1657  }
1658 
1659  TString name("JetCDF");
1660  if (!tracks.IsNull()) { name += "_" + tracks; }
1661  if (!clusters.IsNull()) { name += "_" + clusters; }
1662  if (!cells.IsNull()) { name += "_" + cells; }
1663  if (!suffix.IsNull()) { name += "_" + suffix; }
1664 
1665  AliAnalysisTaskEmcalJetCDF* cdfTask = new AliAnalysisTaskEmcalJetCDF ( name.Data() );
1666  cdfTask->SetVzRange(-10,10);
1667  cdfTask->SetCaloCellsName(cells.Data());
1668 
1669  if ( tracks.EqualTo("mcparticles") ) {
1670  // AliMCParticleContainer* mcpartCont =
1671  cdfTask->AddMCParticleContainer ( tracks.Data() );
1672  }
1673  else
1674  if ( tracks.EqualTo("tracks") || tracks.EqualTo("Tracks") ) {
1675  // AliTrackContainer* trackCont =
1676  cdfTask->AddTrackContainer( tracks.Data() );
1677  }
1678  else
1679  if ( !tracks.IsNull())
1680  { cdfTask->AddParticleContainer(tracks.Data()); }
1681 
1682  cdfTask->AddClusterContainer(clusters.Data());
1683 
1684  //-------------------------------------------------------
1685  // Final settings, pass to manager and set the containers
1686  //-------------------------------------------------------
1687  mgr->AddTask ( cdfTask );
1688 
1689  // Create containers for input/output
1690  AliAnalysisDataContainer *cinput1 = mgr->GetCommonInputContainer();
1691 
1692  TString contname = name + "_histos";
1693  TString outfile (Form("%s", AliAnalysisManager::GetCommonFileName()));
1694  AliAnalysisDataContainer *coutput1 = mgr->CreateContainer ( contname.Data(), TList::Class(), AliAnalysisManager::kOutputContainer, outfile.Data() );
1695 
1696  mgr->ConnectInput ( cdfTask, 0, cinput1 );
1697  mgr->ConnectOutput ( cdfTask, 1, coutput1 );
1698 
1699  return cdfTask;
1700  }
1701 
1711  void NS_AliAnalysisTaskEmcalJetCDF::jetContSetParams ( AliJetContainer* jetCont, Float_t jetptmin, Float_t jetptmax, Float_t jetareacutperc, Int_t leadhadtype, Int_t nLeadJets, Float_t mintrackpt, Float_t maxtrackpt)
1712  {
1713  if (!jetCont) { return; }
1714  jetCont->SetJetPtCut ( jetptmin );
1715  jetCont->SetJetPtCutMax ( jetptmax );
1716  jetCont->SetPercAreaCut ( jetareacutperc );
1717  jetCont->SetLeadingHadronType ( leadhadtype ); // 0 = charged, 1 = neutral, 2 = both
1718  jetCont->SetNLeadingJets(nLeadJets);
1719  jetCont->SetMinTrackPt(mintrackpt);
1720  jetCont->SetMaxTrackPt(maxtrackpt);
1721  }
1722 
1735 TChain* NS_AliAnalysisTaskEmcalJetCDF::CreateChain ( const char* filelist, const char* cTreeNameArg, const char* friends, UInt_t iNumFiles, UInt_t iStartWithFile ) {
1736 TString sTreeNameArg (cTreeNameArg), treeName;
1737 
1738 TFileCollection filecoll ("anachain","File collection for analysis"); // easy manipulation of file collections
1739 Int_t iAddedFiles = filecoll.AddFromFile(filelist,iNumFiles,iStartWithFile);
1740 if ( iAddedFiles < 1 ) { std::cout << "NO Files added to collection !!!" << std::endl; return NULL; }
1741 
1742 // if cTreeNameArg is auto lets try to autodetect what type of tree we have;
1743 // the assuption is that all files are the same and the first one is reprezentative
1744 THashList* list = filecoll.GetList();
1745 if ( sTreeNameArg.EqualTo("auto") ) { // if tree name is not specified
1746  TRegexp tree_regex ("[aod,esd]Tree");
1747  TFileInfo* fileinfo = dynamic_cast<TFileInfo*>(list->At(0)); // get first fileinfo in list
1748  TFile file (fileinfo->GetFirstUrl()->GetFile()); // get the actual TFile
1749  if (file.IsZombie()) { cout << "Should not reach this message!! Error opening file" << endl; return NULL; }
1750 
1751  // lets parse the TFile
1752  TIter next(file.GetListOfKeys());
1753  TKey* key = NULL;
1754  while (( key = dynamic_cast<TKey*>(next()) )) {
1755  TString class_name = key->GetClassName();
1756  if ( ! class_name.EqualTo("TTree") ) { continue; } // searching for first TTree
1757 
1758  TString key_name = key->GetName();
1759  if ( key_name.Contains(tree_regex) ) { treeName = key_name; break;} // that is named either aodTree or esdTree
1760  }
1761  file.Close();
1762  }
1763 else
1764  { treeName = sTreeNameArg ; } // tree name is specified
1765 
1766 TChain* chain = new TChain (treeName.Data(),""); // lets create the chain
1767 if ( chain->AddFileInfoList(list) == 0 ) { return NULL; } // and add file collection (THashList is a Tlist that is a TCollection)
1768 
1769 // start treatment of friends
1770 TChain* chainFriend = NULL;
1771 TString sFriends (friends);
1772 if (!sFriends.IsNull()) {
1773  TString friends_treename, friends_filename;
1774  TObjArray* arr = sFriends.Tokenize("/");
1775  TObjString* strobj_file = dynamic_cast<TObjString*>(arr->At(0));
1776  if (strobj_file) { friends_filename = strobj_file->GetString(); }
1777  TObjString* strobj_tree = dynamic_cast<TObjString*>(arr->At(1));
1778  if (strobj_tree) { friends_treename = strobj_tree->GetString(); }
1779  delete arr;
1780 
1781  if (friends_treename.IsNull()) {
1782  if (treeName.EqualTo("esdTree")) {
1783  friends_treename = "esdFriendTree"; }
1784  else if (treeName.EqualTo("aodTree")) {
1785  friends_treename = "aodTree"; }
1786  else {
1787  cout << "friends argument specified but tree name neither specified nor auto-detected (unknown tree name to associate with known friend tree name)";
1788  return chain; // stop processing of friends, just return the chain created so far
1789  }
1790  }
1791 
1792  chainFriend = new TChain(friends_treename.Data());
1793  TString friendinfo_for_chain = "/" + friends_filename + "/" + friends_treename;
1794 
1795  TIter next(list);
1796  TFileInfo* fileinfo = NULL;
1797  while (( fileinfo = dynamic_cast<TFileInfo*>(next()) )) {
1798  TString dirname = gSystem->DirName(fileinfo->GetFirstUrl()->GetFile());
1799  TString friend_for_chain = dirname + friendinfo_for_chain;
1800  if (chainFriend) { chainFriend->Add(friend_for_chain.Data()); }
1801  }
1802 
1803  if (chainFriend) { chain->AddFriend(chainFriend); }
1804  }
1805 
1806 return chain;
1807 }
1808 
1809 
1810 
1811 // kate: indent-mode none; indent-width 2; replace-tabs on;
1812 
THistManager fHistManager
Histogram manager.
THashList * CreateHistoGroup(const char *groupname)
Create a new group of histograms within a parent group.
double Double_t
Definition: External.C:58
Double_t GetXi(const AliVParticle *trk) const
const AliParticleIterableContainer all() const
void SetLeadingHadronType(Int_t t)
TSystem * gSystem
Declaration of class AliTLorentzVector.
AliClusterContainer * GetClusterContainer() const
Int_t fCentBin
!event centrality bin
std::vector< Int_t > SortTracksPt(AliVEvent *event)
void SetPercAreaCut(Float_t p)
void SetVzRange(Double_t min, Double_t max)
Set pre-configured event cut object.
void jetContSetParams(AliJetContainer *jetCont, Float_t jetptmin=1., Float_t jetptmax=500., Float_t jetareacutperc=0., Int_t leadhadtype=2, Int_t nLeadJets=1, Float_t mintrackpt=0.15, Float_t maxtrackpt=1000.)
AliClusterContainer * AddClusterContainer(const char *n)
Create new cluster container and attach it to the task.
UShort_t GetNumberOfConstituents() const
Definition: AliEmcalJet.h:140
Container for particles within the EMCAL framework.
Int_t TrackAt(Int_t idx) const
Definition: AliEmcalJet.h:160
UShort_t GetNumberOfTracks() const
Definition: AliEmcalJet.h:139
std::ostream & Print(std::ostream &in) const
AliParticleContainer * GetParticleContainer() const
AliEmcalJet * GetLeadingJet(const char *opt="")
void SetMinTrackPt(Float_t b)
TH2 * CreateTH2(const char *name, const char *title, int nbinsx, double xmin, double xmax, int nbinsy, double ymin, double ymax, Option_t *opt="")
Create a new TH2 within the container.
TObject * FindObject(const char *name) const
Find an object inside the container.
int Int_t
Definition: External.C:63
void SetJetPtCut(Float_t cut)
unsigned int UInt_t
Definition: External.C:33
THashList * GetListOfHistograms() const
Get the list of histograms.
Definition: THistManager.h:671
float Float_t
Definition: External.C:68
void SetNLeadingJets(Int_t t)
Double_t GetMinPt() const
TObject * GetHistogram(const char *histName)
Int_t GetNJets() const
AliParticleContainer * AddParticleContainer(const char *n)
Create new particle container and attach it to the task.
TH1 * CreateTH1(const char *name, const char *title, int nbins, double xmin, double xmax, Option_t *opt="")
Create a new TH1 within the container.
Double_t Z_pt(const AliEmcalJet *jet, const AliVParticle *trk)
Definition: External.C:228
Int_t fNcentBins
how many centrality bins
Definition: External.C:212
AliAnalysisTaskEmcalJetCDF * AddTaskEmcalJetCDF(const char *ntracks="usedefault", const char *nclusters="usedefault", const char *ncells="usedefault", const char *tag="CDF")
Int_t GetNAcceptedClusters() const
TClonesArray * GetArray() const
AliMCParticleContainer * AddMCParticleContainer(const char *n)
Create new container for MC particles and attach it to the task.
Double_t GetZ(const Double_t trkPx, const Double_t trkPy, const Double_t trkPz) const
TObjArray fJetCollArray
jet collection array
Double_t DeltaR(const AliVParticle *part) const
Double_t Perp(const AliVParticle &trk1, const AliVParticle &trk2)
Double_t Pt() const
Definition: AliEmcalJet.h:109
const char * GetName() const
AliEmcalList * fOutput
!output list
virtual ~AliAnalysisTaskEmcalJetCDF()
Destructor.
AliTrackContainer * AddTrackContainer(const char *n)
Create new track container and attach it to the task.
TFile * file
TList with histograms for a given trigger.
void SetMakeGeneralHistograms(Bool_t g)
Enable general histograms.
Base task in the EMCAL jet framework.
Represent a jet reconstructed using the EMCal jet framework.
Definition: AliEmcalJet.h:51
unsigned short UShort_t
Definition: External.C:28
functional for sorting pair by first element - descending
const char Option_t
Definition: External.C:48
void UserCreateOutputObjects()
Main initialization function on the worker.
Int_t GetNEntries() const
const AliJetIterableContainer accepted() const
bool Bool_t
Definition: External.C:53
Declaration of class AliAnalysisTaskEmcalJetCDF.
EDataType_t
Switch for the data type.
void SetCaloCellsName(const char *n)
Int_t GetNAcceptedParticles() const
void SetMaxTrackPt(Float_t b)
Container structure for EMCAL clusters.
TChain * CreateChain(const char *filelist="filelist.txt", const char *cTreeNameArg="auto", const char *friends="", UInt_t iNumFiles=-1, UInt_t iStartWithFile=1)
Double_t GetMaxPt() const
Double_t Z_ptot(const AliEmcalJet *jet, const AliVParticle *trk)
Container for jet within the EMCAL jet framework.
void SetJetPtCutMax(Float_t cut)
std::vector< int > GetPtSortedTrackConstituentIndexes(TClonesArray *tracks) const
Analysis of jet shapes and FF of all jets and leading jets.