AliPhysics  cf1a5e2 (cf1a5e2)
AliAnalysisTaskIDFragmentationFunction.cxx
Go to the documentation of this file.
1 // *************************************************************************
2 // * *
3 // * Task for ID Fragmentation Function Analysis *
4 // * *
5 // *************************************************************************
6 
7 
8 /**************************************************************************
9  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
10  * *
11  * Author: The ALICE Off-line Project. *
12  * Contributors are mentioned in the code where appropriate. *
13  * *
14  * Permission to use, copy, modify and distribute this software and its *
15  * documentation strictly for non-commercial purposes is hereby granted *
16  * without fee, provided that the above copyright notice appears in all *
17  * copies and that both the copyright notice and this permission notice *
18  * appear in the supporting documentation. The authors make no claims *
19  * about the suitability of this software for any purpose. It is *
20  * provided "as is" without express or implied warranty. *
21  **************************************************************************/
22 
23 /* $Id: */
24 
25 #include "TList.h"
26 #include "TH1F.h"
27 #include "TH2F.h"
28 #include "TH3F.h"
29 #include "TString.h"
30 #include "THnSparse.h"
31 #include "TProfile.h"
32 #include "TFile.h"
33 #include "TKey.h"
34 #include "TRandom3.h"
35 #include "TAxis.h"
36 #include "TVector3.h"
37 
38 #include "AliAODInputHandler.h"
39 #include "AliAODHandler.h"
40 #include "AliESDEvent.h"
41 #include "AliAODMCParticle.h"
42 #include "AliEmcalJet.h"
43 #include "AliAODJetEventBackground.h"
44 #include "AliGenPythiaEventHeader.h"
45 #include "AliGenHijingEventHeader.h"
46 #include "AliInputEventHandler.h"
47 #include "AliJetContainer.h"
48 
50 #include "AliAnalysisManager.h"
51 #include "AliAnalysisTaskSE.h"
52 #include "AliAnalysisUtils.h"
53 #include "AliVParticle.h"
54 #include "AliVEvent.h"
55 
56 #include "AliAnalysisTaskPID.h"
57 #include "AliPIDResponse.h"
58 
59 
61 using std::cout;
62 using std::endl;
63 using std::cerr;
64 
66 
67 //____________________________________________________________________________
70  ,fESD(0)
71  ,fAOD(0)
72  ,fAODJets(0)
73  ,fAODExtension(0)
74  ,fNonStdFile("")
75  ,fCentralityEstimator("V0M")
76  ,fNameTrackContainer("Tracks")
77  ,fNameTrackContainerEfficiency("")
78  ,fNameMCParticleContainer("")
79  ,fNameJetContainer("")
80  ,fNameMCParticleJetContainer("")
81  ,fUseAODInputJets(kTRUE)
82  ,fUsePhysicsSelection(kTRUE)
83  ,fEvtSelectionMask(0)
84  ,fEventClass(0)
85  ,fMaxVertexZ(10)
86  ,fFFRadius(0)
87  ,fFFMinLTrackPt(-1)
88  ,fFFMaxTrackPt(-1)
89  ,fFFMinnTracks(0)
90  ,fAvgTrials(0)
91  ,fStoreXi(0)
92  ,fStudyTransversalJetStructure(kFALSE)
93  ,fCommonHistList(0)
94  ,fh1EvtSelection(0)
95  ,fh1VtxSelection(0)
96  ,fh1VertexNContributors(0)
97  ,fh1VertexZ(0)
98  ,fh1EvtMult(0)
99  ,fh1EvtCent(0)
100  ,fh1Xsec(0)
101  ,fh1Trials(0)
102  ,fh1PtHard(0)
103  ,fh1PtHardTrials(0)
104  ,fh1EvtsPtHardCut(0)
105  ,fh1nRecJetsCuts(0)
106  ,fh1nGenJets(0)
107  ,fhDCA_XY(0)
108  ,fhDCA_Z(0)
109  ,fhJetPtRefMultEta5(0)
110  ,fhJetPtRefMultEta8(0)
111  ,fhJetPtMultPercent(0)
112 
113  ,fRandom(0)
114 
115  ,fOnlyLeadingJets(kFALSE)
116  ,fMCPtHardCut(-1.)
117 
118  ,fAnaUtils(0)
119 
120  // PID framework
121  ,fNumInclusivePIDtasks(0)
122  ,fNumJetPIDtasks(0)
123  ,fNumJetUEPIDtasks(0)
124  ,fNameInclusivePIDtask(0x0)
125  ,fNameJetPIDtask(0x0)
126  ,fNameJetUEPIDtask(0x0)
127  ,fInclusivePIDtask(0x0)
128  ,fJetPIDtask(0x0)
129  ,fJetUEPIDtask(0x0)
130  ,fUseInclusivePIDtask(kFALSE)
131  ,fUseJetPIDtask(kFALSE)
132  ,fUseJetUEPIDtask(kFALSE)
133  ,fIsPP(kFALSE)
134  ,fFillDCA(kFALSE)
135  ,fRCTrials(1)
136  ,fUEMethods(0x0)
137 {
138  // default constructor
139 
140  if (fFillDCA) {
141  for (Int_t i = 0; i < AliPID::kSPECIES; i++) {
142 
143  fhDCA_XY_prim_MCID[i] = 0x0;
144  fhDCA_Z_prim_MCID[i] = 0x0;
145 
146  fhDCA_XY_sec_MCID[i] = 0x0;
147  fhDCA_Z_sec_MCID[i] = 0x0;
148  }
149  }
150 }
151 
152 //_______________________________________________________________________________________________
154  : AliAnalysisTaskEmcalJet(name, kTRUE)
155  ,fESD(0)
156  ,fAOD(0)
157  ,fAODJets(0)
158  ,fAODExtension(0)
159  ,fNonStdFile("")
160  ,fCentralityEstimator("V0M")
161  ,fNameTrackContainer("Tracks")
162  ,fNameTrackContainerEfficiency("")
163  ,fNameMCParticleContainer("")
164  ,fNameJetContainer("")
165  ,fNameMCParticleJetContainer("")
166  ,fUseAODInputJets(kTRUE)
167  ,fUsePhysicsSelection(kTRUE)
168  ,fEvtSelectionMask(0)
169  ,fEventClass(0)
170  ,fMaxVertexZ(10)
171  ,fFFRadius(0)
172  ,fFFMinLTrackPt(-1)
173  ,fFFMaxTrackPt(-1)
174  ,fFFMinnTracks(0)
175  ,fAvgTrials(0)
176  ,fStoreXi(0)
177  ,fStudyTransversalJetStructure(kFALSE)
178  ,fCommonHistList(0)
179  ,fh1EvtSelection(0)
180  ,fh1VtxSelection(0)
181  ,fh1VertexNContributors(0)
182  ,fh1VertexZ(0)
183  ,fh1EvtMult(0)
184  ,fh1EvtCent(0)
185  ,fh1Xsec(0)
186  ,fh1Trials(0)
187  ,fh1PtHard(0)
188  ,fh1PtHardTrials(0)
189  ,fh1EvtsPtHardCut(0)
190  ,fh1nRecJetsCuts(0)
191  ,fh1nGenJets(0)
192  ,fhDCA_XY(0)
193  ,fhDCA_Z(0)
194  ,fhJetPtRefMultEta5(0)
195  ,fhJetPtRefMultEta8(0)
196  ,fhJetPtMultPercent(0)
197  ,fRandom(0)
198  ,fOnlyLeadingJets(kFALSE)
199  ,fMCPtHardCut(-1.)
200  ,fAnaUtils(0)
201  // PID framework
202  ,fNumInclusivePIDtasks(0)
203  ,fNumJetPIDtasks(0)
204  ,fNumJetUEPIDtasks(0)
205  ,fNameInclusivePIDtask(0x0)
206  ,fNameJetPIDtask(0x0)
207  ,fNameJetUEPIDtask(0x0)
208  ,fInclusivePIDtask(0x0)
209  ,fJetPIDtask(0x0)
210  ,fJetUEPIDtask(0x0)
211  ,fUseInclusivePIDtask(kFALSE)
212  ,fUseJetPIDtask(kFALSE)
213  ,fUseJetUEPIDtask(kFALSE)
214  ,fIsPP(kFALSE)
215  ,fFillDCA(kFALSE)
216  ,fRCTrials(1)
217  ,fUEMethods(0x0)
218 {
219  // constructor
220 
221  if (fFillDCA) {
222  for (Int_t i = 0; i < AliPID::kSPECIES; i++) {
223 
224  fhDCA_XY_prim_MCID[i] = 0x0;
225  fhDCA_Z_prim_MCID[i] = 0x0;
226 
227  fhDCA_XY_sec_MCID[i] = 0x0;
228  fhDCA_Z_sec_MCID[i] = 0x0;
229  }
230  }
231 
232  DefineOutput(1,TList::Class());
233 }
234 
235 //___________________________________________________________________________
237 {
238  // destructor
239 
240  if(fRandom) delete fRandom;
241 
242  delete [] fNameInclusivePIDtask;
243  fNameInclusivePIDtask = 0x0;
244 
245  delete [] fInclusivePIDtask;
246  fInclusivePIDtask = 0x0;
247 
248  delete [] fNameJetPIDtask;
249  fNameJetPIDtask = 0x0;
250 
251  delete [] fJetPIDtask;
252  fJetPIDtask = 0x0;
253 
254  delete [] fNameJetUEPIDtask;
255  fNameJetUEPIDtask = 0x0;
256 
257  delete [] fJetUEPIDtask;
258  fJetUEPIDtask = 0x0;
259 
260  delete fAnaUtils;
261  fAnaUtils = 0x0;
262 }
263 
264 //_________________________________________________________________________________
266 {
267  //
268  // Implemented Notify() to read the cross sections
269  // and number of trials from pyxsec.root
270  // (taken from AliAnalysisTaskJetSpectrum2)
271  //
272  TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
273  Float_t xsection = 0;
274  Float_t ftrials = 1;
275 
276  fAvgTrials = 1;
277  if(tree){
278  TFile *curfile = tree->GetCurrentFile();
279  if (!curfile) {
280  Error("Notify","No current file");
281  return kFALSE;
282  }
283 
284  AliAnalysisHelperJetTasks::PythiaInfoFromFile(curfile->GetName(),xsection,ftrials);
285 
286  if (fUseInclusivePIDtask) {
287  for (Int_t i = 0; i < fNumInclusivePIDtasks; i++)
288  fInclusivePIDtask[i]->FillXsec(xsection);
289  }
290 
291  if (fUseJetPIDtask) {
292  for (Int_t i = 0; i < fNumJetPIDtasks; i++)
293  fJetPIDtask[i]->FillXsec(xsection);
294  }
295 
296  if (fUseJetUEPIDtask) {
297  for (Int_t i = 0; i < fNumJetUEPIDtasks; i++)
298  fJetUEPIDtask[i]->FillXsec(xsection);
299  }
300 
301  if(!fh1Xsec||!fh1Trials){
302  Printf("%s:%d No Histogram fh1Xsec",(char*)__FILE__,__LINE__);
303  return kFALSE;
304  }
305 
306  fh1Xsec->Fill("<#sigma>",xsection);
307  // construct a poor man average trials
308  Float_t nEntries = (Float_t)tree->GetTree()->GetEntries();
309  if(ftrials>=nEntries && nEntries>0.)fAvgTrials = ftrials/nEntries;
310  }
311 
312  // Set seed for backg study
313  delete fRandom;
314  fRandom = new TRandom3();
315  fRandom->SetSeed(0);
316 
317  return kTRUE;
318 }
319 
320 //__________________________________________________________________
322 {
323  // create output objects
324 
325  if(fDebug > 1) Printf("AliAnalysisTaskIDFragmentationFunction::UserCreateOutputObjects()");
326 
327  //
328  // Create histograms / output container
329  //
330 
331  OpenFile(1);
332  fCommonHistList = new TList();
333  fCommonHistList->SetOwner(kTRUE);
334 
335  Bool_t oldStatus = TH1::AddDirectoryStatus();
336  TH1::AddDirectory(kFALSE);
337 
338 
339  // Histograms
340  fh1EvtSelection = new TH1F("fh1EvtSelection", "Event Selection", 6, -0.5, 5.5);
341  fh1EvtSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
342  fh1EvtSelection->GetXaxis()->SetBinLabel(2,"event selection: rejected");
343  fh1EvtSelection->GetXaxis()->SetBinLabel(3,"event class: rejected");
344  fh1EvtSelection->GetXaxis()->SetBinLabel(4,"vertex Ncontr: rejected");
345  fh1EvtSelection->GetXaxis()->SetBinLabel(5,"vertex z: rejected");
346  fh1EvtSelection->GetXaxis()->SetBinLabel(6,"vertex type: rejected");
347 
348  fh1VtxSelection = new TH1F("fh1VtxSelection", "Vertex Selection", 10, -1, 9);
349  fh1VtxSelection->GetXaxis()->SetBinLabel(1,"Undef");
350  fh1VtxSelection->GetXaxis()->SetBinLabel(2,"Primary");
351  fh1VtxSelection->GetXaxis()->SetBinLabel(3,"Kink");
352  fh1VtxSelection->GetXaxis()->SetBinLabel(4,"V0");
353  fh1VtxSelection->GetXaxis()->SetBinLabel(5,"Cascade");
354  fh1VtxSelection->GetXaxis()->SetBinLabel(6,"Multi");
355  fh1VtxSelection->GetXaxis()->SetBinLabel(7,"SPD");
356  fh1VtxSelection->GetXaxis()->SetBinLabel(8,"PileUpSPD");
357  fh1VtxSelection->GetXaxis()->SetBinLabel(9,"PileUpTracks");
358  fh1VtxSelection->GetXaxis()->SetBinLabel(10,"TPC");
359 
360  fh1VertexNContributors = new TH1F("fh1VertexNContributors", "Vertex N contributors", 2500,-.5, 2499.5);
361  fh1VertexZ = new TH1F("fh1VertexZ", "Vertex z distribution", 30, -15., 15.);
362  fh1EvtMult = new TH1F("fh1EvtMult","Event multiplicity, track pT cut > 150 MeV/c, |#eta| < 0.9",120,0.,12000.);
363  fh1EvtCent = new TH1F("fh1EvtCent","centrality",100,0.,100.);
364 
365  fh1Xsec = new TProfile("fh1Xsec","xsec from pyxsec.root",1,0,1);
366  fh1Xsec->GetXaxis()->SetBinLabel(1,"<#sigma>");
367  fh1Trials = new TH1F("fh1Trials","trials from pyxsec.root",1,0,1);
368  fh1Trials->GetXaxis()->SetBinLabel(1,"#sum{ntrials}");
369  fh1PtHard = new TH1F("fh1PtHard","PYTHIA Pt hard;p_{T,hard}",350,-.5,349.5);
370  fh1PtHardTrials = new TH1F("fh1PtHardTrials","PYTHIA Pt hard weight with trials;p_{T,hard}",350,-.5,349.5);
371 
372  fh1EvtsPtHardCut = new TH1F("fh1EvtsPtHardCut", "#events before and after MC #it{p}_{T,hard} cut;;Events",2,0,2);
373  fh1EvtsPtHardCut->GetXaxis()->SetBinLabel(1, "All");
374  fh1EvtsPtHardCut->GetXaxis()->SetBinLabel(2, "#it{p}_{T,hard}");
375 
376 
377  fh1nRecJetsCuts = new TH1F("fh1nRecJetsCuts","reconstructed jets per event",10,-0.5,9.5);
378  fh1nGenJets = new TH1F("fh1nGenJets","generated jets per event",10,-0.5,9.5);
379 
380 
381  // ____________ define histograms ___________________
382 
390  fCommonHistList->Add(fh1Xsec);
395 
397 
398  // Default analysis utils
399  fAnaUtils = new AliAnalysisUtils();
400 
401  // Not used yet, but to be save, forward vertex z cut to analysis utils object
402  fAnaUtils->SetMaxVtxZ(fMaxVertexZ);
403 
404  // Load PID framework if desired
405  if(fDebug > 1) Printf("AliAnalysisTaskIDFragmentationFunction::UserCreateOutputObjects() -> Loading PID framework");
406 
408  TObjArray* tasks = AliAnalysisManager::GetAnalysisManager()->GetTasks();
409  if (!tasks) {
410  Printf("ERROR loading PID tasks: Failed to retrieve tasks from analysis manager!\n");
411 
412  fUseInclusivePIDtask = kFALSE;
413  fUseJetPIDtask = kFALSE;
414  fUseJetUEPIDtask = kFALSE;
415  }
416 
417  if (fUseInclusivePIDtask) {
418  delete [] fInclusivePIDtask;
419  fInclusivePIDtask = 0x0;
420 
421  if (fNumInclusivePIDtasks > 0) {
423 
424  for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
425  fInclusivePIDtask[i] = (AliAnalysisTaskPID*)tasks->FindObject(fNameInclusivePIDtask[i].Data());
426 
427  if (!fInclusivePIDtask[i]) {
428  Printf("ERROR: Failed to load inclusive pid task!\n");
429  fUseInclusivePIDtask = kFALSE;
430  }
431  }
432  }
433  else {
434  Printf("WARNING: zero inclusive pid tasks!\n");
435  fUseInclusivePIDtask = kFALSE;
436  }
437  }
438 
439  if (fUseJetPIDtask) {
440  delete [] fJetPIDtask;
441  fJetPIDtask = 0x0;
442 
443  if (fNumJetPIDtasks > 0) {
445 
446  for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
447  fJetPIDtask[i] = (AliAnalysisTaskPID*)tasks->FindObject(fNameJetPIDtask[i].Data());
448 
449  if (!fJetPIDtask[i]) {
450  Printf("ERROR: Failed to load jet pid task!\n");
451  fUseJetPIDtask = kFALSE;
452  }
453  }
454  }
455  else {
456  Printf("WARNING: zero jet pid tasks!\n");
457  fUseJetPIDtask = kFALSE;
458  }
459  }
460 
461  if (fUseJetUEPIDtask) {
462  delete [] fJetUEPIDtask;
463  fJetUEPIDtask = 0x0;
464 
465  if (fNumJetUEPIDtasks > 0) {
467 
468  for (Int_t i = 0; i < fNumJetUEPIDtasks; i++) {
469  fJetUEPIDtask[i] = (AliAnalysisTaskPID*)tasks->FindObject(fNameJetUEPIDtask[i].Data());
470 
471  if (!fJetUEPIDtask[i]) {
472  Printf("ERROR: Failed to load jet underlying event pid task!\n");
473  fUseJetUEPIDtask = kFALSE;
474  }
475  }
476  }
477  else {
478  Printf("WARNING: zero jet underlying event pid tasks!\n");
479  fUseJetUEPIDtask = kFALSE;
480  }
481  }
482 
483  }
484 
485  const Int_t nRefMultBins = 100;
486  const Double_t refMultUp = 100.;
487  const Double_t refMultDown = 0.;
488 
489  const Int_t nJetPtBins = 20;
490  const Double_t jetPtUp = 100.;
491  const Double_t jetPtDown = 0.;
492 
493  const Int_t nCentBins = 12;
494  const Double_t binsCentpp[nCentBins+1] = { 0, 0.01, 0.1, 1, 5, 10, 15, 20, 30, 40, 50, 70, 100};
495 
496  fhJetPtRefMultEta5 = new TH2F("fhJetPtRefMultEta5",
497  "Correlation between jet energy and event multiplicity (|#eta| < 0.5);Ref. mult. (|#eta| < 0.5);#it{p}_{T, jet}^{ch} (GeV/#it{c})",
498  nRefMultBins, refMultDown, refMultUp, nJetPtBins, jetPtDown, jetPtUp);
499  fhJetPtRefMultEta8 = new TH2F("fhJetPtRefMultEta8",
500  "Correlation between jet energy and event multiplicity (|#eta| < 0.8);Ref. mult. (|#eta| < 0.8);#it{p}_{T, jet}^{ch} (GeV/#it{c})",
501  nRefMultBins, refMultDown, refMultUp, nJetPtBins, jetPtDown, jetPtUp);
502  fhJetPtMultPercent = new TH2F("fhJetPtMultPercent",
503  "Correlation between jet energy and event multiplicity percentile (V0M);Multiplicity Percentile (V0M);#it{p}_{T, jet}^{ch} (GeV/#it{c})",
504  nCentBins, binsCentpp, nJetPtBins, jetPtDown, jetPtUp);
508 
509  if (fUseJetPIDtask && fFillDCA) {
510  const Int_t nPtBins = 68;
511  Double_t binsPt[nPtBins+1] = {0. , 0.05, 0.1, 0.15, 0.2, 0.25, 0.3, 0.35, 0.4, 0.45,
512  0.5, 0.55, 0.6, 0.65, 0.7, 0.75, 0.8, 0.85, 0.9, 0.95,
513  1.0, 1.1 , 1.2, 1.3 , 1.4, 1.5 , 1.6, 1.7 , 1.8, 1.9 ,
514  2.0, 2.2 , 2.4, 2.6 , 2.8, 3.0 , 3.2, 3.4 , 3.6, 3.8 ,
515  4.0, 4.5 , 5.0, 5.5 , 6.0, 6.5 , 7.0, 8.0 , 9.0, 10.0,
516  11.0, 12.0, 13.0, 14.0, 15.0, 16.0, 18.0, 20.0, 22.0, 24.0,
517  26.0, 28.0, 30.0, 32.0, 34.0, 36.0, 40.0, 45.0, 50.0 };
518 
519  const Int_t DCAbins = 320;
520  const Double_t DCA_Z_max = 3.5;
521  const Double_t DCA_XY_max = 2.5;
522 
523  fhDCA_XY = new TH2F("fhDCA_XY", "All rec. tracks;#it{p}_{T} (GeV/#it{c});DCA_{XY}", nPtBins, binsPt, DCAbins, -DCA_XY_max, DCA_XY_max);
524  fhDCA_Z = new TH2F("fhDCA_Z", "All rec. tracks;#it{p}_{T} (GeV/#it{c});DCA_{Z}", nPtBins, binsPt, DCAbins, -DCA_Z_max, DCA_Z_max);
526  fCommonHistList->Add(fhDCA_Z);
527 
528  for (Int_t i = 0; i < AliPID::kSPECIES; i++) {
529  fhDCA_XY_prim_MCID[i] = new TH2F(Form("fhDCA_XY_prim_MCID_%s", AliPID::ParticleShortName(i)),
530  Form("Rec. %s (prim.);#it{p}_{T} (GeV/#it{c});DCA_{XY}", AliPID::ParticleLatexName(i)),
531  nPtBins, binsPt, DCAbins, -DCA_XY_max, DCA_XY_max);
532  fhDCA_Z_prim_MCID[i] = new TH2F(Form("fhDCA_Z_prim_MCID_%s", AliPID::ParticleShortName(i)),
533  Form("Rec. %s (prim.);#it{p}_{T} (GeV/#it{c});DCA_{Z}", AliPID::ParticleLatexName(i)),
534  nPtBins, binsPt, DCAbins, -DCA_Z_max, DCA_Z_max);
537 
538  fhDCA_XY_sec_MCID[i] = new TH2F(Form("fhDCA_XY_sec_MCID_%s", AliPID::ParticleShortName(i)),
539  Form("Rec. %s (sec.);#it{p}_{T} (GeV/#it{c});DCA_{XY}", AliPID::ParticleLatexName(i)),
540  nPtBins, binsPt, DCAbins, -DCA_XY_max, DCA_XY_max);
541  fhDCA_Z_sec_MCID[i] = new TH2F(Form("fhDCA_Z_sec_MCID_%s", AliPID::ParticleShortName(i)),
542  Form("Rec. %s (sec.);#it{p}_{T} (GeV/#it{c});DCA_{Z}", AliPID::ParticleLatexName(i)),
543  nPtBins, binsPt, DCAbins, -DCA_Z_max, DCA_Z_max);
546  }
547  }
548 
549 
550  // =========== Switch on Sumw2 for all histos ===========
551  for (Int_t i=0; i<fCommonHistList->GetEntries(); ++i){
552  TH1 *h1 = dynamic_cast<TH1*>(fCommonHistList->At(i));
553  if (h1) h1->Sumw2();
554  else{
555  THnSparse *hnSparse = dynamic_cast<THnSparse*>(fCommonHistList->At(i));
556  if(hnSparse) hnSparse->Sumw2();
557  }
558  }
559 
560  TH1::AddDirectory(oldStatus);
561 
562  if(fDebug > 2) Printf("AliAnalysisTaskIDFragmentationFunction::UserCreateOutputObjects() -> Posting Output");
563 
564  PostData(1, fCommonHistList);
565 
566  if(fDebug > 2) Printf("AliAnalysisTaskIDFragmentationFunction::UserCreateOutputObjects() -> Done");
567 }
568 
569 //_______________________________________________
571 {
572  // Initialization
573  if(fDebug > 1) Printf("AliAnalysisTaskIDFragmentationFunction::Init()");
574 
575 }
576 
577 //_____________________________________________________________
579 {
580 
581  if(fDebug > 1) Printf("AliAnalysisTaskIDFragmentationFunction::FillHistograms()");
582 
583 
584  if(fDebug > 1) Printf("Analysis event #%5d", (Int_t) fEntry);
585 
586  fMCEvent = MCEvent();
587  if(!fMCEvent){
588  if(fDebug>3) Printf("%s:%d MCEvent not found in the input", (char*)__FILE__,__LINE__);
589  }
590 
591 
592  // Extract pThard and nTrials in case of MC.
593 
594  Double_t ptHard = 0.;
595  Double_t nTrials = 1; // trials for MC trigger weight for real data
596  Bool_t pythiaGenHeaderFound = kFALSE;
597 
598  if(fMCEvent){
599  AliGenEventHeader* genHeader = fMCEvent->GenEventHeader();
600 
601  if(genHeader){
602  AliGenPythiaEventHeader* pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
603  AliGenHijingEventHeader* hijingGenHeader = 0x0;
604 
605  if(pythiaGenHeader){
606  if(fDebug>3) Printf("%s:%d pythiaGenHeader found", (char*)__FILE__,__LINE__);
607  pythiaGenHeaderFound = kTRUE;
608  nTrials = pythiaGenHeader->Trials();
609  ptHard = pythiaGenHeader->GetPtHard();
610  } else { // no pythia, hijing?
611 
612  if(fDebug>3) Printf("%s:%d no pythiaGenHeader found", (char*)__FILE__,__LINE__);
613 
614  hijingGenHeader = dynamic_cast<AliGenHijingEventHeader*>(genHeader);
615  if(!hijingGenHeader){
616  Printf("%s:%d no pythiaGenHeader or hjingGenHeader found", (char*)__FILE__,__LINE__);
617  } else {
618  if(fDebug>3) Printf("%s:%d hijingGenHeader found", (char*)__FILE__,__LINE__);
619  }
620  }
621 
622  //fh1Trials->Fill("#sum{ntrials}",fAvgTrials);
623  }
624  }
625 
626 
627  // Cut on pThard if fMCEvent and pThard >= 0 and fill histo with #evt before and after the cut
628  if (fMCEvent) {
629  // Before cut
630  fh1EvtsPtHardCut->Fill(0.);
631 
632  if (fUseInclusivePIDtask) {
633  for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
635  }
636  }
637 
638  if (fUseJetPIDtask) {
639  for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
641  }
642  }
643 
644  if (fUseJetUEPIDtask) {
645  for (Int_t i = 0; i < fNumJetUEPIDtasks; i++) {
647  }
648  }
649 
650  // Cut
651  if (fMCPtHardCut >= 0. && ptHard >= fMCPtHardCut) {
652  if (fDebug>3) Printf("%s:%d skipping event with pThard %f (>= %f)", (char*)__FILE__,__LINE__, ptHard, fMCPtHardCut);
653  PostData(1, fCommonHistList);
654  return kFALSE;
655  }
656 
657  // After cut
658  fh1EvtsPtHardCut->Fill(1.);
659 
660  if (fUseInclusivePIDtask) {
661  for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
663  }
664  }
665 
666  if (fUseJetPIDtask) {
667  for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
669  }
670  }
671 
672  if (fUseJetUEPIDtask) {
673  for (Int_t i = 0; i < fNumJetUEPIDtasks; i++) {
675  }
676  }
677 
678  }
679 
680  // Trigger selection
681  AliInputEventHandler* inputHandler = (AliInputEventHandler*)
682  ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
683 
684  if(!(inputHandler->IsEventSelected() & fEvtSelectionMask)){
685  fh1EvtSelection->Fill(1.);
686  if (fDebug > 1 ) Printf(" Trigger Selection: event REJECTED ... ");
687  PostData(1, fCommonHistList);
688  return kFALSE;
689  }
690 
691  fESD = dynamic_cast<AliESDEvent*>(InputEvent());
692  if(!fESD){
693  if(fDebug>3) Printf("%s:%d ESDEvent not found in the input", (char*)__FILE__,__LINE__);
694  }
695 
696  // get AOD event from input/ouput
697  TObject* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
698  if( handler && handler->InheritsFrom("AliAODInputHandler") ) {
699  fAOD = ((AliAODInputHandler*)handler)->GetEvent();
701  if (fDebug > 1) Printf("%s:%d AOD event from input", (char*)__FILE__,__LINE__);
702  }
703  else {
704  handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
705  if( handler && handler->InheritsFrom("AliAODHandler") ) {
706  fAOD = ((AliAODHandler*)handler)->GetAOD();
707  fAODJets = fAOD;
708  if (fDebug > 1) Printf("%s:%d AOD event from output", (char*)__FILE__,__LINE__);
709  }
710  }
711 
712  if(!fAODJets && !fUseAODInputJets){ // case we have AOD in input & output and want jets from output
713  TObject* outHandler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
714  if( outHandler && outHandler->InheritsFrom("AliAODHandler") ) {
715  fAODJets = ((AliAODHandler*)outHandler)->GetAOD();
716  if (fDebug > 1) Printf("%s:%d jets from output AOD", (char*)__FILE__,__LINE__);
717  }
718  }
719 
720  if(fNonStdFile.Length()!=0){
721  // case we have an AOD extension - fetch the jets from the extended output
722 
723  AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
724  fAODExtension = (aodH?aodH->GetExtension(fNonStdFile.Data()):0);
725  if(!fAODExtension){
726  if(fDebug>1)Printf("AODExtension not found for %s",fNonStdFile.Data());
727  }
728  }
729 
730  if(!fAOD){
731  Printf("%s:%d AODEvent not found", (char*)__FILE__,__LINE__);
732  return kFALSE;
733  }
734  if(!fAODJets){
735  Printf("%s:%d AODEvent with jet branch not found", (char*)__FILE__,__LINE__);
736  return kFALSE;
737  }
738 
739 
740  // event selection **************************************************
741  // *** event class ***
742  AliVEvent* evtForCentDetermination = handler->InheritsFrom("AliAODInputHandler") ? fAOD : InputEvent();
743 
744  Double_t centPercent = -1;
745 
746  if(fEventClass>0){
747  Int_t cl = 0;
748  if(handler->InheritsFrom("AliAODInputHandler")){
749  // since it is not supported by the helper task define own classes
750  centPercent = ((AliAODHeader*)fAOD->GetHeader())->GetCentrality();
751  cl = 1;
752  if(centPercent>10) cl = 2;
753  if(centPercent>30) cl = 3;
754  if(centPercent>50) cl = 4;
755  }
756  else {
758  if(fESD) centPercent = fESD->GetCentrality()->GetCentralityPercentile(fCentralityEstimator.Data());
759  }
760 
761  if(cl!=fEventClass){
762  // event not in selected event class, reject event
763  if (fDebug > 1) Printf("%s:%d event not in selected event class: event REJECTED ...",(char*)__FILE__,__LINE__);
764  fh1EvtSelection->Fill(2.);
765  PostData(1, fCommonHistList);
766  return kFALSE;
767  }
768  }
769 
770  if (fCentralityEstimator.Contains("NoCentrality",TString::kIgnoreCase)) {
771  centPercent = -1;
772  }
773  else {
774  if (!fIsPP) centPercent = evtForCentDetermination->GetCentrality()->GetCentralityPercentile(fCentralityEstimator.Data());
775  }
776 
777  AliPIDResponse *fPIDResponse = inputHandler->GetPIDResponse();
778  if (!fPIDResponse) {
779  AliError("PIDResponse object was not created");
780  }
781 
782  fPIDResponse->SetCurrentCentrality(centPercent);
783 
784  // Retrieve reference multiplicities in |eta|<0.8 and <0.5
785  const Int_t refMult5 = ((AliAODHeader*)fAOD->GetHeader())->GetRefMultiplicityComb05();
786  const Int_t refMult8 = ((AliAODHeader*)fAOD->GetHeader())->GetRefMultiplicityComb08();
787  const Double_t centPercentPP = fAnaUtils->GetMultiplicityPercentile(fAOD, "V0M");
788 
789 
790  // Count events with trigger selection, note: Set centrality percentile fix to -1 for pp for PID framework
791  if (fUseInclusivePIDtask) {
792  for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
794  }
795  }
796 
797  if (fUseJetPIDtask) {
798  for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
800  }
801  }
802 
803  if (fUseJetUEPIDtask) {
804  for (Int_t i = 0; i < fNumJetUEPIDtasks; i++) {
806  }
807  }
808 
809  // *** vertex cut ***
810  AliAODVertex* primVtx = fAOD->GetPrimaryVertex();
811  if (!primVtx) {
812  Printf("%s:%d Primary vertex not found", (char*)__FILE__,__LINE__);
813  return kFALSE;
814  }
815 
816  Int_t nTracksPrim = primVtx->GetNContributors();
817  fh1VertexNContributors->Fill(nTracksPrim);
818 
819 
820  if (fDebug > 1) Printf("%s:%d primary vertex selection: %d", (char*)__FILE__,__LINE__,nTracksPrim);
821  if(nTracksPrim <= 0) {
822  if (fDebug > 1) Printf("%s:%d primary vertex selection: event REJECTED...",(char*)__FILE__,__LINE__);
823  fh1EvtSelection->Fill(3.);
824  PostData(1, fCommonHistList);
825  return kFALSE;
826  }
827 
828  TString primVtxName(primVtx->GetName());
829 
830  if(primVtxName.CompareTo("TPCVertex",TString::kIgnoreCase) == 1){
831  if (fDebug > 1) Printf("%s:%d primary vertex selection: TPC vertex, event REJECTED...",(char*)__FILE__,__LINE__);
832  fh1EvtSelection->Fill(5.);
833  PostData(1, fCommonHistList);
834  return kFALSE;
835  }
836 
837  // Count events with trigger selection and vtx cut, note: Set centrality percentile fix to -1 for pp for PID framework
838  if (fUseInclusivePIDtask) {
839  for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
841  }
842  }
843 
844  if (fUseJetPIDtask) {
845  for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
847  }
848  }
849 
850  if (fUseJetUEPIDtask) {
851  for (Int_t i = 0; i < fNumJetUEPIDtasks; i++) {
853  }
854  }
855 
856 
857  fh1VertexZ->Fill(primVtx->GetZ());
858 
859  if(TMath::Abs(primVtx->GetZ())>fMaxVertexZ){
860  if (fDebug > 1) Printf("%s:%d primary vertex z = %f: event REJECTED...",(char*)__FILE__,__LINE__,primVtx->GetZ());
861  fh1EvtSelection->Fill(4.);
862  PostData(1, fCommonHistList);
863  return kFALSE;
864  }
865 
866  // Count events with trigger selection, vtx cut and z vtx cut, note: Set centrality percentile fix to -1 for pp for PID framework
867  if (fUseInclusivePIDtask) {
868  for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
870  }
871  }
872 
873  if (fUseJetPIDtask) {
874  for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
876  }
877  }
878 
879  if (fUseJetUEPIDtask) {
880  for (Int_t i = 0; i < fNumJetUEPIDtasks; i++) {
882  }
883  }
884 
885 
886  // Store for each task, whether this task would tag this event as pile-up or not
887  const Int_t arrSizeInclusive = TMath::Max(1, fNumInclusivePIDtasks);
888  const Int_t arrSizeJet = TMath::Max(1, fNumJetPIDtasks);
889  const Int_t arrSizeJetUE = TMath::Max(1, fNumJetUEPIDtasks);
890  Bool_t isPileUpInclusivePIDtask[arrSizeInclusive];
891  Bool_t isPileUpJetPIDtask[arrSizeJet];
892  Bool_t isPileUpJetUEPIDtask[arrSizeJetUE];
893 
894  for (Int_t i = 0; i < arrSizeInclusive; i++)
895  isPileUpInclusivePIDtask[i] = kFALSE;
896 
897  for (Int_t i = 0; i < arrSizeJet; i++)
898  isPileUpJetPIDtask[i] = kFALSE;
899 
900  for (Int_t i = 0; i < arrSizeJetUE; i++)
901  isPileUpJetUEPIDtask[i] = kFALSE;
902 
903  // Check whether there is at least one task that does not reject the event (saves processing time in the following code)
904  Bool_t isPileUpForAllInclusivePIDTasks = kTRUE;
905  Bool_t isPileUpForAllJetPIDTasks = kTRUE;
906  Bool_t isPileUpForAllJetUEPIDTasks = kTRUE;
907 
908  // Count events with trigger selection, vtx cut, z vtx cut and after pile-up rejection (if enabled in that task)
909  // Note: Set centrality percentile fix to -1 for pp for PID framework
910  if (fUseInclusivePIDtask) {
911  for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
912  isPileUpInclusivePIDtask[i] = fInclusivePIDtask[i]->GetIsPileUp(fAOD, fInclusivePIDtask[i]->GetPileUpRejectionType());
913  if (!isPileUpInclusivePIDtask[i])
915 
916  isPileUpForAllInclusivePIDTasks = isPileUpForAllInclusivePIDTasks && isPileUpInclusivePIDtask[i];
917  }
918  }
919 
920  if (fUseJetPIDtask) {
921  for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
922  isPileUpJetPIDtask[i] = fJetPIDtask[i]->GetIsPileUp(fAOD, fJetPIDtask[i]->GetPileUpRejectionType());
923  if (!isPileUpJetPIDtask[i])
925 
926  isPileUpForAllJetPIDTasks = isPileUpForAllJetPIDTasks && isPileUpJetPIDtask[i];
927  }
928  }
929 
930  if (fUseJetUEPIDtask) {
931  for (Int_t i = 0; i < fNumJetUEPIDtasks; i++) {
932  isPileUpJetUEPIDtask[i] = fJetUEPIDtask[i]->GetIsPileUp(fAOD, fJetUEPIDtask[i]->GetPileUpRejectionType());
933  if (!isPileUpJetUEPIDtask[i])
935 
936  isPileUpForAllJetUEPIDTasks = isPileUpForAllJetUEPIDTasks && isPileUpJetUEPIDtask[i];
937  }
938  }
939 
940 
941 
942  if (fDebug > 1) Printf("%s:%d event ACCEPTED ...",(char*)__FILE__,__LINE__);
943  fh1EvtSelection->Fill(0.);
944  fh1VtxSelection->Fill(primVtx->GetType());
945  fh1EvtCent->Fill(centPercent);
946 
947  // Set centrality percentile fix to -1 for pp to be used for the PID framework
948  if (fIsPP)
949  centPercent = -1;
950 
951  // Call ConfigureTaskForCurrentEvent of PID tasks to ensure that everything is set up properly for the current event
952  // (e.g. run/period dependence of max eta variation map)
953  if (fUseInclusivePIDtask) {
954  for (Int_t i = 0; i < fNumInclusivePIDtasks; i++)
955  if (!isPileUpInclusivePIDtask[i])
957  }
958 
959  if (fUseJetPIDtask) {
960  for (Int_t i = 0; i < fNumJetPIDtasks; i++)
961  if (!isPileUpJetPIDtask[i])
963  }
964 
965  if (fUseJetUEPIDtask) {
966  for (Int_t i = 0; i < fNumJetUEPIDtasks; i++)
967  if (!isPileUpJetUEPIDtask[i])
969  }
970 
971 
972  //___ fill MC information __________________________________________________________________
973 
974  fh1Trials->Fill("#sum{ntrials}",fAvgTrials);
975 
976  if (fUseInclusivePIDtask) {
977  for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
978  if (!isPileUpInclusivePIDtask[i])
980  }
981  }
982 
983  if (fUseJetPIDtask) {
984  for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
985  if (!isPileUpJetPIDtask[i])
987  }
988  }
989 
990  if (fUseJetUEPIDtask) {
991  for (Int_t i = 0; i < fNumJetUEPIDtasks; i++) {
992  if (!isPileUpJetUEPIDtask[i])
994  }
995  }
996 
997  if (fMCEvent) {
998  fh1PtHard->Fill(ptHard);
999  fh1PtHardTrials->Fill(ptHard,nTrials);
1000  }
1001 
1002  AliPIDResponse* pidResponse = 0x0;
1003  Bool_t tuneOnDataTPC = kFALSE;
1005  if (!inputHandler) {
1006  AliFatal("Input handler needed");
1007  return kFALSE;
1008  }
1009  else {
1010  // PID response object
1011  pidResponse = inputHandler->GetPIDResponse();
1012  if (!pidResponse) {
1013  AliFatal("PIDResponse object was not created");
1014  return kFALSE;
1015  }
1016  else {
1017  tuneOnDataTPC = pidResponse->IsTunedOnData() &&
1018  ((pidResponse->GetTunedOnDataMask() & AliPIDResponse::kDetTPC) == AliPIDResponse::kDetTPC);
1019  }
1020  }
1021  }
1022 
1023  if(fDebug>2)Printf("%s:%d Starting processing...",(char*)__FILE__,__LINE__);
1024 
1025  //____ analysis, fill histos ___________________________________________________
1026 
1029  Double_t trackEtaMin = -0.9;
1030  Double_t trackEtaMax = 0.9;
1031  if (trackContainer) {
1032  trackEtaMin = trackContainer->GetParticleEtaMin();
1033  trackEtaMax = trackContainer->GetParticleEtaMax();
1034  }
1035 
1036  // Fill efficiency for generated primaries and also fill histos for generated yields (primaries + all)
1037  // Efficiency, inclusive - particle level
1038  if(fDebug>2)Printf("%s:%d Starting Inclusive Efficiency...",(char*)__FILE__,__LINE__);
1039  if (fUseInclusivePIDtask && mcParticleContainer && !isPileUpForAllInclusivePIDTasks) {
1040  for (auto part : mcParticleContainer->accepted()) {
1041 
1042  if (!part)
1043  continue;
1044 
1045  // Define clean MC sample with corresponding particle level track cuts:
1046  // - MC-track must be in desired eta range
1047  // - MC-track must be physical primary
1048  // - Species must be one of those in question
1049 
1050  if (part->Eta() > trackEtaMax || part->Eta() < trackEtaMin)
1051  continue;
1052 
1053  Int_t mcID = AliAnalysisTaskPID::PDGtoMCID(part->GetPdgCode());
1054 
1055  // Following lines are not needed - just keep other species (like casecades) - will end up in overflow bin
1056  // and only affect the efficiencies for all (i.e. not identified) what is desired!
1057  //if (mcID == AliPID::kUnknown)
1058  // continue;
1059 
1060  if (!part->IsPhysicalPrimary())
1061  continue;
1062  /*
1063  Int_t iMother = part->GetMother();
1064  if (iMother >= 0)
1065  continue; // Not a physical primary
1066  */
1067 
1068  Double_t pT = part->Pt();
1069 
1070  // AliAODMCParticle->Charge() calls TParticlePDG->Charge(), which returns the charge in units of e0 / 3
1071  Double_t chargeMC = part->Charge() / 3.;
1072 
1073  if (TMath::Abs(chargeMC) < 0.01)
1074  continue; // Reject neutral particles (only relevant, if mcID is not used)
1075 
1076  Double_t valuesGenYield[AliAnalysisTaskPID::kGenYieldNumAxes] = { static_cast<Double_t>(mcID), pT, centPercent,
1077  -1, -1, -1, -1, -1, -1 };
1078 
1079  for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
1080  if (!isPileUpInclusivePIDtask[i] && fInclusivePIDtask[i]->IsInAcceptedEtaRange(TMath::Abs(part->Eta()))) {
1081  valuesGenYield[fInclusivePIDtask[i]->GetIndexOfChargeAxisGenYield()] = fInclusivePIDtask[i]->GetStoreCharge() ? chargeMC : -2;
1082  fInclusivePIDtask[i]->FillGeneratedYield(valuesGenYield);
1083  }
1084  }
1085 
1086  // Always store the charge for the efficiency (needed for geant-fluka correction)
1087  Double_t valuesEff[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), pT, part->Eta(), chargeMC,
1088  centPercent, -1, -1, -1, -1, -1 };// no jet pT etc since inclusive spectrum
1089  for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
1090  if (!isPileUpInclusivePIDtask[i])
1092  }
1093  }
1094  }
1095 
1097  if (fUseInclusivePIDtask && mcParticleContainer && !isPileUpForAllInclusivePIDTasks) {
1098  //Efficiency, inclusive - detector level
1099  for(auto inclusiveaod : trackContainerEfficiency->accepted()) {
1100  // fill inclusive tracks XXX, they have the same track cuts!
1101 
1102  if(inclusiveaod){
1103  Double_t dEdxTPC = tuneOnDataTPC ? pidResponse->GetTPCsignalTunedOnData(inclusiveaod)
1104  : inclusiveaod->GetTPCsignal();
1105 
1106  if (dEdxTPC <= 0)
1107  continue;
1108 
1109  Bool_t survivedTPCCutMIGeo = AliAnalysisTaskPID::TPCCutMIGeo(inclusiveaod, InputEvent());
1110  Bool_t survivedTPCnclCut = AliAnalysisTaskPID::TPCnclCut(inclusiveaod);
1111 
1112  Int_t label = TMath::Abs(inclusiveaod->GetLabel());
1113 
1114  // find MC track in our list, if available
1115  AliAODMCParticle* gentrack = mcParticleContainer->GetMCParticleWithLabel(label);
1116  Int_t pdg = 0;
1117 
1118  if (gentrack)
1119  pdg = gentrack->GetPdgCode();
1120 
1121  // For efficiency: Reconstructed track has survived all cuts on the detector level (no cut on eta acceptance)
1122  // and has an associated MC track
1123  // -> Check whether associated MC track belongs to the clean MC sample defined above,
1124  // i.e. survives the particle level track cuts
1125  if (gentrack) {
1127 
1128  // Following lines are not needed - just keep other species (like casecades) - will end up in overflow bin
1129  // and only affect the efficiencies for all (i.e. not identified) what is desired!
1130  //if (mcID == AliPID::kUnknown)
1131  // continue;
1132 
1133  // Fill efficiency for reconstructed primaries
1134  if (!gentrack->IsPhysicalPrimary())
1135  continue;
1136  /*
1137  * Int_t iMother = gentrack->GetMother();
1138  * if (iMother >= 0)
1139  * continue; // Not a physical primary
1140  */
1141 
1142  if (gentrack->Eta() > trackEtaMax || gentrack->Eta() < trackEtaMin)
1143  continue;
1144 
1145  // AliAODMCParticle->Charge() calls TParticlePDG->Charge(), which returns the charge in units of e0 / 3
1146  // Always store the charge for the efficiency (needed for geant-fluka correction)
1147  Double_t value[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), gentrack->Pt(), gentrack->Eta(),
1148  gentrack->Charge() / 3., centPercent, -1, -1, -1, -1, -1 };// no jet pT etc since inclusive spectrum
1149  for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
1150  if (!isPileUpInclusivePIDtask[i] &&
1151  ((!fInclusivePIDtask[i]->GetUseTPCCutMIGeo() && !fInclusivePIDtask[i]->GetUseTPCnclCut()) ||
1152  (survivedTPCCutMIGeo && fInclusivePIDtask[i]->GetUseTPCCutMIGeo()) ||
1153  (survivedTPCnclCut && fInclusivePIDtask[i]->GetUseTPCnclCut())))
1155  }
1156 
1157  Double_t valueMeas[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), inclusiveaod->Pt(), inclusiveaod->Eta(),
1158  static_cast<Double_t>(inclusiveaod->Charge()), centPercent,
1159  -1, -1, -1, -1, -1 };// no jet pT etc since inclusive spectrum
1160  for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
1161  if (!isPileUpInclusivePIDtask[i] &&
1162  ((!fInclusivePIDtask[i]->GetUseTPCCutMIGeo() && !fInclusivePIDtask[i]->GetUseTPCnclCut()) ||
1163  (survivedTPCCutMIGeo && fInclusivePIDtask[i]->GetUseTPCCutMIGeo()) ||
1164  (survivedTPCnclCut && fInclusivePIDtask[i]->GetUseTPCnclCut())))
1165  fInclusivePIDtask[i]->FillEfficiencyContainer(valueMeas, AliAnalysisTaskPID::kStepRecWithGenCutsMeasuredObs);
1166  }
1167  }
1168  }
1169  }
1170  }
1171 
1172  if(fDebug>2)Printf("%s:%d Process inclusive tracks...",(char*)__FILE__,__LINE__);
1173  if (fUseInclusivePIDtask && trackContainer && !isPileUpForAllInclusivePIDTasks) {
1174  for(auto part : trackContainer->accepted()) {
1175  if (!part)
1176  continue;
1177 
1178  Double_t dEdxTPC = tuneOnDataTPC ? pidResponse->GetTPCsignalTunedOnData(part)
1179  : part->GetTPCsignal();
1180 
1181  if (dEdxTPC <= 0)
1182  continue;
1183 
1184  Bool_t survivedTPCCutMIGeo = AliAnalysisTaskPID::TPCCutMIGeo(part, InputEvent());
1185  Bool_t survivedTPCnclCut = AliAnalysisTaskPID::TPCnclCut(part);
1186 
1187  Int_t label = TMath::Abs(part->GetLabel());
1188 
1189  // find MC track in our list, if available
1190  AliAODMCParticle* gentrack = mcParticleContainer ? mcParticleContainer->GetMCParticleWithLabel(label) : 0x0;
1191  Int_t pdg = 0;
1192 
1193  if (gentrack)
1194  pdg = gentrack->GetPdgCode();
1195 
1196  for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
1197  if (!isPileUpInclusivePIDtask[i] &&
1198  ((!fInclusivePIDtask[i]->GetUseTPCCutMIGeo() && !fInclusivePIDtask[i]->GetUseTPCnclCut()) ||
1199  (survivedTPCCutMIGeo && fInclusivePIDtask[i]->GetUseTPCCutMIGeo()) ||
1200  (survivedTPCnclCut && fInclusivePIDtask[i]->GetUseTPCnclCut()))) {
1201  if (fInclusivePIDtask[i]->IsInAcceptedEtaRange(TMath::Abs(part->Eta())))
1202  fInclusivePIDtask[i]->ProcessTrack(part, pdg, centPercent, -1, kFALSE, kTRUE, fStoreXi, -1, -1); // no jet pT etc since inclusive spectrum
1203  }
1204  }
1205 
1206  if (gentrack) {
1208  // Always store the charge for the efficiency (needed for geant-fluka correction)
1209  Double_t valueRecAllCuts[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), part->Pt(), part->Eta(),
1210  static_cast<Double_t>(part->Charge()), centPercent,
1211  -1, -1, -1, -1, -1 };
1212  for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
1213  if (!isPileUpInclusivePIDtask[i] &&
1214  ((!fInclusivePIDtask[i]->GetUseTPCCutMIGeo() && !fInclusivePIDtask[i]->GetUseTPCnclCut()) ||
1215  (survivedTPCCutMIGeo && fInclusivePIDtask[i]->GetUseTPCCutMIGeo()) ||
1216  (survivedTPCnclCut && fInclusivePIDtask[i]->GetUseTPCnclCut())))
1217  fInclusivePIDtask[i]->FillEfficiencyContainer(valueRecAllCuts, AliAnalysisTaskPID::kStepRecWithRecCutsMeasuredObs);
1218  }
1219 
1220  Double_t weight = IsSecondaryWithStrangeMotherMC(gentrack, mcParticleContainer)
1221  ? GetMCStrangenessFactorCMS(gentrack, mcParticleContainer) : 1.0;
1222  for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
1223  if (!isPileUpInclusivePIDtask[i] &&
1224  ((!fInclusivePIDtask[i]->GetUseTPCCutMIGeo() && !fInclusivePIDtask[i]->GetUseTPCnclCut()) ||
1225  (survivedTPCCutMIGeo && fInclusivePIDtask[i]->GetUseTPCCutMIGeo()) ||
1226  (survivedTPCnclCut && fInclusivePIDtask[i]->GetUseTPCnclCut())))
1227  fInclusivePIDtask[i]->FillEfficiencyContainer(valueRecAllCuts,
1229  weight);
1230  }
1231 
1232  if (gentrack->IsPhysicalPrimary()) {
1233  // AliAODMCParticle->Charge() calls TParticlePDG->Charge(), which returns the charge in units of e0 / 3
1234  // Always store the charge for the efficiency (needed for geant-fluka correction)
1235  Double_t valueGenAllCuts[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), gentrack->Pt(), gentrack->Eta(),
1236  gentrack->Charge() / 3., centPercent, -1, -1, -1, -1, -1 };
1237 
1238  Double_t valuePtResolution[AliAnalysisTaskPID::kPtResNumAxes] = { -1, gentrack->Pt(), part->Pt(),
1239  gentrack->Charge() / 3., centPercent };
1240 
1241  for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
1242  if (!isPileUpInclusivePIDtask[i] &&
1243  ((!fInclusivePIDtask[i]->GetUseTPCCutMIGeo() && !fInclusivePIDtask[i]->GetUseTPCnclCut()) ||
1244  (survivedTPCCutMIGeo && fInclusivePIDtask[i]->GetUseTPCCutMIGeo()) ||
1245  (survivedTPCnclCut && fInclusivePIDtask[i]->GetUseTPCnclCut()))) {
1246  fInclusivePIDtask[i]->FillEfficiencyContainer(valueRecAllCuts,
1248  fInclusivePIDtask[i]->FillEfficiencyContainer(valueGenAllCuts,
1250 
1251  fInclusivePIDtask[i]->FillPtResolution(mcID, valuePtResolution);
1252  }
1253  }
1254  }
1255  }
1256  }
1257  }
1258 
1259  if(fDebug>2)Printf("%s:%d Process Jets - efficiency, particle level..",(char*)__FILE__,__LINE__);
1261  if (fOnlyLeadingJets)
1262  mcJetContainer->SortArray();
1263 
1264  if (fUseJetPIDtask && mcJetContainer && !isPileUpForAllJetPIDTasks) {
1265 
1266  for (Int_t i=0;i<fNumJetPIDtasks;i++) {
1267  if (!isPileUpJetPIDtask[i]) {
1268 
1269  for(auto jet : mcJetContainer->accepted()) {
1270 
1271  if(!jet)
1272  continue;
1273 
1274  Float_t jetPt = jet->Pt();
1275  if (mcJetContainer->GetRhoParameter())
1276  jetPt = jetPt - mcJetContainer->GetRhoVal() * jet->Area();
1277 
1278  fJetPIDtask[i]->FillGenJets(centPercent, jetPt);
1279 
1280  TClonesArray *arrayWithTracks = mcParticleContainer ? mcParticleContainer->GetArray() : 0x0;
1281 
1282  for(Int_t it=0; it<jet->GetNumberOfTracks(); ++it) {
1283 
1284  AliVParticle* trackVP = dynamic_cast<AliVParticle*>(jet->TrackAt(it,arrayWithTracks));
1285  PerformJetMonteCarloAnalysisGeneratedYield(jet, trackVP, fJetPIDtask[i], centPercent);
1286  }
1287 
1288  if (fOnlyLeadingJets)
1289  break;
1290  }
1291  }
1292  }
1293  }
1294 
1295  if(fDebug>2)Printf("%s:%d Process jets...",(char*)__FILE__,__LINE__);
1297  if (fOnlyLeadingJets)
1298  jetContainer->SortArray();
1299 
1300  if (fUseJetPIDtask && jetContainer && !isPileUpForAllJetPIDTasks) {
1301 
1302  for (auto jet : jetContainer->accepted()) {
1303  if (!jet)
1304  continue;
1305 
1306  Float_t jetPt = jet->Pt();
1307  if (jetContainer->GetRhoParameter())
1308  jetPt = jetPt - jetContainer->GetRhoVal() * jet->Area();
1309 
1310  for (Int_t i=0;i<fNumJetPIDtasks;++i) {
1311  if (!isPileUpJetPIDtask[i])
1312  fJetPIDtask[i]->FillRecJets(centPercent, jetPt);
1313  }
1314 
1315  fhJetPtRefMultEta5->Fill(refMult5, jetPt);
1316  fhJetPtRefMultEta8->Fill(refMult8, jetPt);
1317  fhJetPtMultPercent->Fill(centPercentPP, jetPt);
1318 
1319  TClonesArray *arrayWithTracks = trackContainer->GetArray();
1320 
1321  for(Int_t j=0; j<jet->GetNumberOfTracks(); ++j) {
1322 
1323  AliVTrack *track = dynamic_cast<AliVTrack*>(jet->TrackAt(j,arrayWithTracks));
1324  Bool_t trackRejectedByTask[arrSizeJet];
1325  Bool_t trackRejectedByAllTasks = kTRUE;
1326  if (track) {
1327  Bool_t survivedTPCCutMIGeo = AliAnalysisTaskPID::TPCCutMIGeo(track, InputEvent());
1328  Bool_t survivedTPCnclCut = AliAnalysisTaskPID::TPCnclCut(track);
1329  Double_t dEdxTPC = tuneOnDataTPC ? pidResponse->GetTPCsignalTunedOnData(track) : track->GetTPCsignal();
1330  for (Int_t i=0;i<fNumJetPIDtasks;++i) {
1331  trackRejectedByTask[i] = isPileUpJetPIDtask[i] || !((!fJetPIDtask[i]->GetUseTPCCutMIGeo() && !fJetPIDtask[i]->GetUseTPCnclCut()) || (survivedTPCCutMIGeo && fJetPIDtask[i]->GetUseTPCCutMIGeo()) || (survivedTPCnclCut && fJetPIDtask[i]->GetUseTPCnclCut()) && fJetPIDtask[i]->IsInAcceptedEtaRange(TMath::Abs(track->Eta())) && (dEdxTPC > 0.0));
1332  trackRejectedByAllTasks = trackRejectedByAllTasks && trackRejectedByTask[i];
1333  }
1334 
1335  if (!trackRejectedByAllTasks) {
1336  AnalyseJetTrack(track, jet, 0x0, trackRejectedByTask, centPercent, mcParticleContainer);
1337  FillDCA(track, mcParticleContainer);
1338  }
1339  }
1340 
1341  if (fOnlyLeadingJets)
1342  break;
1343  }
1344  }
1345  }
1346 
1347  if(fDebug>2)Printf("%s:%d Process Underlying event...",(char*)__FILE__,__LINE__);
1348  if (fUseJetUEPIDtask && jetContainer && !isPileUpForAllJetUEPIDTasks) {
1349  for (Int_t i=0;i<fNumJetUEPIDtasks;++i) {
1350  if (!isPileUpJetUEPIDtask[i]) {
1351  TString method = fUEMethods[i];
1352  AliAnalysisTaskPID* task = fJetUEPIDtask[i];
1353 
1354  TList *jetUElist = 0x0;
1355  TList* mcJetUElist = 0x0;
1356 
1357  if (mcJetContainer) {
1358  if (fOnlyLeadingJets)
1359  mcJetContainer->SortArray();
1360 
1361  if (method.Contains("RC",TString::kIgnoreCase) || method.Contains("Random",TString::kIgnoreCase)) {
1362  mcJetUElist = GetUEJetsWithRandomConeMethod(mcJetContainer, TMath::Abs(GetFFRadius()), task->GetEtaAbsCutUp());
1363  }
1364  else if (method.Contains("PC",TString::kIgnoreCase) || method.Contains("Perpendicular",TString::kIgnoreCase)) {
1365  mcJetUElist = GetUEJetsWithPerpendicularConeMethod(mcJetContainer);
1366  }
1367  }
1368 
1369  if (jetContainer) {
1370  //Get Underlying event jets
1371  if (fOnlyLeadingJets)
1372  jetContainer->SortArray();
1373  //Random Cones______________________________________________
1374  if (method.Contains("RC",TString::kIgnoreCase) || method.Contains("Random",TString::kIgnoreCase)) {
1375  jetUElist = GetUEJetsWithRandomConeMethod(jetContainer, TMath::Abs(GetFFRadius()), task->GetEtaAbsCutUp());
1376  }
1377  //Perpendicular Cones______________________________________________
1378  else if (method.Contains("PC",TString::kIgnoreCase) || method.Contains("Perpendicular",TString::kIgnoreCase)) {
1379  jetUElist = GetUEJetsWithPerpendicularConeMethod(jetContainer);
1380  }
1381  }
1382  //Process Tracks
1383  if (jetUElist) {
1384  for (Int_t i=0;i<jetUElist->GetEntries();++i) {
1385  AliEmcalJet* jet = (AliEmcalJet*)(jetUElist->At(i));
1386  Double_t UEPtDensity = 0.0;
1387  Float_t jetPt = jet->Pt();
1388  if (jetContainer->GetRhoParameter())
1389  jetPt = jetPt - jetContainer->GetRhoVal() * jet->Area();
1390  task->FillRecJets(centPercent, jetPt);
1391  task->FillJetArea(centPercent, jet->Area());
1392 
1393  TList *jetUEtracklist = GetTracksInCone(jet, trackContainer);
1394  if (!jetUEtracklist) {
1395  cout << "No track list for an underlying event jet. Check Code!" << endl;
1396  continue;
1397  }
1398 
1399  jetUEtracklist->SetOwner(kFALSE);
1400 
1401  for (Int_t j=0;j<jetUEtracklist->GetEntries();++j) {
1402  AliVTrack* UEtrack = (AliVTrack*)(jetUEtracklist->At(j));
1403 
1404  if (UEtrack) {
1405  Bool_t survivedTPCCutMIGeo = AliAnalysisTaskPID::TPCCutMIGeo(UEtrack, InputEvent());
1406  Bool_t survivedTPCnclCut = AliAnalysisTaskPID::TPCnclCut(UEtrack);
1407  Double_t dEdxTPC = tuneOnDataTPC ? pidResponse->GetTPCsignalTunedOnData(UEtrack) : UEtrack->GetTPCsignal();
1408 
1409  if ((!task->GetUseTPCCutMIGeo() && !task->GetUseTPCnclCut()) || (survivedTPCCutMIGeo && task->GetUseTPCCutMIGeo()) || (survivedTPCnclCut && task->GetUseTPCnclCut()) && task->IsInAcceptedEtaRange(TMath::Abs(UEtrack->Eta())) && (dEdxTPC > 0.0))
1410  AnalyseJetTrack(UEtrack, jet, task, 0x0, centPercent, mcParticleContainer);
1411 
1412  UEPtDensity += UEtrack->Pt();
1413  }
1414  }
1415  delete jetUEtracklist;
1416  jetUEtracklist = 0x0;
1417  UEPtDensity = UEPtDensity/(TMath::Abs(GetFFRadius()) * TMath::Abs(GetFFRadius()) * TMath::Pi());
1418  task->FillUEDensity(centPercent, UEPtDensity);
1419  }
1420  delete jetUElist;
1421  jetUElist = 0x0;
1422  }
1423  if (mcJetUElist) {
1424  for (Int_t i=0;i<mcJetUElist->GetEntries();++i) {
1425  AliEmcalJet* jet = (AliEmcalJet*)(mcJetUElist->At(i));
1426 
1427  task->FillGenJets(centPercent, jet->Pt());
1428  TList* mcJetUEtracklist = GetTracksInCone(jet, mcParticleContainer);
1429  if (!mcJetUEtracklist) {
1430  cout << "No track list for MC Tracks in Underlying event. Check Code!" << endl;
1431  continue;
1432  }
1433  mcJetUEtracklist->SetOwner(kFALSE);
1434  for (Int_t j=0;j<mcJetUEtracklist->GetEntries();++j) {
1435  AliVParticle* mcUEParticle = dynamic_cast<AliVParticle*>(mcJetUEtracklist->At(j));
1436  PerformJetMonteCarloAnalysisGeneratedYield(jet, mcUEParticle, task, centPercent);
1437  }
1438  }
1439  }
1440  }
1441  }
1442  }
1443  //End of underlying event
1444 
1445  //___________________
1446 
1447  if(fDebug > 2) Printf("%s:%d Processing done, posting output data now...",(char*)__FILE__,__LINE__);
1448 
1449  //Post output data.
1450  PostData(1, fCommonHistList);
1451 
1452  if (fUseInclusivePIDtask) {
1453  for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
1455  }
1456  }
1457 
1458  if (fUseJetPIDtask) {
1459  for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
1461  }
1462  }
1463 
1464  if (fUseJetUEPIDtask) {
1465  for (Int_t i = 0; i < fNumJetUEPIDtasks; i++) {
1467  }
1468  }
1469 
1470  if(fDebug > 2) Printf("%s:%d Done",(char*)__FILE__,__LINE__);
1471 
1472  return kTRUE;
1473 }
1474 
1475 //______________________________________________________________
1477 {
1478  // terminated
1479  if (fUseJetUEPIDtask) {
1480  for (Int_t i=0;i<fNumJetUEPIDtasks;++i) {
1481  fJetUEPIDtask[i]->NormalizeJetArea(TMath::Abs(GetFFRadius()));
1482  }
1483  }
1484 
1485  if(fDebug > 1) printf("AliAnalysisTaskIDFragmentationFunction::Terminate() \n");
1486 }
1487 
1488 // _________________________________________________________________________________________________________
1489 void AliAnalysisTaskIDFragmentationFunction::SetProperties(THnSparse* h, Int_t dim, const char** labels)
1490 {
1491  // Set properties of THnSparse
1492 
1493  for(Int_t i=0; i<dim; i++){
1494  h->GetAxis(i)->SetTitle(labels[i]);
1495  h->GetAxis(i)->SetTitleColor(1);
1496  }
1497 }
1498 
1499 // __________________________________________________________________________________________
1500 void AliAnalysisTaskIDFragmentationFunction::SetProperties(TH1* h,const char* x, const char* y)
1501 {
1502  //Set properties of histos (x and y title)
1503 
1504  h->SetXTitle(x);
1505  h->SetYTitle(y);
1506  h->GetXaxis()->SetTitleColor(1);
1507  h->GetYaxis()->SetTitleColor(1);
1508 }
1509 
1510 // _________________________________________________________________________________________________________
1511 void AliAnalysisTaskIDFragmentationFunction::SetProperties(TH1* h,const char* x, const char* y, const char* z)
1512 {
1513  //Set properties of histos (x,y and z title)
1514 
1515  h->SetXTitle(x);
1516  h->SetYTitle(y);
1517  h->SetZTitle(z);
1518  h->GetXaxis()->SetTitleColor(1);
1519  h->GetYaxis()->SetTitleColor(1);
1520  h->GetZaxis()->SetTitleColor(1);
1521 }
1522 
1523 
1524 //_______________________________________________________
1526  Bool_t useInternalFlag) const
1527 {
1528  // Calculate the distance between jet and track in the eta-phi-plane. If useInternalFlag is on and fStudyTransversalJetStructure is kFALSE,
1529  // the function just returns -1
1530 
1531  if (useInternalFlag && !fStudyTransversalJetStructure)
1532  return -1;
1533 
1534  if (!jet || !track)
1535  return -1.;
1536 
1537  TVector3 v_track(track->Px(), track->Py(), track->Pz());
1538  TVector3 v_jet(jet->Px(), jet->Py(), jet->Pz());
1539 
1540  return v_track.DeltaR(v_jet);
1541 }
1542 
1543 
1544 //_______________________________________________________
1546 {
1547  // factor strangeness data/MC as function of pt from UE analysis (Sara Vallero)
1548 
1549  Double_t alpha = 1;
1550 
1551  if(0.150<pt && pt<0.200) alpha = 3.639;
1552  if(0.200<pt && pt<0.250) alpha = 2.097;
1553  if(0.250<pt && pt<0.300) alpha = 1.930;
1554  if(0.300<pt && pt<0.350) alpha = 1.932;
1555  if(0.350<pt && pt<0.400) alpha = 1.943;
1556  if(0.400<pt && pt<0.450) alpha = 1.993;
1557  if(0.450<pt && pt<0.500) alpha = 1.989;
1558  if(0.500<pt && pt<0.600) alpha = 1.963;
1559  if(0.600<pt && pt<0.700) alpha = 1.917;
1560  if(0.700<pt && pt<0.800) alpha = 1.861;
1561  if(0.800<pt && pt<0.900) alpha = 1.820;
1562  if(0.900<pt && pt<1.000) alpha = 1.741;
1563  if(1.000<pt && pt<1.500) alpha = 0.878;
1564 
1565  return alpha;
1566 }
1567 
1568 //__________________________________________________________________________________________________
1570  AliMCParticleContainer* mcParticleContainer) const
1571 {
1572  // strangeness ratio MC/data as function of mother pt from CMS data in |eta|<2.0
1573 
1574  if(!mcParticleContainer) return 1;
1575 
1576  AliAODMCParticle* currentMother = daughter;
1577  AliAODMCParticle* currentDaughter = daughter;
1578 
1579 
1580  // find first primary mother K0s, Lambda or Xi
1581  while(1){
1582 
1583  Int_t daughterPDG = currentDaughter->GetPdgCode();
1584 
1585  Int_t motherLabel = currentDaughter->GetMother();
1586 
1587  currentMother = (AliAODMCParticle*) mcParticleContainer->GetMCParticleWithLabel(motherLabel);
1588 
1589  if(!currentMother){
1590  currentMother = currentDaughter;
1591  break;
1592  }
1593 
1594  Int_t motherPDG = currentMother->GetPdgCode();
1595 
1596  // phys. primary found ?
1597  if(currentMother->IsPhysicalPrimary()) break;
1598 
1599  if(TMath::Abs(daughterPDG) == 321){ // K+/K- e.g. from phi (ref data not feeddown corrected)
1600  currentMother = currentDaughter; break;
1601  }
1602  if(TMath::Abs(motherPDG) == 310 ){ // K0s e.g. from phi (ref data not feeddown corrected)
1603  break;
1604  }
1605  if(TMath::Abs(motherPDG) == 3212 && TMath::Abs(daughterPDG) == 3122){ // mother Sigma0, daughter Lambda (this case not included in feeddown corr.)
1606  currentMother = currentDaughter; break;
1607  }
1608 
1609  currentDaughter = currentMother;
1610  }
1611 
1612 
1613  Int_t motherPDG = currentMother->GetPdgCode();
1614  Double_t motherGenPt = currentMother->Pt();
1615 
1616  return AliAnalysisTaskPID::GetMCStrangenessFactorCMS(motherPDG, motherGenPt);
1617 }
1618 
1619 
1620 //_______________________________________________________
1622  Bool_t useInternalFlag) const
1623 {
1624  // Calculate the momentum of the track perpendicular to the jet axis. If useInternalFlag is on and fStudyTransversalJetStructure is kFALSE,
1625  // the function just returns -1
1626 
1627  if (useInternalFlag && !fStudyTransversalJetStructure)
1628  return -1;
1629 
1630  if (!jet || !track)
1631  return -1.;
1632 
1633  TVector3 v_track(track->Px(), track->Py(), track->Pz());
1634  TVector3 v_jet(jet->Px(), jet->Py(), jet->Pz());
1635 
1636  // Important that the jet is the argument, otherwise one would get the transverse jet momentum w.r.t. the track direction
1637  return v_track.Perp(v_jet);
1638 }
1639 
1640 
1641 // _________________________________________________________________________________
1643 {
1644  // Check whether particle is a secondary with strange mother, i.e. returns kTRUE if a strange mother is found
1645  // and the particle is NOT a physical primary. In all other cases kFALSE is returned
1646 
1647  if (!mcParticleContainer || !part)
1648  return kFALSE;
1649 
1650  if (part->IsPhysicalPrimary())
1651  return kFALSE;
1652 
1653  Int_t iMother = part->GetMother();
1654  if (iMother < 0)
1655  return kFALSE;
1656 
1657 
1658  AliAODMCParticle* partM = dynamic_cast<AliAODMCParticle*>(mcParticleContainer->GetMCParticleWithLabel(iMother));
1659  if (!partM)
1660  return kFALSE;
1661 
1662  Int_t codeM = TMath::Abs(partM->GetPdgCode());
1663  Int_t mfl = Int_t(codeM / TMath::Power(10, Int_t(TMath::Log10(codeM))));
1664  if (mfl == 3 && codeM != 3) // codeM = 3 is for s quark
1665  return kTRUE;
1666 
1667  return kFALSE;
1668 }
1669 
1671  TList *jetUElist = new TList();
1672  jetUElist->SetOwner(kTRUE);
1673  if (jetContainer) {
1674  AliEmcalJet* jetRC;
1675  for (auto jet : jetContainer->accepted()) {
1676  jetRC = GetRandomCone(jet, maxEtaTrack - coneRadius, 2.0 * coneRadius);
1677  if (jetRC) {
1678  jetUElist->Add(jetRC);
1679  }
1680  if (fOnlyLeadingJets)
1681  break;
1682  }
1683  }
1684  return jetUElist;
1685 }
1686 
1688  TList *jetUElist = new TList();
1689  jetUElist->SetOwner(kTRUE);
1690  if (jetContainer) {
1691  Double_t perpAngle = TMath::Pi()/2.0;
1692  for (auto jet : jetContainer->accepted()) {
1693  AliEmcalJet* jetPC = 0x0;
1694  jetPC = GetPerpendicularCone(jet,perpAngle);
1695  if (jetPC) {
1696  jetUElist->Add(jetPC);
1697  }
1698  jetPC = GetPerpendicularCone(jet,-perpAngle);
1699  if (jetPC) {
1700  jetUElist->Add(jetPC);
1701  }
1702  if (fOnlyLeadingJets)
1703  break;
1704  }
1705  }
1706  return jetUElist;
1707 }
1708 
1710  if (!processedJet)
1711  return 0x0;
1712 
1713  TLorentzVector vecRdCone;
1714  AliEmcalJet* jetRC = 0x0; //random cone candidate
1715  Double_t dEta, dPhi; //random eta and phi value for RC
1716  UInt_t i = 0;
1717  Double_t jetPt = processedJet->Pt(); //ATTENTION Do NOT subtract rho here - we want to have an exact copy of the processed jet.
1718  Double_t jetM = processedJet->M();
1719 
1720  //Loop is running until the number of trials equals the maximum number of trials or the created random jet is not overlapping with any jet in the jetlist
1721  do {
1722  delete jetRC; //Delete previously created jet
1723  dEta = dEtaConeMax * (2 * fRandom->Rndm() - 1.0); //random eta value in range: [-dEtaConeMax,+dEtaConeMax]
1724  dPhi = TMath::TwoPi() * fRandom->Rndm(); //random phi value in range: [0,2*Pi]
1725  jetRC = new AliEmcalJet(jetPt, dEta, dPhi, jetM); //new RC candidate
1726  ++i;
1727  }
1728  while (i<fRCTrials && IsRCJCOverlap(jetRC,dDistance));
1729 
1730  //If the last jet is also overlapping, delete it and set pointer to zero
1731  if(IsRCJCOverlap(jetRC,dDistance)) {
1732  delete jetRC;
1733  jetRC = 0x0;
1734  }
1735 
1736  if (jetRC)
1737  jetRC->SetArea(processedJet->Area());
1738 
1739  return jetRC;
1740 }
1741 
1743  if (!processedJet)
1744  return 0x0;
1745 
1746  Double_t dDistance = 2.0 * TMath::Abs(GetFFRadius());
1747 
1748  AliEmcalJet* perpJet = new AliEmcalJet(processedJet->Pt(), processedJet->Eta(), processedJet->Phi()+perpAngle, processedJet->M());
1749 
1750  perpJet->SetArea(processedJet->Area());
1751 
1752  if (!IsParticleInCone(processedJet,perpJet,dDistance))
1753  return perpJet;
1754  else
1755  return 0x0;
1756 }
1757 
1759 
1760  if(!part) return kFALSE;
1761  if(!dDistance) return kFALSE;
1762 
1764 
1765  for(auto jet : jetContainer->accepted()){ //loop over all reconstructed jets in events
1766  if(!jet){
1767  if (fDebug>2) std::cout << "AliAnalysisTaskJetChem::IsRCJCOverlap jet pointer invalid!" << std::endl;
1768  continue;
1769  }
1770  if(IsParticleInCone(jet, part, dDistance) == kTRUE) return kTRUE;//RC and JC are overlapping
1771 
1772  }//end loop testing RC-JC overlap
1773  return kFALSE;//RC and JC are not overlapping -> good!
1774 }
1775 
1776 Bool_t AliAnalysisTaskIDFragmentationFunction::IsParticleInCone(const AliVParticle* part1, const AliVParticle* part2, Double_t dRMax) const {
1777 // decides whether a particle is inside a jet cone, or has a minimum distance to a second track/axis
1778  if (!part1 || !part2)
1779  return kFALSE;
1780 
1781  TVector3 vecMom2(part2->Px(),part2->Py(),part2->Pz());
1782  TVector3 vecMom1(part1->Px(),part1->Py(),part1->Pz());
1783  Double_t dR = vecMom2.DeltaR(vecMom1); // = sqrt(dEta*dEta+dPhi*dPhi)
1784  if(dR < dRMax) // momentum vectors of part1 and part2 are closer than dRMax
1785  return kTRUE;
1786  return kFALSE;
1787 }
1788 
1790  TList* jetTrackList = new TList();
1791  jetTrackList->SetOwner(kFALSE);
1792  if (!jet)
1793  return jetTrackList;
1794 
1795  if (!particleContainer)
1796  particleContainer = GetTrackContainer(GetNameTrackContainer());
1797 
1798  for (auto track : particleContainer->accepted()) {
1799 
1800  if(!track)
1801  continue;
1802 
1803  if (IsParticleInCone(track,jet,TMath::Abs(GetFFRadius()))) {
1804  jetTrackList->Add(track);
1805  }
1806  }
1807  return jetTrackList;
1808 }
1809 
1810 //End of underlying event calculations
1811 
1813  if (!jet || !trackVP || !task)
1814  return;
1815 
1816  if (!mcJetContainer)
1817  mcJetContainer = GetJetContainer(GetNameMCParticleJetContainer());
1818 
1819  if (!mcJetContainer) {
1820  std::cout << "Monte-Carlo Jet Container not found." << std::endl;
1821  }
1822 
1823  Float_t jetPt = jet->Pt();
1824  if (mcJetContainer->GetRhoParameter())
1825  jetPt = jetPt - mcJetContainer->GetRhoVal() * jet->Area();
1826  Float_t trackPt = trackVP->Pt();
1827 
1828  // Efficiency, jets - particle level
1829  AliAODMCParticle* part = dynamic_cast<AliAODMCParticle*>(trackVP);
1830  if (!part) {
1831  AliError("expected ref track not found ");
1832  return;
1833  }
1834 
1835  Double_t trackEtaMin = -1.0 * (task->GetEtaAbsCutUp());
1836  Double_t trackEtaMax = task->GetEtaAbsCutUp();
1837 
1839  if (mcParticleContainer) {
1840  trackEtaMin = mcParticleContainer->GetParticleEtaMin();
1841  trackEtaMax = mcParticleContainer->GetParticleEtaMax();
1842  }
1843 
1844  // Fill efficiency for generated primaries and also fill histos for generated yields (primaries + all)
1845  if (part->Eta() > trackEtaMax || part->Eta() < trackEtaMin)
1846  return;
1847 
1848  Int_t mcID = AliAnalysisTaskPID::PDGtoMCID(part->GetPdgCode());
1849 
1850  // Following lines are not needed - just keep other species (like casecades) - will end up in overflow bin
1851  // and only affect the efficiencies for all (i.e. not identified) what is desired!
1852  //if (mcID == AliPID::kUnknown)
1853  // continue;
1854 
1855  if (!part->IsPhysicalPrimary())
1856  return;
1857  //
1858  // Int_t iMother = part->GetMother();
1859  // if (iMother >= 0)
1860  // continue; // Not a physical primary
1861  //
1862 
1863  Double_t z = -1., xi = -1.;
1864  AliAnalysisTaskPID::GetJetTrackObservables(trackPt, jetPt, z, xi, fStoreXi);
1865  Double_t distance = GetDistanceJetTrack(jet, trackVP);
1866  Double_t jT = GetPerpendicularMomentumTrackJet(jet, trackVP);
1867 
1868  // AliAODMCParticle->Charge() calls TParticlePDG->Charge(), which returns the charge in units of e0 / 3
1869  Double_t chargeMC = part->Charge() / 3.;
1870 
1871  if (TMath::Abs(chargeMC) < 0.01)
1872  return; // Reject neutral particles (only relevant, if mcID is not used)
1873 
1874  Double_t valuesGenYield[AliAnalysisTaskPID::kGenYieldNumAxes] = { static_cast<Double_t>(mcID), trackPt, centPercent, jetPt, z, xi,
1875  chargeMC, distance, jT };
1876 
1877  if (task->IsInAcceptedEtaRange(TMath::Abs(part->Eta()))) {
1878  valuesGenYield[task->GetIndexOfChargeAxisGenYield()] = task->GetStoreCharge() ? chargeMC : -2;
1879  task->FillGeneratedYield(valuesGenYield);
1880  }
1881 
1882  // Always store the charge for the efficiency (needed for geant-fluka correction)
1883  Double_t valuesEff[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), trackPt, part->Eta(), chargeMC,
1884  centPercent, jetPt, z, xi, distance, jT };
1886 }
1887 
1888 void AliAnalysisTaskIDFragmentationFunction::AnalyseJetTrack(AliVTrack* track, AliEmcalJet* jet, AliAnalysisTaskPID* task, const Bool_t* trackRejectedByTask, Double_t centPercent, AliMCParticleContainer* mcParticleContainer) {
1889 
1890  if (!track || (!task && !trackRejectedByTask)) {
1891  std::cout << "Cannot analyse track" << std::endl;
1892  std::cout << "Track: " << track << std::endl;
1893  std::cout << "Task: " << task << std::endl;
1894  std::cout << "Rejection Array: " << trackRejectedByTask << std::endl;
1895  return;
1896  }
1897 
1899 
1900  Float_t jetPt = jet->Pt();
1901  if (jetContainer->GetRhoParameter())
1902  jetPt = jetPt - jetContainer->GetRhoVal() * jet->Area();
1903 
1904  Double_t trackPt = track->Pt();
1905 
1906  Int_t label = TMath::Abs(track->GetLabel());
1907 
1908  Int_t pdg = 0;
1909  Int_t mcID = AliPID::kUnknown;
1910 
1911  Double_t z = -1., xi = -1.;
1912  AliAnalysisTaskPID::GetJetTrackObservables(trackPt, jetPt, z, xi, fStoreXi);
1913  Double_t distance = GetDistanceJetTrack(jet, track);
1914  Double_t jT = GetPerpendicularMomentumTrackJet(jet, track);
1915 
1916  if (!mcParticleContainer)
1917  mcParticleContainer = GetMCParticleContainer(GetNameMCParticleContainer());
1918 
1919  AliAODMCParticle* gentrack = mcParticleContainer ? mcParticleContainer->GetMCParticleWithLabel(label) : 0x0;
1920 
1921  //Get PDG Code
1922  if (gentrack) {
1923  pdg = gentrack->GetPdgCode();
1924  }
1925 
1926  if (task)
1927  task->ProcessTrack(track, pdg, centPercent, jetPt, kFALSE, kTRUE, fStoreXi, distance, jT);
1928  else {
1929  for (Int_t i=0;i<fNumJetPIDtasks;++i) {
1930  if (!trackRejectedByTask[i])
1931  fJetPIDtask[i]->ProcessTrack(track, pdg, centPercent, jetPt, kFALSE, kTRUE, fStoreXi, distance, jT);
1932  }
1933  }
1934 
1935  //Process
1936  if (gentrack) {
1937 
1938  // Secondaries, jets
1939  mcID = AliAnalysisTaskPID::PDGtoMCID(pdg);
1940 
1941  Double_t valueRecAllCuts[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), trackPt, track->Eta(),
1942  static_cast<Double_t>(track->Charge()),
1943  centPercent, jetPt, z, xi, distance, jT };
1944 
1945  Double_t weight = IsSecondaryWithStrangeMotherMC(gentrack, mcParticleContainer)
1946  ? GetMCStrangenessFactorCMS(gentrack, mcParticleContainer) : 1.0;
1947 
1948  if (task) {
1951  weight);
1952  }
1953  else {
1954  for (Int_t i=0;i<fNumJetPIDtasks;++i) {
1955  if (!trackRejectedByTask[i])
1958  weight);
1959  }
1960  }
1961 
1962  if (gentrack->IsPhysicalPrimary()) {
1963  Double_t genPt = gentrack->Pt();
1964  Double_t genZ = -1., genXi = -1.;
1965  AliAnalysisTaskPID::GetJetTrackObservables(genPt, jetPt, genZ, genXi, fStoreXi);
1966  Double_t genDistance = GetDistanceJetTrack(jet, gentrack);
1967  Double_t genJt = GetPerpendicularMomentumTrackJet(jet, gentrack);
1968 
1969  // Always store the charge for the efficiency (needed for geant-fluka correction)
1970  Double_t valueGenAllCuts[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), genPt, gentrack->Eta(),
1971  gentrack->Charge() / 3., centPercent, jetPt, genZ,
1972  genXi, genDistance, genJt };
1973 
1974  Double_t valuePtResolution[AliAnalysisTaskPID::kPtResNumAxes] = { jetPt, genPt, trackPt, gentrack->Charge() / 3., centPercent };
1975 
1976  if (task) {
1979  task->FillPtResolution(mcID, valuePtResolution);
1980  }
1981  else {
1982  for (Int_t i=0;i<fNumJetPIDtasks;++i) {
1983  if (!trackRejectedByTask[i])
1986  fJetPIDtask[i]->FillPtResolution(mcID, valuePtResolution);
1987  }
1988  }
1989 
1990  //Efficiency, jets - detector level
1991  Double_t trackEtaMin, trackEtaMax;
1992  if (task) {
1993  trackEtaMin = -1.0 * (task->GetEtaAbsCutUp());
1994  trackEtaMax = task->GetEtaAbsCutUp();
1995  }
1996  else {
1997  trackEtaMin = -1.0 * (fJetPIDtask[0]->GetEtaAbsCutUp());
1998  trackEtaMax = fJetPIDtask[0]->GetEtaAbsCutUp();
1999  }
2001  if (mcParticleContainer) {
2002  trackEtaMin = mcParticleContainer->GetParticleEtaMin();
2003  trackEtaMax = mcParticleContainer->GetParticleEtaMax();
2004  }
2005 
2006  if (gentrack->Eta() <= trackEtaMax || gentrack->Eta() >= trackEtaMin) {
2007 
2008  Double_t genZ = -1., genXi = -1.;
2009  Double_t genPt = gentrack->Pt();
2010  AliAnalysisTaskPID::GetJetTrackObservables(genPt, jetPt, genZ, genXi, fStoreXi);
2011  Double_t genDistance = GetDistanceJetTrack(jet, gentrack);
2012  Double_t genJt = GetPerpendicularMomentumTrackJet(jet, gentrack);
2013 
2014  Double_t measZ = -1., measXi = -1.;
2015  AliAnalysisTaskPID::GetJetTrackObservables(trackPt, jetPt, measZ, measXi, fStoreXi);
2016  Double_t measDistance = GetDistanceJetTrack(jet, track);
2017  Double_t measJt = GetPerpendicularMomentumTrackJet(jet, track);
2018 
2019  // AliAODMCParticle->Charge() calls TParticlePDG->Charge(), which returns the charge in units of e0 / 3
2020  // Always store the charge for the efficiency (needed for geant-fluka correction)
2021  Double_t value[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), genPt, gentrack->Eta(), gentrack->Charge() / 3., centPercent, jetPt, genZ, genXi, genDistance, genJt };
2022 
2023  Double_t valueMeas[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), trackPt, track->Eta(),
2024  static_cast<Double_t>(track->Charge()),
2025  centPercent, jetPt, measZ, measXi, measDistance, measJt };
2026 
2027  if (task) {
2030  }
2031  else {
2032  for (Int_t i=0;i<fNumJetPIDtasks;++i) {
2033  if (!trackRejectedByTask[i])
2036  }
2037  }
2038  }
2039  }
2040  }
2041  return;
2042 }
2043 
2044 void AliAnalysisTaskIDFragmentationFunction::FillDCA(AliVTrack* track, AliMCParticleContainer* mcParticleContainer) {
2045 
2046  if (!fFillDCA)
2047  return;
2048 
2049  Double_t trackPt = track->Pt();
2050 
2051  Int_t label = TMath::Abs(track->GetLabel());
2052 
2053  if (!mcParticleContainer)
2054  mcParticleContainer = GetMCParticleContainer(GetNameMCParticleContainer());
2055 
2056  Int_t mcID = AliPID::kUnknown;
2057 
2058  AliAODMCParticle* gentrack = mcParticleContainer ? mcParticleContainer->GetMCParticleWithLabel(label) : 0x0;
2059 
2060  if (gentrack) {
2061  Int_t pdg = gentrack->GetPdgCode();
2062  mcID = AliAnalysisTaskPID::PDGtoMCID(pdg);
2063  }
2064 
2065  Double_t dca[2] = {0., 0.}; // 0: xy; 1: z
2066 
2067  Double_t v[3] = {0, };
2068  Double_t pos[3] = {0, };
2069  AliAODVertex* primVtx = fAOD->GetPrimaryVertex();
2070  if (!primVtx) {
2071  std::cout << "Primary Vertex not found, do not fill DCA" << std::endl;
2072  return;
2073  }
2074  primVtx->GetXYZ(v);
2075  track->GetXYZ(pos);
2076  dca[0] = pos[0] - v[0];
2077  dca[1] = pos[1] - v[1];
2078 
2079  // "Unidentified" for data and MC
2080  fhDCA_XY->Fill(trackPt, dca[0]);
2081  fhDCA_Z->Fill(trackPt, dca[1]);
2082 
2083  // "Identified" for MC
2084  if (gentrack && mcID != AliPID::kUnknown) {
2085  // MC
2086  if (gentrack->IsPhysicalPrimary()) {
2087  fhDCA_XY_prim_MCID[mcID]->Fill(trackPt, dca[0]);
2088  fhDCA_Z_prim_MCID[mcID]->Fill(trackPt, dca[1]);
2089  }
2090  else {
2091  fhDCA_XY_sec_MCID[mcID]->Fill(trackPt, dca[0]);
2092  fhDCA_Z_sec_MCID[mcID]->Fill(trackPt, dca[1]);
2093  }
2094  }
2095 }
TString fNonStdFile
where we take the jets from can be input or output AOD
Int_t pdg
Double_t Area() const
Definition: AliEmcalJet.h:130
TH1F * fh1EvtCent
number of reconstructed tracks after cuts
Double_t GetRhoVal() const
double Double_t
Definition: External.C:58
void NormalizeJetArea(Double_t jetParameter)
static void SetProperties(TH1 *h, const char *x, const char *y)
TH1F * fh1nGenJets
number of jets from reconstructed tracks per event
Definition: External.C:236
Bool_t FillEfficiencyContainer(const Double_t *values, EffSteps step, Double_t weight=1.0)
AliJetContainer * GetJetContainer(Int_t i=0) const
Double_t Eta() const
Definition: AliEmcalJet.h:121
virtual void AnalyseJetTrack(AliVTrack *track, AliEmcalJet *jet, AliAnalysisTaskPID *task, const Bool_t *trackRejectedByTask, Double_t centPercent, AliMCParticleContainer *mcParticleContainer=0x0)
static Bool_t TPCCutMIGeo(const AliVTrack *track, const AliVEvent *evt, TTreeStream *streamer=0x0)
Bool_t FillGenJets(Double_t centralityPercentile, Double_t jetPt, Double_t norm=-1.)
TH2F * fhDCA_XY_sec_MCID[AliPID::kSPECIES]
DCA Z for all rec. prim. particles sorted by MC-ID.
Double_t Py() const
Definition: AliEmcalJet.h:107
const AliTrackIterableContainer accepted() const
Double_t Phi() const
Definition: AliEmcalJet.h:117
Double_t GetEtaAbsCutUp() const
Container with name, TClonesArray and cuts for particles.
AliAnalysisTaskPID ** fJetUEPIDtask
Pointer to tasks for jet PID spectra.
Int_t fNumInclusivePIDtasks
Object to use analysis utils like pile-up rejection.
void SetArea(Double_t a)
Definition: AliEmcalJet.h:256
static Double_t GetMCStrangenessFactorCMS(Int_t motherPDG, Double_t motherGenPt)
virtual Bool_t GetUseTPCCutMIGeo() const
virtual Bool_t GetUseTPCnclCut() const
Bool_t IsInAcceptedEtaRange(Double_t etaAbs) const
Int_t nCentBins
virtual Bool_t GetIsPileUp(AliVEvent *event, PileUpRejectionType pileUpRejection=kPileUpRejectionClass) const
Double_t GetDistanceJetTrack(const AliEmcalJet *jet, const AliVParticle *track, Bool_t useInternalFlag=kTRUE) const
Double_t GetParticleEtaMin() const
Bool_t IsSecondaryWithStrangeMotherMC(AliAODMCParticle *part, AliMCParticleContainer *mcParticleContainer)
TString part
use mixed event to constrain combinatorial background
Definition: InvMassFit.C:52
Container for particles within the EMCAL framework.
virtual TList * GetUEJetsWithRandomConeMethod(AliJetContainer *jetContainer, Double_t coneRadius, Double_t maxEtaTrack)
Double_t Px() const
Definition: AliEmcalJet.h:106
Bool_t FillPythiaTrials(Double_t avgTrials)
Bool_t FillGeneratedYield(const Double_t *values, Double_t weight=1.0)
const Int_t nPtBins
Double_t GetMCStrangenessFactorCMS(AliAODMCParticle *daughter, AliMCParticleContainer *mcParticleContainer) const
virtual Bool_t IsRCJCOverlap(const AliVParticle *part, Double_t dDistance) const
static void GetJetTrackObservables(Double_t trackPt, Double_t jetPt, Double_t &z, Double_t &xi, Bool_t storeXi=kTRUE)
int Int_t
Definition: External.C:63
static Int_t EventClass(Bool_t bSet=kFALSE, Int_t iNew=0)
TH2F * fhDCA_Z_sec_MCID[AliPID::kSPECIES]
DCA XY for all rec. sec. particles sorted by MC-ID.
unsigned int UInt_t
Definition: External.C:33
float Float_t
Definition: External.C:68
Bool_t GetStoreCharge() const
Int_t method
Bool_t ProcessTrack(const AliVTrack *track, Int_t particlePDGcode, Double_t centralityPercentile, Double_t jetPt, Bool_t isMBSelected=kFALSE, Bool_t isMultSelected=kTRUE, Bool_t storeXi=kTRUE, Double_t radialDistanceToJet=-1, Double_t jT=-1)
void FillUEDensity(Double_t cent, Double_t UEpt)
Bool_t fUseInclusivePIDtask
Pointer to tasks for jet UE PID spectra.
TRandom3 * fRandom
DCA Z for all rec. sec. particles sorted by MC-ID.
Int_t GetIndexOfChargeAxisGenYield() const
virtual void ConfigureTaskForCurrentEvent(AliVEvent *event)
static Bool_t TPCnclCut(const AliVTrack *track)
Double_t GetPerpendicularMomentumTrackJet(const AliEmcalJet *jet, const AliVParticle *track, Bool_t useInternalFlag=kTRUE) const
AliAnalysisTaskPID ** fJetPIDtask
Pointer to tasks for inclusive PID spectra.
virtual void FillDCA(AliVTrack *track, AliMCParticleContainer *mcParticleContainer)
TH2F * fhJetPtRefMultEta8
Jet pT vs. reference multiplicity (|eta|<0.5)
AliRhoParameter * GetRhoParameter()
Double_t Pt() const
Definition: AliEmcalJet.h:109
virtual void PerformJetMonteCarloAnalysisGeneratedYield(AliEmcalJet *jet, AliVParticle *trackVP, AliAnalysisTaskPID *task, Double_t centPercent, AliJetContainer *mcJetContainer=0x0)
virtual AliEmcalJet * GetPerpendicularCone(AliEmcalJet *processedJet, Double_t perpAngle) const
TH2F * fhJetPtMultPercent
Jet pT vs. reference multiplicity (|eta|<0.8)
AliMCParticleContainer * GetMCParticleContainer(Int_t i=0) const
TH2F * fhDCA_XY_prim_MCID[AliPID::kSPECIES]
Jet pT vs. multiplicity percentile (usually V0M)
virtual TList * GetUEJetsWithPerpendicularConeMethod(AliJetContainer *jetContainer)
TH1F * fh1nRecJetsCuts
Number events before and after the cut on MC pT hard.
Double_t GetParticleEtaMax() const
AliTrackContainer * GetTrackContainer(Int_t i=0) const
AliAODExtension * fAODExtension
AOD event with jet branch (case we have AOD both in input and output)
Base task in the EMCAL jet framework.
Bool_t IsParticleInCone(const AliVParticle *part1, const AliVParticle *part2, Double_t dRMax) const
Represent a jet reconstructed using the EMCal jet framework.
Definition: AliEmcalJet.h:51
static Int_t PDGtoMCID(Int_t pdg)
static Bool_t PythiaInfoFromFile(const char *currFile, Float_t &fXsec, Float_t &fTrials)
Double_t Pz() const
Definition: AliEmcalJet.h:108
void FillJetArea(Double_t cent, Double_t area)
virtual AliEmcalJet * GetRandomCone(AliEmcalJet *processedJet, Double_t dEtaConeMax, Double_t dDistance) const
const AliParticleIterableContainer accepted() const
const char Option_t
Definition: External.C:48
virtual TList * GetTracksInCone(const AliEmcalJet *jet, AliParticleContainer *particleContainer=0x0) const
const AliJetIterableContainer accepted() const
bool Bool_t
Definition: External.C:53
Bool_t fOnlyLeadingJets
TRandom3 for background estimation.
TH2F * fhDCA_XY
number of jets from generated tracks per event
TH2F * fhDCA_Z_prim_MCID[AliPID::kSPECIES]
DCA XY for all rec. prim. particles sorted by MC-ID.
Container for MC-true particles within the EMCAL framework.
virtual AliAODMCParticle * GetMCParticleWithLabel(Int_t lab) const
Double_t M() const
Definition: AliEmcalJet.h:120
Bool_t FillRecJets(Double_t centralityPercentile, Double_t jetPt, Double_t norm=-1.)
const AliMCParticleIterableContainer accepted() const
Container for jet within the EMCAL jet framework.
Definition: External.C:196
Bool_t FillCutHisto(Double_t value, CutHistoType type)
Bool_t IncrementEventCounter(Double_t centralityPercentile, EventCounterType type)
TList * OpenFile(const char *fname)
Definition: DrawAnaELoss.C:65
Bool_t FillPtResolution(Int_t mcID, const Double_t *values)