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