AliPhysics  608b256 (608b256)
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 #include <iostream>
18 
19 #include <TSystem.h>
20 #include <TObject.h>
21 #include <TObjArray.h>
22 #include <TObjString.h>
23 #include <TString.h>
24 #include <TFile.h>
25 #include <TKey.h>
26 #include <TChain.h>
27 #include <TFileCollection.h>
28 #include <TCollection.h>
29 #include <THashList.h>
30 #include <TRegexp.h>
31 #include <TFileInfo.h>
32 
33 #include <TClonesArray.h>
34 #include <TH1D.h>
35 #include <TH2D.h>
36 #include <TArrayD.h>
37 #include <TString.h>
38 
39 
40 #include <AliVCluster.h>
41 #include <AliVParticle.h>
42 #include <AliLog.h>
43 
44 #include "AliTLorentzVector.h"
45 #include "AliEmcalJet.h"
46 #include "AliRhoParameter.h"
47 #include "AliJetContainer.h"
48 #include "AliParticleContainer.h"
49 #include "AliClusterContainer.h"
50 #include "AliAnalysisManager.h"
51 #include "AliVEventHandler.h"
52 #include "AliAnalysisDataContainer.h"
53 
54 #include "AliEmcalList.h"
55 
57 
58 using std::cout;
59 using std::endl;
60 
62 ClassImp ( AliAnalysisTaskEmcalJetCDF );
64 
70  fHistManager()
71 {}
72 
79  AliAnalysisTaskEmcalJet ( name, kTRUE ),
80  fHistManager(name)
81  {
82  // Standard constructor.
83  SetMakeGeneralHistograms ( kTRUE );
84  }
85 
88 
97  {
98  return kTRUE;
99  }
100 
101 //________________________________________________________________________
103  {
104  TH1::SetDefaultSumw2(kTRUE);
105  TH2::SetDefaultSumw2(kTRUE);
106 
108  TString histname = "", groupname = "", fullgroupname = "";
109 
110  AliJetContainer* jetCont = NULL;
111  TIter next(&fJetCollArray);
112  while ( (jetCont = static_cast<AliJetContainer*>(next())) ) {
113  //##### EARLY VALIDITY CHECKS - BAIL OUT FAST
114  // get particles connected to jets
115  AliParticleContainer* fTracksCont = dynamic_cast<AliParticleContainer*>(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->GetNAcceptedJets();
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  if ( fDebug > 1 ) { std::cout << "fNJets = " << fNJets_accepted << " ; fNPart = " << fNaccPart << std::endl; }
131 
132 
133  // Only fill the embedding qa plots if:
134  // - We are using the embedding helper
135  // - The class has been initialized
136  // - Both jet collections are available
137 // if (fEmbeddingQA.IsInitialized()) { fEmbeddingQA.RecordEmbeddedEventProperties(); }
138 
139 //######################################################################################################
140  groupname = jetCont->GetName();
141 
142  Double_t jet_pt_min = jetCont->GetMinPt();
143  Double_t jet_pt_max = jetCont->GetMaxPt();
144 
145  TString jetstrmin = TString::Itoa((Int_t)jet_pt_min,10);
146  TString jetstrmax = TString::Itoa((Int_t)jet_pt_max,10);
147 
148  // add to groupname the min,max pt cuts of jets in the container
149  groupname = groupname + "_" + "ptbin" + "_" + jetstrmin + "_" + jetstrmax;
150 
151 //######################################################################################################
152 // Get histo pointers from Hist Manager
153  histname = TString::Format("%s/histo1_%d", groupname.Data(), fCentBin);
154  TH1D* fH1 = (TH1D*)GetHistogram(histname.Data());
155 
156  histname = TString::Format("%s/histo2_%d", groupname.Data(), fCentBin);
157  TH1D* fH2 = (TH1D*)GetHistogram(histname.Data());
158 
159  histname = TString::Format("%s/histo3_%d", groupname.Data(), fCentBin);
160  TH1D* fH3 = (TH1D*)GetHistogram(histname.Data());
161 
162  histname = TString::Format("%s/histo4_%d", groupname.Data(), fCentBin);
163  TH1D* fH4 = (TH1D*)GetHistogram(histname.Data());
164 
165  histname = TString::Format("%s/histo4c_%d", groupname.Data(), fCentBin);
166  TH1D* fH4c = (TH1D*)GetHistogram(histname.Data());
167 
168  histname = TString::Format("%s/histo5_%d", groupname.Data(), fCentBin);
169  TH1D* fH5 = (TH1D*)GetHistogram(histname.Data());
170 
171  histname = TString::Format("%s/histo7all_%d", groupname.Data(), fCentBin);
172  TH2D* fH7all = (TH2D*)GetHistogram(histname.Data());
173 //######################################################################################################
174  histname = TString::Format("%s/histo8_all_pt_%d", groupname.Data(), fCentBin);
175  TH1D* fH8_all_pt = (TH1D*)GetHistogram(histname.Data());
176 
177  histname = TString::Format("%s/histo8xi_all_pt_%d", groupname.Data(), fCentBin);
178  TH1D* fH8xi_all_pt = (TH1D*)GetHistogram(histname.Data());
179 //######################################################################################################
180 
181  histname = TString::Format("%s/histo15all_%d", groupname.Data(), fCentBin);
182  TH2D* fH15all = (TH2D*)GetHistogram(histname.Data());
183 
184  histname = TString::Format("%s/histo15all_n80_%d", groupname.Data(), fCentBin);
185  TH2D* fH15all_n80 = (TH2D*)GetHistogram(histname.Data());
186 
187  histname = TString::Format("%s/histo15all_n90_%d", groupname.Data(), fCentBin);
188  TH2D* fH15all_n90 = (TH2D*)GetHistogram(histname.Data());
189 
190  histname = TString::Format("%s/histo15all_pt80_%d", groupname.Data(), fCentBin);
191  TH2D* fH15all_pt80 = (TH2D*)GetHistogram(histname.Data());
192 
193  histname = TString::Format("%s/histo15all_pt90_%d", groupname.Data(), fCentBin);
194  TH2D* fH15all_pt90 = (TH2D*)GetHistogram(histname.Data());
195 //######################################################################################################
196 
197  histname = TString::Format("%s/histo20all_%d", groupname.Data(), fCentBin);
198  TH1D* fH20all = (TH1D*)GetHistogram(histname.Data());
199 
200  histname = TString::Format("%s/histo20all_n80_%d", groupname.Data(), fCentBin);
201  TH1D* fH20all_n80 = (TH1D*)GetHistogram(histname.Data());
202 
203  histname = TString::Format("%s/histo20all_n90_%d", groupname.Data(), fCentBin);
204  TH1D* fH20all_n90 = (TH1D*)GetHistogram(histname.Data());
205 
206  histname = TString::Format("%s/histo20all_pt80_%d", groupname.Data(), fCentBin);
207  TH1D* fH20all_pt80 = (TH1D*)GetHistogram(histname.Data());
208 
209  histname = TString::Format("%s/histo20all_pt90_%d", groupname.Data(), fCentBin);
210  TH1D* fH20all_pt90 = (TH1D*)GetHistogram(histname.Data());
211 //######################################################################################################
212 
213  histname = TString::Format("%s/histo_g_%d", groupname.Data(), fCentBin);
214  TH1D* fHg = (TH1D*)GetHistogram(histname.Data());
215 
216  histname = TString::Format("%s/histo_g_n80_%d", groupname.Data(), fCentBin);
217  TH1D* fHg_n80 = (TH1D*)GetHistogram(histname.Data());
218 
219  histname = TString::Format("%s/histo_g_n90_%d", groupname.Data(), fCentBin);
220  TH1D* fHg_n90 = (TH1D*)GetHistogram(histname.Data());
221 
222  histname = TString::Format("%s/histo_g_pt80_%d", groupname.Data(), fCentBin);
223  TH1D* fHg_pt80 = (TH1D*)GetHistogram(histname.Data());
224 
225  histname = TString::Format("%s/histo_g_pt90_%d", groupname.Data(), fCentBin);
226  TH1D* fHg_pt90 = (TH1D*)GetHistogram(histname.Data());
227 //######################################################################################################
228 
229  histname = TString::Format("%s/histo_ptd_%d", groupname.Data(), fCentBin);
230  TH1D* fHptd = (TH1D*)GetHistogram(histname.Data());
231 
232  histname = TString::Format("%s/histo_ptd_n80_%d", groupname.Data(), fCentBin);
233  TH1D* fHptd_n80 = (TH1D*)GetHistogram(histname.Data());
234 
235  histname = TString::Format("%s/histo_ptd_n90_%d", groupname.Data(), fCentBin);
236  TH1D* fHptd_n90 = (TH1D*)GetHistogram(histname.Data());
237 
238  histname = TString::Format("%s/histo_ptd_pt80_%d", groupname.Data(), fCentBin);
239  TH1D* fHptd_pt80 = (TH1D*)GetHistogram(histname.Data());
240 
241  histname = TString::Format("%s/histo_ptd_pt90_%d", groupname.Data(), fCentBin);
242  TH1D* fHptd_pt90 = (TH1D*)GetHistogram(histname.Data());
243 //######################################################################################################
244 
245  histname = TString::Format("%s/histo_Rjt_%d", groupname.Data(), fCentBin);
246  TH2D* fH_Rjt = (TH2D*)GetHistogram(histname.Data());
247 
248  histname = TString::Format("%s/histo_Rjt_n80_%d", groupname.Data(), fCentBin);
249  TH2D* fH_Rjt_n80 = (TH2D*)GetHistogram(histname.Data());
250 
251  histname = TString::Format("%s/histo_Rjt_n90_%d", groupname.Data(), fCentBin);
252  TH2D* fH_Rjt_n90 = (TH2D*)GetHistogram(histname.Data());
253 
254  histname = TString::Format("%s/histo_Rjt_pt80_%d", groupname.Data(), fCentBin);
255  TH2D* fH_Rjt_pt80 = (TH2D*)GetHistogram(histname.Data());
256 
257  histname = TString::Format("%s/histo_Rjt_pt90_%d", groupname.Data(), fCentBin);
258  TH2D* fH_Rjt_pt90 = (TH2D*)GetHistogram(histname.Data());
259 //######################################################################################################
260 
261  histname = TString::Format("%s/histo_jt_%d", groupname.Data(), fCentBin);
262  TH1D* fH_jt = (TH1D*)GetHistogram(histname.Data());
263 
264  histname = TString::Format("%s/histo_jt_n80_%d", groupname.Data(), fCentBin);
265  TH1D* fH_jt_n80 = (TH1D*)GetHistogram(histname.Data());
266 
267  histname = TString::Format("%s/histo_jt_n90_%d", groupname.Data(), fCentBin);
268  TH1D* fH_jt_n90 = (TH1D*)GetHistogram(histname.Data());
269 
270  histname = TString::Format("%s/histo_jt_pt80_%d", groupname.Data(), fCentBin);
271  TH1D* fH_jt_pt80 = (TH1D*)GetHistogram(histname.Data());
272 
273  histname = TString::Format("%s/histo_jt_pt90_%d", groupname.Data(), fCentBin);
274  TH1D* fH_jt_pt90 = (TH1D*)GetHistogram(histname.Data());
275 
276 
277 //######################################################################################################
278  // get clusters connected to jets
279  AliClusterContainer* fCaloClustersCont = jetCont->GetClusterContainer();
280  // accepted clusters in cluster container
281  Int_t fNaccClus = -1;
282  if (fCaloClustersCont) { fNaccClus = fCaloClustersCont->GetNAcceptedClusters(); }
283 
284  fH5->Fill ( fNJets_accepted ); // Distribution of jets in events;
285 
286  UShort_t counter_part = 0; Double_t counter_pt = 0.; // counter for npart and pt recording
287  UShort_t jet_n90 = -99 ; Double_t jet_pt90 = -99.99 ;
288  UShort_t jet_n80 = -99 ; Double_t jet_pt80 = -99.99 ;
289 
290  // variables used to compute g and ptD
291  Double_t g_tot = 0.; Double_t sum_part_pt_tot = 0.; Double_t sum_part_pt2_tot = 0.;
292  Double_t g_n90 = 0.; Double_t sum_part_pt_n90 = 0.; Double_t sum_part_pt2_n90 = 0.;
293  Double_t g_n80 = 0.; Double_t sum_part_pt_n80 = 0.; Double_t sum_part_pt2_n80 = 0.;
294  Double_t g_pt90 = 0.; Double_t sum_part_pt_pt90 = 0.; Double_t sum_part_pt2_pt90 = 0.;
295  Double_t g_pt80 = 0.; Double_t sum_part_pt_pt80 = 0.; Double_t sum_part_pt2_pt80 = 0.;
296 
297 
298  // **************************************************************
299  // ALL JETS
300  // **************************************************************
301  Double_t jet_pt = 0. ; UShort_t jet_npart = 0; UShort_t jet_nconst = 0;
302 
303  // loop over all jets
304  for( auto jet : jetCont->accepted()) {
305  if (!jet) { continue; }
306 
307  // vector of sorted indexes of particles in jet
308  std::vector< int > jet_sorted_idxvec ;
309 
310  Int_t track_idx = -999; // index variable for tracks
311  jet_pt = 0. ; jet_npart = 0; jet_nconst = 0;
312  counter_part = 0; counter_pt = 0.; // reset counters
313 
314  // jet : Sorting by p_T jet constituents
315  jet_sorted_idxvec.clear();
316  jet_sorted_idxvec = jet->GetPtSortedTrackConstituentIndexes ( fTracksContArray );
317 
318  jet_pt = jet->Pt();
319  jet_npart = jet->GetNumberOfTracks();
320  jet_nconst = jet->GetNumberOfConstituents();
321 
322  // variables for g and pdt computations
323  g_tot = 0.; sum_part_pt_tot = 0.; sum_part_pt2_tot = 0.;
324  g_n90 = 0.; sum_part_pt_n90 = 0.; sum_part_pt2_n90 = 0.;
325  g_n80 = 0.; sum_part_pt_n80 = 0.; sum_part_pt2_n80 = 0.;
326  g_pt90 = 0.; sum_part_pt_pt90 = 0.; sum_part_pt2_pt90 = 0.;
327  g_pt80 = 0.; sum_part_pt_pt80 = 0.; sum_part_pt2_pt80 = 0.;
328 
329  // sentinels for the pt tail cuts
330  jet_n90 = ( UShort_t ) ( 0.9 * jet_npart );
331  jet_pt90 = 0.9 * jet_pt;
332 
333  jet_n80 = ( UShort_t ) ( 0.8 * jet_npart );
334  jet_pt80 = 0.8 * jet_pt;
335 
336  fH1->Fill ( jet_pt ); // Pt distribution of jets
337  fH2->Fill ( jet->Eta() ); // Eta distribution of jets
338  fH3->Fill ( jet->Phi() ); // Phi distribution of jets
339  fH4->Fill ( jet_npart ); // Multiplicity of jets
340  fH4c->Fill ( jet_nconst ); // Multiplicity of jets - all constituents
341  fH7all->Fill ( jet_pt, jet_nconst ); // N(jet) vs P_{T} - all jets
342 
343  // parsing all jet tracks
344  for (std::size_t i = 0; i < jet_npart; i++ ) {
345  track_idx = jet_sorted_idxvec.at (i);
346  AliVParticle* track = jet->TrackAt ( track_idx, fTracksContArray );
347  if (!track) { std::cout << "Parsing tracks of jets :: track not defined but it should!!!" << std::endl; continue; }
348 
349  Double_t dpart = jet->DeltaR ( track );
350  Double_t track_pt = track->Pt();
351  Double_t jt = CDF::Perp (*track, *jet);
352 
353  // https://arxiv.org/abs/1209.6086 pag 4, (1)
354  Double_t z_pt = CDF::Z_pt(jet,track);
355  fH8_all_pt->Fill ( z_pt ); // Momentum distribution for jets (FF)
356  fH8xi_all_pt->Fill ( CDF::Xi (z_pt) ); // Momentum distribution for jets (FF) xi
357 
358  fH15all->Fill ( dpart, track_pt, track_pt ); // p_{T} track vs the Distance R from jet
359  fH20all->Fill ( dpart ); // Distribution of R in leading jet
360 
361  fH_Rjt->Fill ( dpart, jt, jt ); // jt track vs dR weighted with jt
362  fH_jt->Fill ( jt ); // jt track vs dR
363 
364  // computing components for g and ptD in the jet tracks loop
365  g_tot += (track_pt * dpart)/jet_pt;
366  sum_part_pt_tot += TMath::Abs(track_pt);
367  sum_part_pt2_tot += TMath::Power( track_pt, 2 );
368 
369  //#############################################################################################
370  if ( counter_part <= jet_n90 ) {
371  fH15all_n90->Fill ( dpart, track_pt ); // p_{T} track vs the Distance R from Jet - 80% of particles
372  fH20all_n90->Fill ( dpart );
373  fH_Rjt_n90->Fill ( dpart, jt, jt );
374  fH_jt_n90->Fill ( jt ); // jt track vs dR
375 
376  // computing components for g and ptD in the jet tracks loop
377  g_n90 += (track_pt * dpart)/jet_pt;
378  sum_part_pt_n90 += track_pt;
379  sum_part_pt2_n90 += TMath::Power( track_pt, 2 );
380  }
381 
382  if ( counter_pt <= jet_pt90 ) {
383  fH15all_pt90->Fill ( dpart, track_pt ); // p_{T} track vs the Distance R from Jet - 80% of pt
384  fH20all_pt90->Fill ( dpart );
385  fH_Rjt_pt90->Fill ( dpart, jt, jt );
386  fH_jt_pt90->Fill ( jt ); // jt track vs dR
387 
388  // computing components for g and ptD in the jet tracks loop
389  g_pt90 += (track_pt * dpart)/jet_pt;
390  sum_part_pt_pt90 += track_pt;
391  sum_part_pt2_pt90 += TMath::Power( track_pt, 2 );
392  }
393 
394  //#############################################################################################
395  if ( counter_part <= jet_n80 ) {
396  fH15all_n80->Fill ( dpart, track_pt ); // p_{T} track vs the Distance R from Jet - 80% of particles
397  fH20all_n80->Fill ( dpart );
398  fH_Rjt_n80->Fill ( dpart, jt, jt );
399  fH_jt_n80->Fill ( jt ); // jt track vs dR
400 
401  // computing components for g and ptD in the jet tracks loop
402  g_n80 += (track_pt * dpart)/jet_pt;
403  sum_part_pt_n80 += track_pt;
404  sum_part_pt2_n80 += TMath::Power( track_pt, 2 );
405  }
406 
407  if ( counter_pt <= jet_pt80 ) {
408  fH15all_pt80->Fill ( dpart, track_pt ); // p_{T} track vs the Distance R from Jet - 80% of pt
409  fH20all_pt80->Fill ( dpart );
410  fH_Rjt_pt80->Fill ( dpart, jt, jt );
411  fH_jt_pt80->Fill ( jt ); // jt track vs dR
412 
413  // computing components for g and ptD in the jet tracks loop
414  g_pt80 += (track_pt * dpart)/jet_pt;
415  sum_part_pt_pt80 += track_pt;
416  sum_part_pt2_pt80 += TMath::Power( track_pt, 2 );
417  }
418 
419  ++counter_part; counter_pt += track_pt;
420  } // end of loop over jet tracks
421 
422  fHg->Fill ( g_tot );
423  fHg_n80->Fill ( g_n80 ); fHg_pt80->Fill ( g_pt80 );
424  fHg_n90->Fill ( g_n90 ); fHg_pt90->Fill ( g_pt90 );
425 
426  if ( sum_part_pt_tot > 1e-8 )
427  { fHptd->Fill( TMath::Sqrt(sum_part_pt2_tot)/sum_part_pt_tot ); }
428  else
429  { if ( fDebug > 2 ) cout << "sum_part_pt_tot aprox 0!!!!" << endl; }
430 
431  if ( sum_part_pt_n80 > 1e-8 )
432  { fHptd_n80->Fill( TMath::Sqrt(sum_part_pt2_n80)/sum_part_pt_n80 ); }
433  else
434  { if ( fDebug > 2 ) cout << "sum_part_pt_n80 aprox 0!!!!" << endl; }
435 
436  if ( sum_part_pt_n90 > 1e-8 )
437  { fHptd_n90->Fill( TMath::Sqrt(sum_part_pt2_n90)/sum_part_pt_n90 ); }
438  else
439  { if ( fDebug > 2 ) cout << "sum_part_pt_n90 aprox 0!!!!" << endl; }
440 
441  if ( sum_part_pt_pt80 > 1e-8 )
442  { fHptd_pt80->Fill( TMath::Sqrt(sum_part_pt2_pt80)/sum_part_pt_pt80 ); }
443  else
444  { if ( fDebug > 2 ) cout << "sum_part_pt_pt80 aprox 0!!!!" << endl; }
445 
446  if ( sum_part_pt_pt90 > 1e-8 )
447  { fHptd_pt90->Fill( TMath::Sqrt(sum_part_pt2_pt90)/sum_part_pt_pt90 ); }
448  else
449  { if ( fDebug > 2 ) cout << "sum_part_pt_pt90 aprox 0!!!!" << endl; }
450 
451  } // end of loopt over all jets
452  } // end of loop over jet container collection
453 
454  // post data at every processing
455  PostData ( 1, fOutput ); // Post data for ALL output slots > 0 here.
456  return kTRUE;
457  }
458 
459 //________________________________________________________________________
461  {
462  // Create user output.
464 
465  TString histname = "", histtitle = "", groupname = "", fullgroupname = "";
466  AliJetContainer* jetCont = 0;
467  TIter next(&fJetCollArray);
468  while ((jetCont = static_cast<AliJetContainer*>(next()))) {
469  groupname = jetCont->GetName();
470 
471  Double_t jet_pt_min = jetCont->GetMinPt();
472  Double_t jet_pt_max = jetCont->GetMaxPt();
473 
474  TString jetstrmin = TString::Itoa((Int_t)jet_pt_min,10);
475  TString jetstrmax = TString::Itoa((Int_t)jet_pt_max,10);
476 
477  // add to groupname the min,max pt cuts of jets in the container
478  groupname = groupname + "_" + "ptbin" + "_" + jetstrmin + "_" + jetstrmax;
479 
480  fHistManager.CreateHistoGroup(groupname);
481  for (Int_t cent = 0; cent < fNcentBins; cent++) {
482  //=====================================================================================
483  Int_t h1_nbin = 200; Double_t h1_binwidth = 1; Double_t h1_low = 0;
484  Double_t h1_high = h1_low + h1_binwidth * h1_nbin; // 1GeV/bin
485  histname = TString::Format("%s/histo1_%d", groupname.Data(), cent);
486  histtitle = TString::Format("%s;#it{p}_{T,jet} (GeV/#it{c}) (accepted);Jets", histname.Data()); // Pt distro of jets
487  fHistManager.CreateTH1(histname, histtitle, h1_nbin, h1_low, h1_high);
488 
489  //=====================================================================================
490  Int_t h2_nbin = 200; Double_t h2_binwidth = 0.01; Double_t h2_low = -1;
491  Double_t h2_high = h2_low + h2_binwidth * h2_nbin;
492  histname = TString::Format("%s/histo2_%d", groupname.Data(), cent);
493  histtitle = TString::Format("%s;#it{#eta}_{jet};Jets", histname.Data()); // Eta distro of jets
494  fHistManager.CreateTH1(histname, histtitle, h2_nbin, h2_low, h2_high); // 1 unit of rapidity / 100 bin
495 
496  //=====================================================================================
497  Int_t h3_nbin = 126; Double_t h3_binwidth = 0.05; Double_t h3_low = 0.;
498  Double_t h3_high = h3_low + h3_binwidth * h3_nbin;
499  histname = TString::Format("%s/histo3_%d", groupname.Data(), cent);
500  histtitle = TString::Format("%s;#it{#phi}_{jet};Jets", histname.Data()); // Phi distro of jets
501  fHistManager.CreateTH1(histname, histtitle, h3_nbin, h3_low, h3_high);
502 
503  //=====================================================================================
504  Int_t h4_nbin = 100; Double_t h4_binwidth = 1; Double_t h4_low = 0;
505  Double_t h4_high = h4_low + h4_binwidth * h4_nbin; // 1 unit of multiplicity /bin
506  histname = TString::Format("%s/histo4_%d", groupname.Data(), cent);
507  histtitle = TString::Format("%s;N_{tracks}(jet);Jets", histname.Data()); // Multiplicity of all jets; chg tracks
508  fHistManager.CreateTH1(histname, histtitle, h4_nbin, h4_low, h4_high);
509 
510  histname = TString::Format("%s/histo4c_%d", groupname.Data(), cent);
511  histtitle = TString::Format("%s;N_{tracks}(jet);Jets", histname.Data());
512  fHistManager.CreateTH1(histname, histtitle, h4_nbin, h4_low, h4_high); // Multiplicity of all jets; all tracks
513 
514  //#####################################
515 
516  //=====================================================================================
517  Int_t h5_nbin = 100; Double_t h5_binwidth = 1; Double_t h5_low = 0;
518  Double_t h5_high = h5_low + h5_binwidth * h5_nbin;
519  histname = TString::Format("%s/histo5_%d", groupname.Data(), cent);
520  histtitle = TString::Format("%s;N_{jets};Events", histname.Data()); // Distribution of jets in events
521  fHistManager.CreateTH1(histname, histtitle, h5_nbin, h5_low, h5_high);
522 
523  //=====================================================================================
524  Int_t h7_xnbin = 300; Double_t h7_xbinwidth = 1; Double_t h7_xlow = 0;
525  Double_t h7_xhigh = h7_xlow + h7_xbinwidth * h7_xnbin;
526  Int_t h7_ynbin = 100; Double_t h7_ybinwidth = 1; Double_t h7_ylow = 0;
527  Double_t h7_yhigh = h7_ylow + h7_ybinwidth * h7_ynbin;
528 
529  histname = TString::Format("%s/histo7all_%d", groupname.Data(), cent);
530  histtitle = TString::Format("%s;#it{p}_{T,jet} (GeV/c);N_{tracks}(jets)", histname.Data()); // N vs pt all jets
531  fHistManager.CreateTH2(histname, histtitle, h7_xnbin, h7_xlow, h7_xhigh, h7_ynbin, h7_ylow, h7_yhigh);
532 
533  //=====================================================================================
534  Int_t h8_nbin = 101; Double_t h8_binwidth = 0.01; Double_t h8_low = 0;
535  Double_t h8_high = h8_low + h8_binwidth * h8_nbin;
536 
537  // Pt implementation of Z
538  histname = TString::Format("%s/histo8_all_pt_%d", groupname.Data(), cent);
539  histtitle = TString::Format("%s - all jets Pt;z = p_{T,track}/p_{T,jet1};F(Z) = 1/N_{jets} dN/dz", histname.Data());
540  fHistManager.CreateTH1(histname, histtitle, h8_nbin, h8_low, h8_high);
541 
542  //=====================================================================================
543  Int_t h8xi_nbin = 140; Double_t h8xi_binwidth = 0.05; Double_t h8xi_low = 0;
544  Double_t h8xi_high = h8xi_low + h8xi_binwidth * h8xi_nbin;
545 
546  histname = TString::Format("%s/histo8xi_all_pt_%d", groupname.Data(), cent);
547  histtitle = TString::Format("%s - all jets Pt;#xi = ln(1/z);D(#xi) = 1/N_{jets} dN/d#xi", histname.Data());
548  fHistManager.CreateTH1(histname, histtitle, h8xi_nbin, h8xi_low, h8xi_high);
549 
550  //=====================================================================================
551  Int_t h15_xnbin = 60; Double_t h15_xbinwidth = 0.01; Double_t h15_xlow = 0.;
552  Double_t h15_xhigh = h15_xlow + h15_xbinwidth * h15_xnbin;
553  Int_t h15_ynbin = 300; Double_t h15_ybinwidth = 1.; Double_t h15_ylow = 0.;
554  Double_t h15_yhigh = h15_ylow + h15_ybinwidth * h15_ynbin;
555 
556  histname = TString::Format("%s/histo15all_%d", groupname.Data(), cent);
557  histtitle = TString::Format("%s - all jets;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
558  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
559 
560  //########################################################
561  histname = TString::Format("%s/histo15all_n80_%d", groupname.Data(), cent);
562  histtitle = TString::Format("%s - all jets;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
563  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
564 
565  histname = TString::Format("%s/histo15all_n90_%d", groupname.Data(), cent);
566  histtitle = TString::Format("%s - all jets;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
567  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
568 
569  //########################################################
570  histname = TString::Format("%s/histo15all_pt80_%d", groupname.Data(), cent);
571  histtitle = TString::Format("%s - all jets;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
572  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
573 
574  histname = TString::Format("%s/histo15all_pt90_%d", groupname.Data(), cent);
575  histtitle = TString::Format("%s - all jets;dR;#it{p}_{T,track} (GeV/c)", histname.Data()); // p_T track vs dR
576  fHistManager.CreateTH2(histname, histtitle, h15_xnbin, h15_xlow, h15_xhigh, h15_ynbin, h15_ylow, h15_yhigh);
577 
578  //=====================================================================================
579  Int_t h20_nbin = 60; Double_t h20_binwidth = 0.01; Double_t h20_low = 0.;
580  Double_t h20_high = h20_low + h20_binwidth * h20_nbin;
581 
582  //########################################################
583  histname = TString::Format("%s/histo20all_%d", groupname.Data(), cent);
584  histtitle = TString::Format("%s - all jets;R_{tracks};dN/dR", histname.Data());
585  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
586 
587  //########################################################
588  histname = TString::Format("%s/histo20all_n80_%d", groupname.Data(), cent);
589  histtitle = TString::Format("%s - all jets;R_{tracks};dN/dR", histname.Data());
590  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
591 
592  histname = TString::Format("%s/histo20all_n90_%d", groupname.Data(), cent);
593  histtitle = TString::Format("%s - all jets;R_{tracks};dN/dR", histname.Data());
594  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
595 
596  //########################################################
597  histname = TString::Format("%s/histo20all_pt80_%d", groupname.Data(), cent);
598  histtitle = TString::Format("%s - all jets;R_{tracks};dN/dR", histname.Data());
599  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
600 
601  histname = TString::Format("%s/histo20all_pt90_%d", groupname.Data(), cent);
602  histtitle = TString::Format("%s - all jets;R_{tracks};dN/dR", histname.Data());
603  fHistManager.CreateTH1(histname, histtitle, h20_nbin, h20_low, h20_high);
604 
605  //=====================================================================================
606  // Distribution of girth (radial girth) g = sum_jet_parts ( r_i * ( pt_i/pt_jet ) )
607  Int_t hg_nbin = 40; Double_t hg_binwidth = 0.005; Double_t hg_low = 0.;
608  Double_t hg_high = hg_low + hg_binwidth * hg_nbin;
609 
610  //########################################################
611  histname = TString::Format("%s/histo_g_%d", groupname.Data(), cent);
612  histtitle = TString::Format("%s - all jets;g;1/N_{jets} dN/dg", histname.Data());
613  fHistManager.CreateTH1(histname, histtitle, hg_nbin, hg_low, hg_high);
614 
615  //########################################################
616  histname = TString::Format("%s/histo_g_n80_%d", groupname.Data(), cent);
617  histtitle = TString::Format("%s - all jets;g;1/N_{jets} dN/dg", histname.Data());
618  fHistManager.CreateTH1(histname, histtitle, hg_nbin, hg_low, hg_high);
619 
620  histname = TString::Format("%s/histo_g_n90_%d", groupname.Data(), cent);
621  histtitle = TString::Format("%s - all jets;g;1/N_{jets} dN/dg", histname.Data());
622  fHistManager.CreateTH1(histname, histtitle, hg_nbin, hg_low, hg_high);
623 
624  //########################################################
625  histname = TString::Format("%s/histo_g_pt80_%d", groupname.Data(), cent);
626  histtitle = TString::Format("%s - all jets;g;1/N_{jets} dN/dg", histname.Data());
627  fHistManager.CreateTH1(histname, histtitle, hg_nbin, hg_low, hg_high);
628 
629  histname = TString::Format("%s/histo_g_pt90_%d", groupname.Data(), cent);
630  histtitle = TString::Format("%s - all jets;g;1/N_{jets} dN/dg", histname.Data());
631  fHistManager.CreateTH1(histname, histtitle, hg_nbin, hg_low, hg_high);
632 
633  //=====================================================================================
634  // Distribution of dispersion d pt_D = sqrt ( sum (pt_i^2) )/sum (pt_i)
635  Int_t hptd_nbin = 110; Double_t hptd_binwidth = 0.01; Double_t hptd_low = 0.;
636  Double_t hptd_high = hptd_low + hptd_binwidth * hptd_nbin;
637 
638  //########################################################
639  histname = TString::Format("%s/histo_ptd_%d", groupname.Data(), cent);
640  histtitle = TString::Format("%s - all jets;ptd;1/N_{jets} dN/dp_{T}D", histname.Data());
641  fHistManager.CreateTH1(histname, histtitle, hptd_nbin, hptd_low, hptd_high);
642 
643  //########################################################
644  histname = TString::Format("%s/histo_ptd_n80_%d", groupname.Data(), cent);
645  histtitle = TString::Format("%s - all jets;ptd;1/N_{jets} dN/dp_{T}D", histname.Data());
646  fHistManager.CreateTH1(histname, histtitle, hptd_nbin, hptd_low, hptd_high);
647 
648  histname = TString::Format("%s/histo_ptd_n90_%d", groupname.Data(), cent);
649  histtitle = TString::Format("%s - all jets;ptd;1/N_{jets} dN/dp_{T}D", histname.Data());
650  fHistManager.CreateTH1(histname, histtitle, hptd_nbin, hptd_low, hptd_high);
651  //########################################################
652 
653  histname = TString::Format("%s/histo_ptd_pt80_%d", groupname.Data(), cent);
654  histtitle = TString::Format("%s - all jets;ptd;1/N_{jets} dN/dp_{T}D", histname.Data());
655  fHistManager.CreateTH1(histname, histtitle, hptd_nbin, hptd_low, hptd_high);
656 
657  histname = TString::Format("%s/histo_ptd_pt90_%d", groupname.Data(), cent);
658  histtitle = TString::Format("%s - all jets;ptd;1/N_{jets} dN/dp_{T}D", histname.Data());
659  fHistManager.CreateTH1(histname, histtitle, hptd_nbin, hptd_low, hptd_high);
660 
661  //=====================================================================================
662  Int_t h_Rjt_xnbin = 60; Double_t h_Rjt_xbinwidth = 0.01; Double_t h_Rjt_xlow = 0.;
663  Double_t h_Rjt_xhigh = h_Rjt_xlow + h_Rjt_xbinwidth * h_Rjt_xnbin;
664  Int_t h_Rjt_ynbin = 500; Double_t h_Rjt_ybinwidth = 0.01; Double_t h_Rjt_ylow = 0.;
665  Double_t h_Rjt_yhigh = h_Rjt_ylow + h_Rjt_ybinwidth * h_Rjt_ynbin;
666 
667  //########################################################
668  histname = TString::Format("%s/histo_Rjt_%d", groupname.Data(), cent);
669  histtitle = TString::Format("%s ;dR;j_{T} (GeV/c);", histname.Data()); // j_T track vs dR
670  fHistManager.CreateTH2(histname, histtitle, h_Rjt_xnbin, h_Rjt_xlow, h_Rjt_xhigh, h_Rjt_ynbin, h_Rjt_ylow, h_Rjt_yhigh);
671 
672  //########################################################
673  histname = TString::Format("%s/histo_Rjt_n80_%d", groupname.Data(), cent);
674  histtitle = TString::Format("%s ;dR;j_{T} (GeV/c);", histname.Data()); // j_T track vs dR
675  fHistManager.CreateTH2(histname, histtitle, h_Rjt_xnbin, h_Rjt_xlow, h_Rjt_xhigh, h_Rjt_ynbin, h_Rjt_ylow, h_Rjt_yhigh);
676 
677  histname = TString::Format("%s/histo_Rjt_n90_%d", groupname.Data(), cent);
678  histtitle = TString::Format("%s ;dR;j_{T} (GeV/c);", histname.Data()); // j_T track vs dR
679  fHistManager.CreateTH2(histname, histtitle, h_Rjt_xnbin, h_Rjt_xlow, h_Rjt_xhigh, h_Rjt_ynbin, h_Rjt_ylow, h_Rjt_yhigh);
680 
681  //########################################################
682  histname = TString::Format("%s/histo_Rjt_pt80_%d", groupname.Data(), cent);
683  histtitle = TString::Format("%s ;dR;j_{T} (GeV/c);", histname.Data()); // j_T track vs dR
684  fHistManager.CreateTH2(histname, histtitle, h_Rjt_xnbin, h_Rjt_xlow, h_Rjt_xhigh, h_Rjt_ynbin, h_Rjt_ylow, h_Rjt_yhigh);
685 
686  histname = TString::Format("%s/histo_Rjt_pt90_%d", groupname.Data(), cent);
687  histtitle = TString::Format("%s ;dR;j_{T} (GeV/c);", histname.Data()); // j_T track vs dR
688  fHistManager.CreateTH2(histname, histtitle, h_Rjt_xnbin, h_Rjt_xlow, h_Rjt_xhigh, h_Rjt_ynbin, h_Rjt_ylow, h_Rjt_yhigh);
689  //########################################################
690 
691  //=====================================================================================
692  Int_t h_jt_xnbin = 500; Double_t h_jt_xbinwidth = 0.01; Double_t h_jt_xlow = 0.;
693  Double_t h_jt_xhigh = h_jt_xlow + h_jt_xbinwidth * h_jt_xnbin;
694 
695  //########################################################
696  histname = TString::Format("%s/histo_jt_%d", groupname.Data(), cent);
697  histtitle = TString::Format("%s ;j_{T} (GeV/c);1/N_{jets} dN/dj_{T};", histname.Data()); // j_T track vs dR
698  fHistManager.CreateTH1(histname, histtitle, h_jt_xnbin, h_jt_xlow, h_jt_xhigh);
699 
700  //########################################################
701  histname = TString::Format("%s/histo_jt_n80_%d", groupname.Data(), cent);
702  histtitle = TString::Format("%s ;j_{T} (GeV/c);1/N_{jets} dN/dj_{T};", histname.Data()); // j_T track vs dR
703  fHistManager.CreateTH1(histname, histtitle, h_jt_xnbin, h_jt_xlow, h_jt_xhigh);
704 
705  histname = TString::Format("%s/histo_jt_n90_%d", groupname.Data(), cent);
706  histtitle = TString::Format("%s ;j_{T} (GeV/c);1/N_{jets} dN/dj_{T};", histname.Data()); // j_T track vs dR
707  fHistManager.CreateTH1(histname, histtitle, h_jt_xnbin, h_jt_xlow, h_jt_xhigh);
708 
709  //########################################################
710  histname = TString::Format("%s/histo_jt_pt80_%d", groupname.Data(), cent);
711  histtitle = TString::Format("%s ;j_{T} (GeV/c);1/N_{jets} dN/dj_{T};", histname.Data()); // j_T track vs dR
712  fHistManager.CreateTH1(histname, histtitle, h_jt_xnbin, h_jt_xlow, h_jt_xhigh);
713 
714  histname = TString::Format("%s/histo_jt_pt90_%d", groupname.Data(), cent);
715  histtitle = TString::Format("%s ;j_{T} (GeV/c);1/N_{jets} dN/dj_{T};", histname.Data()); // j_T track vs dR
716  fHistManager.CreateTH1(histname, histtitle, h_jt_xnbin, h_jt_xlow, h_jt_xhigh);
717  //########################################################
718 
719  } //end of loop over fNcentBins
720  }// end of loop over jet containers
721 
722 
723  // =========== Switch on Sumw2 for all histos ===========
724  TH1::SetDefaultSumw2(kTRUE);
725  TH2::SetDefaultSumw2(kTRUE);
726 
727  // add all fHistManager content to fOutput
728  TIter nexthist(fHistManager.GetListOfHistograms());
729  TObject* obj = NULL;
730  while ((obj = nexthist())) { fOutput->Add(obj); }
731 
732  PostData ( 1, fOutput ); // Post data for ALL output slots > 0 here.
733  }
734 
735 /*
736  * This function is executed automatically for the first event.
737  * Some extra initialization can be performed here.
738  */
741  }
742 
743 /*
744  * This function is called once at the end of the analysis.
745  */
747  }
748 
749 //________________________________________________________________________
751  return fHistManager.FindObject(histName);
752 }
753 
754 // This function (overloading the base class) uses AliEventCuts to perform event selection.
755 //________________________________________________________________________
757  if (fUseAliEventCuts) {
758  if (!fEventCuts.AcceptEvent(InputEvent())) {
759  PostData(1, fOutput);
760  return kFALSE;
761  }
762  }
763  else {
765  return answer;
766  }
767  return kTRUE;
768 }
769 
770 
771 //########################################################################
772 // Namespace AliAnalysisTaskEmcalJetCDF
773 //########################################################################
774 
775 //__________________________________________________________________________________________________
777  {
778  // Sorting by p_T (decreasing) event tracks
779  Int_t entries = event->GetNumberOfTracks();
780 
781  // Create vector for Pt sorting
782  std::vector<ptidx_pair> pair_list;
783  pair_list.reserve ( entries );
784 
785  for ( Int_t i_entry = 0; i_entry < entries; i_entry++ )
786  {
787  AliVParticle* track = event->GetTrack ( i_entry );
788  if ( !track ) { std::cout << Form ("Unable to find track %d in collection %s", i_entry, event->GetName()) << std::endl ; continue; }
789 
790  pair_list.push_back ( std::make_pair ( track->Pt(), i_entry ) );
791  }
792 
793  std::stable_sort ( pair_list.begin(), pair_list.end(), sort_descend() );
794 
795  // return an vector of indexes of constituents (sorted descending by pt)
796  std::vector<Int_t> index_sorted_list;
797  index_sorted_list.reserve ( entries );
798 
799  // populating the return object with indexes of sorted tracks
800  for (auto it : pair_list) { index_sorted_list.push_back(it.second); }
801 
802  return index_sorted_list;
803  }
804 
805 //__________________________________________________________________________________________________
807  {
808  // Sorting by p_T (decreasing) event tracks
809  Int_t entries = trackscont->GetNEntries();
810 
811  // Create vector for Pt sorting
812  std::vector<ptidx_pair> pair_list;
813  pair_list.reserve ( entries );
814 
815  UInt_t i_entry = 0;
816  AliParticleContainer* partCont = 0;
817  for(auto part : partCont->all()) {
818  if (!part) {continue;}
819  pair_list.push_back ( std::make_pair ( part->Pt(), i_entry++ ) );
820  }
821 
822  std::stable_sort ( pair_list.begin(), pair_list.end(), sort_descend() );
823 
824  // return an vector of indexes of constituents (sorted descending by pt)
825  std::vector<Int_t> index_sorted_list;
826  index_sorted_list.reserve ( i_entry );
827 
828  // populating the return object with indexes of sorted tracks
829  for (auto it : pair_list) { index_sorted_list.push_back(it.second); }
830 
831  return index_sorted_list;
832  }
833 
834 //__________________________________________________________________________________________________
835 AliAnalysisTaskEmcalJetCDF* PWGJE::EMCALJetTasks::AliAnalysisTaskEmcalJetCDF_NS::AddTaskEmcalJetCDF ( const char* ntracks, const char* nclusters, const char* ncells, const char* ntracks_mc, const char* tag)
836  {
837  // Get the pointer to the existing analysis manager via the static access method.
838  //==============================================================================
839  AliAnalysisManager* mgr = AliAnalysisManager::GetAnalysisManager();
840  if ( !mgr ) { ::Error ( "AddTaskEmcalJetCDF", "No analysis manager to connect to." ); return NULL; }
841 
842  // Check the analysis type using the event handlers connected to the analysis manager.
843  //==============================================================================
844  AliVEventHandler* handler = mgr->GetInputEventHandler();
845  if (!handler) { ::Error ( "AddTaskEmcalJetCDF", "This task requires an input event handler" ); return NULL; }
846 
847  enum EDataType_t { kUnknown, kESD, kAOD }; EDataType_t dataType = kUnknown;
848 
849  if (handler->InheritsFrom("AliESDInputHandler")) { dataType = kESD; }
850  else
851  if (handler->InheritsFrom("AliAODInputHandler")) { dataType = kAOD; }
852 
853  //-------------------------------------------------------
854  // Init the task and do settings
855  //-------------------------------------------------------
856  TString suffix (tag);
857  TString tracks (ntracks);
858  TString clusters (nclusters);
859  TString tracks_mc (ntracks_mc);
860  TString cells (ncells);
861 
862  if ( tracks.EqualTo("usedefault") ) {
863  if ( dataType == kESD ) { tracks = "Tracks"; }
864  else
865  if ( dataType == kAOD ) { tracks = "tracks"; }
866  else
867  { tracks = ""; }
868  }
869 
870  if ( clusters.EqualTo("usedefault") ) {
871  if ( dataType == kESD ) { clusters = "CaloClusters"; }
872  else
873  if ( dataType == kAOD ) { clusters = "caloClusters"; }
874  else
875  { clusters = ""; }
876  }
877 
878  if ( cells.EqualTo("usedefault") ) {
879  if (dataType == kESD) { cells = "EMCALCells"; }
880  else
881  if (dataType == kAOD) { cells = "emcalCells"; }
882  else
883  { cells = ""; }
884  }
885 
886  TString name("JetCDF");
887  if (!tracks.IsNull()) { name += "_" + tracks; }
888  if (!clusters.IsNull()) { name += "_" + clusters; }
889  if (!cells.IsNull()) { name += "_" + cells; }
890  if (!suffix.IsNull()) { name += "_" + suffix; }
891 
892  AliAnalysisTaskEmcalJetCDF* cdfTask = new AliAnalysisTaskEmcalJetCDF ( name.Data() );
893  cdfTask->SetVzRange(-10,10);
894  cdfTask->SetCaloCellsName(cells.Data());
895 
896  if ( tracks.EqualTo("mcparticles") ) {
897  cdfTask->AddMCParticleContainer(tracks.Data()); }
898  else if ( tracks.EqualTo("tracks") || tracks.EqualTo("Tracks") ) {
899  cdfTask->AddTrackContainer(tracks.Data()); }
900  else if ( !tracks.IsNull()) {
901  cdfTask->AddParticleContainer(tracks.Data()); }
902 
903  // Add the generator-level container, if specified
904  if ( !tracks_mc.IsNull() ) {
905  AliMCParticleContainer* mcpartCont = cdfTask->AddMCParticleContainer(tracks_mc.Data());
906  mcpartCont->SelectPhysicalPrimaries(kTRUE);
907  }
908 
909  if (!clusters.IsNull()) { cdfTask->AddClusterContainer(clusters.Data()); }
910 
911  //-------------------------------------------------------
912  // Final settings, pass to manager and set the containers
913  //-------------------------------------------------------
914  mgr->AddTask ( cdfTask );
915 
916  // Create containers for input/output
917  AliAnalysisDataContainer *cinput1 = mgr->GetCommonInputContainer();
918 
919  TString contname = name + "_histos";
920  TString outfile (Form("%s", AliAnalysisManager::GetCommonFileName()));
921  AliAnalysisDataContainer *coutput1 = mgr->CreateContainer ( contname.Data(), TList::Class(), AliAnalysisManager::kOutputContainer, outfile.Data() );
922 
923  mgr->ConnectInput ( cdfTask, 0, cinput1 );
924  mgr->ConnectOutput ( cdfTask, 1, coutput1 );
925 
926  return cdfTask;
927  }
928 
929 //__________________________________________________________________________________________________
930  AliJetContainer* PWGJE::EMCALJetTasks::AliAnalysisTaskEmcalJetCDF_NS::jetContSetParams ( AliJetContainer* jetCont, Float_t jetptmin, Float_t jetptmax, Float_t jetareacutperc, Int_t leadhadtype, Int_t nLeadJets, Float_t mintrackpt, Float_t maxtrackpt)
931  {
932  if (!jetCont) { return NULL; }
933  jetCont->SetJetPtCut ( jetptmin );
934  jetCont->SetJetPtCutMax ( jetptmax );
935  jetCont->SetPercAreaCut ( jetareacutperc );
936  jetCont->SetLeadingHadronType ( leadhadtype ); // 0 = charged, 1 = neutral, 2 = both
937  jetCont->SetNLeadingJets(nLeadJets);
938  jetCont->SetMinTrackPt(mintrackpt);
939  jetCont->SetMaxTrackPt(maxtrackpt);
940  return jetCont;
941  }
942 
943 //__________________________________________________________________________________________________
944 TChain* PWGJE::EMCALJetTasks::AliAnalysisTaskEmcalJetCDF_NS::CreateChain ( const char* filelist, const char* cTreeNameArg, const char* friends, UInt_t iNumFiles, UInt_t iStartWithFile ) {
945 TString sTreeNameArg (cTreeNameArg), treeName;
946 
947 TFileCollection filecoll ("anachain","File collection for analysis"); // easy manipulation of file collections
948 Int_t iAddedFiles = filecoll.AddFromFile(filelist,iNumFiles,iStartWithFile);
949 if ( iAddedFiles < 1 ) { std::cout << "NO Files added to collection !!!" << std::endl; return NULL; }
950 
951 // if cTreeNameArg is auto lets try to autodetect what type of tree we have;
952 // the assuption is that all files are the same and the first one is reprezentative
953 THashList* list = filecoll.GetList();
954 if ( sTreeNameArg.EqualTo("auto") ) { // if tree name is not specified
955  TRegexp tree_regex ("[aod,esd]Tree");
956  TFileInfo* fileinfo = dynamic_cast<TFileInfo*>(list->At(0)); // get first fileinfo in list
957  TFile file (fileinfo->GetFirstUrl()->GetFile()); // get the actual TFile
958  if (file.IsZombie()) { cout << "Should not reach this message!! Error opening file" << endl; return NULL; }
959 
960  // lets parse the TFile
961  TIter next(file.GetListOfKeys());
962  TKey* key = NULL;
963  while (( key = dynamic_cast<TKey*>(next()) )) {
964  TString class_name = key->GetClassName();
965  if ( ! class_name.EqualTo("TTree") ) { continue; } // searching for first TTree
966 
967  TString key_name = key->GetName();
968  if ( key_name.Contains(tree_regex) ) { treeName = key_name; break;} // that is named either aodTree or esdTree
969  }
970  file.Close();
971  }
972 else
973  { treeName = sTreeNameArg ; } // tree name is specified
974 
975 TChain* chain = new TChain (treeName.Data(),""); // lets create the chain
976 if ( chain->AddFileInfoList(list) == 0 ) { return NULL; } // and add file collection (THashList is a Tlist that is a TCollection)
977 
978 // start treatment of friends
979 TChain* chainFriend = NULL;
980 TString sFriends (friends);
981 if (!sFriends.IsNull()) {
982  TString friends_treename, friends_filename;
983  TObjArray* arr = sFriends.Tokenize("/");
984  TObjString* strobj_file = dynamic_cast<TObjString*>(arr->At(0));
985  if (strobj_file) { friends_filename = strobj_file->GetString(); }
986  TObjString* strobj_tree = dynamic_cast<TObjString*>(arr->At(1));
987  if (strobj_tree) { friends_treename = strobj_tree->GetString(); }
988  delete arr;
989 
990  if (friends_treename.IsNull()) {
991  if (treeName.EqualTo("esdTree")) {
992  friends_treename = "esdFriendTree"; }
993  else if (treeName.EqualTo("aodTree")) {
994  friends_treename = "aodTree"; }
995  else {
996  cout << "friends argument specified but tree name neither specified nor auto-detected (unknown tree name to associate with known friend tree name)";
997  return chain; // stop processing of friends, just return the chain created so far
998  }
999  }
1000 
1001  chainFriend = new TChain(friends_treename.Data());
1002  TString friendinfo_for_chain = "/" + friends_filename + "/" + friends_treename;
1003 
1004  TIter next(list);
1005  TFileInfo* fileinfo = NULL;
1006  while (( fileinfo = dynamic_cast<TFileInfo*>(next()) )) {
1007  TString dirname = gSystem->DirName(fileinfo->GetFirstUrl()->GetFile());
1008  TString friend_for_chain = dirname + friendinfo_for_chain;
1009  if (chainFriend) { chainFriend->Add(friend_for_chain.Data()); }
1010  }
1011 
1012  if (chainFriend) { chain->AddFriend(chainFriend); }
1013  }
1014 
1015 return chain;
1016 }
1017 
1018 // kate: indent-mode none; indent-width 2; replace-tabs on;
1019 
THistManager fHistManager
Histogram manager.
THashList * CreateHistoGroup(const char *groupname)
Create a new group of histograms within a parent group.
Bool_t fUseAliEventCuts
Flag to use AliEventCuts (otherwise AliAnalysisTaskEmcal will be used)
double Double_t
Definition: External.C:58
Bool_t IsEventSelected()
Performing event selection.
const AliParticleIterableContainer all() const
void SetLeadingHadronType(Int_t t)
TSystem * gSystem
Declaration of class AliTLorentzVector.
void ExecOnce()
Perform steps needed to initialize the analysis.
AliClusterContainer * GetClusterContainer() const
Int_t fCentBin
!event centrality bin
void SetPercAreaCut(Float_t p)
void SetVzRange(Double_t min, Double_t max)
Set pre-configured event cut object.
AliClusterContainer * AddClusterContainer(const char *n)
Create new cluster container and attach it to the task.
Container for particles within the EMCAL framework.
AliParticleContainer * GetParticleContainer() const
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)
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.
Definition: External.C:228
Int_t fNcentBins
how many centrality bins
Definition: External.C:212
Int_t GetNAcceptedClusters() const
TClonesArray * GetArray() const
Double_t Z_pt(const AliEmcalJet *jet, const AliVParticle *trk)
AliEventCuts fEventCuts
event selection utility
AliMCParticleContainer * AddMCParticleContainer(const char *n)
Create new container for MC particles and attach it to the task.
TObjArray fJetCollArray
jet collection array
functional for sorting pair by first element - descending
Double_t Perp(const AliVParticle &trk1, const AliVParticle &trk2)
AliJetContainer * 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.)
const char * GetName() const
virtual Bool_t IsEventSelected()
Performing event selection.
AliEmcalList * fOutput
!output list
virtual ~AliAnalysisTaskEmcalJetCDF()
Destructor.
AliTrackContainer * AddTrackContainer(const char *n)
Create new track container and attach it to the task.
AliAnalysisTaskEmcalJetCDF * AddTaskEmcalJetCDF(const char *ntracks="usedefault", const char *nclusters="usedefault", const char *ncells="usedefault", const char *ntracks_mc="", const char *tag="CDF")
void SelectPhysicalPrimaries(Bool_t s)
TFile * file
TList with histograms for a given trigger.
void SetMakeGeneralHistograms(Bool_t g)
Enable general histograms.
Base task in the EMCAL jet framework.
unsigned short UShort_t
Definition: External.C:28
TChain * CreateChain(const char *filelist="filelist.txt", const char *cTreeNameArg="auto", const char *friends="", UInt_t iNumFiles=-1, UInt_t iStartWithFile=1)
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.
Double_t GetMaxPt() const
Container for MC-true particles within the EMCAL framework.
Container for jet within the EMCAL jet framework.
void SetJetPtCutMax(Float_t cut)
Analysis of jet shapes and FF of all jets and leading jets.