AliPhysics  a6017e1 (a6017e1)
AliAnalysisTaskJetExtractor.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-2018, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Author: R. Haake. *
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 <iostream>
17 #include <cstdlib>
18 #include <algorithm>
19 #include <vector>
20 #include <TClonesArray.h>
21 #include <TH1F.h>
22 #include <TH2F.h>
23 #include <TH3F.h>
24 #include <TProfile.h>
25 #include <TTree.h>
26 #include <TList.h>
27 #include <TLorentzVector.h>
28 #include <TNamed.h>
29 
30 #include "AliMCEvent.h"
31 #include "AliESDVertex.h"
32 #include "AliAODVertex.h"
33 #include "AliAODPid.h"
34 
35 #include "AliRhoParameter.h"
36 
37 #include "AliVTrack.h"
38 #include "AliVHeader.h"
39 #include "AliEmcalJet.h"
40 #include "AliLog.h"
41 #include "AliJetContainer.h"
42 #include "AliParticleContainer.h"
43 #include "AliAODTrack.h"
44 #include "AliVParticle.h"
45 #include "TRandom3.h"
46 #include "AliEmcalPythiaInfo.h"
48 #include "AliAnalysisManager.h"
49 #include "AliYAMLConfiguration.h"
50 
51 #include "AliGenHepMCEventHeader.h"
52 #include "AliGenHijingEventHeader.h"
53 #include "AliHFJetsTaggingVertex.h"
54 #include "AliRDHFJetsCutsVertex.h"
55 
57 
59 ClassImp(AliEmcalJetTree)
61 
65 
66 //________________________________________________________________________
67 AliEmcalJetTree::AliEmcalJetTree() : TNamed("CustomTree", "CustomTree"), fJetTree(0), fInitialized(0), fSaveEventProperties(0), fSaveConstituents(0), fSaveConstituentsIP(0), fSaveConstituentPID(0), fSaveJetShapes(0), fSaveMCInformation(0), fSaveSecondaryVertices(0), fSaveTriggerTracks(0), fExtractionPercentages(), fExtractionPercentagePtBins(), fExtractionJetTypes_HM(), fExtractionJetTypes_PM()
68 {
69  // For these arrays, we need to reserve memory
70  fBuffer_Const_Pt = new Float_t[kMaxNumConstituents];
71  fBuffer_Const_Eta = new Float_t[kMaxNumConstituents];
72  fBuffer_Const_Phi = new Float_t[kMaxNumConstituents];
73  fBuffer_Const_Charge = new Float_t[kMaxNumConstituents];
74  fBuffer_Const_Label = new Int_t [kMaxNumConstituents];
75  fBuffer_Const_ProdVtx_X = new Float_t[kMaxNumConstituents];
76  fBuffer_Const_ProdVtx_Y = new Float_t[kMaxNumConstituents];
77  fBuffer_Const_ProdVtx_Z = new Float_t[kMaxNumConstituents];
78 
79 }
80 
81 //________________________________________________________________________
82 AliEmcalJetTree::AliEmcalJetTree(const char* name) : TNamed(name, name), fJetTree(0), fInitialized(0), fSaveEventProperties(0), fSaveConstituents(0), fSaveConstituentsIP(0), fSaveConstituentPID(0), fSaveJetShapes(0), fSaveMCInformation(0), fSaveSecondaryVertices(0), fSaveTriggerTracks(0), fExtractionPercentages(), fExtractionPercentagePtBins(), fExtractionJetTypes_HM(), fExtractionJetTypes_PM()
83 {
84  // For these arrays, we need to reserve memory
93 }
94 
95 
96 //________________________________________________________________________
97 Bool_t AliEmcalJetTree::AddJetToTree(AliEmcalJet* jet, Float_t bgrdDensity, Float_t vertexX, Float_t vertexY, Float_t vertexZ, Float_t centrality, Long64_t eventID, Float_t magField,
98  AliParticleContainer* trackCont, Int_t motherParton, Int_t motherHadron, Int_t partonInitialCollision, Float_t matchedPt, Float_t truePtFraction, Float_t ptHard, Float_t eventWeight, Float_t impactParameter,
99  Float_t* trackPID_ITS, Float_t* trackPID_TPC, Float_t* trackPID_TOF, Float_t* trackPID_TRD, Short_t* trackPID_Reco, Int_t* trackPID_Truth,
100  Int_t numTriggerTracks, Float_t* triggerTrackPt, Float_t* triggerTrackDeltaEta, Float_t* triggerTrackDeltaPhi,
101  Float_t* trackIP_d0, Float_t* trackIP_z0, Float_t* trackIP_d0cov, Float_t* trackIP_z0cov,
102  Int_t numSecVertices, Float_t* secVtx_X, Float_t* secVtx_Y, Float_t* secVtx_Z, Float_t* secVtx_Mass, Float_t* secVtx_Lxy, Float_t* secVtx_SigmaLxy, Float_t* secVtx_Chi2, Float_t* secVtx_Dispersion)
103 {
104  // ############
105  // Approach: Fill buffer and then add to tree
106  //
107  if(!fInitialized)
108  AliFatal("Tree is not initialized.");
109 
110  fBuffer_JetPt = jet->Pt() - bgrdDensity*jet->Area();
111 
112  // Check if jet type is contained in extraction list
113  if( (fExtractionJetTypes_PM.size() || fExtractionJetTypes_HM.size()) &&
114  (std::find(fExtractionJetTypes_PM.begin(), fExtractionJetTypes_PM.end(), motherParton) == fExtractionJetTypes_PM.end()) &&
115  (std::find(fExtractionJetTypes_HM.begin(), fExtractionJetTypes_HM.end(), motherHadron) == fExtractionJetTypes_HM.end()) )
116  return kFALSE;
117 
118  // Check acceptance percentage for the given jet and discard statistically on demand
119  Bool_t inPtRange = kFALSE;
120  for(size_t i = 0; i<fExtractionPercentages.size(); i++)
121  {
123  {
124  inPtRange = kTRUE;
125  if(gRandom->Rndm() >= fExtractionPercentages[i])
126  return kFALSE;
127  break;
128  }
129  }
130  if(!inPtRange && fExtractionPercentagePtBins.size()) // only discard jet if fExtractionPercentagePtBins was given
131  return kFALSE;
132 
133  // Set basic properties
134  fBuffer_JetEta = jet->Eta();
135  fBuffer_JetPhi = jet->Phi();
136  fBuffer_JetArea = jet->Area();
137 
138  // Set event properties
140  {
141  fBuffer_Event_BackgroundDensity = bgrdDensity;
142  fBuffer_Event_Vertex_X = vertexX;
143  fBuffer_Event_Vertex_Y = vertexY;
144  fBuffer_Event_Vertex_Z = vertexZ;
146  fBuffer_Event_ID = eventID;
147  fBuffer_Event_MagneticField = magField;
149  {
150  fBuffer_Event_PtHard = ptHard;
151  fBuffer_Event_Weight = eventWeight;
152  fBuffer_Event_ImpactParameter = impactParameter;
153  }
154  }
155 
156  // Extract basic constituent properties directly from AliEmcalJet object
158  if(trackCont && (fSaveConstituents || fSaveConstituentsIP))
159  for(Int_t i = 0; i < jet->GetNumberOfTracks(); i++)
160  {
161  AliVParticle* particle = static_cast<AliVParticle*>(jet->TrackAt(i, trackCont->GetArray()));
162  if(!particle) continue;
163 
165  {
166  fBuffer_Const_Pt[fBuffer_NumConstituents] = particle->Pt();
167  fBuffer_Const_Eta[fBuffer_NumConstituents] = particle->Eta();
168  fBuffer_Const_Phi[fBuffer_NumConstituents] = particle->Phi();
169  fBuffer_Const_Charge[fBuffer_NumConstituents] = particle->Charge();
170  fBuffer_Const_Label[fBuffer_NumConstituents] = particle->GetLabel();
171  }
173  {
177  }
179  }
180 
181 
182  // Set constituent arrays for impact parameters
184  {
185  fJetTree->SetBranchAddress("Jet_Const_CovIPd", trackIP_d0cov);
186  fJetTree->SetBranchAddress("Jet_Const_CovIPz", trackIP_z0cov);
187  fJetTree->SetBranchAddress("Jet_Const_IPd", trackIP_d0);
188  fJetTree->SetBranchAddress("Jet_Const_IPz", trackIP_z0);
189  }
190  // Set constituent arrays for PID parameters
192  {
193  fJetTree->SetBranchAddress("Jet_Const_PID_ITS", trackPID_ITS);
194  fJetTree->SetBranchAddress("Jet_Const_PID_TPC", trackPID_TPC);
195  fJetTree->SetBranchAddress("Jet_Const_PID_TOF", trackPID_TOF);
196  fJetTree->SetBranchAddress("Jet_Const_PID_TRD", trackPID_TRD);
197  fJetTree->SetBranchAddress("Jet_Const_PID_Reconstructed", trackPID_Reco);
199  fJetTree->SetBranchAddress("Jet_Const_PID_Truth", trackPID_Truth);
200  }
201 
202  // Set secondary vertices arrays
204  {
205  fBuffer_NumSecVertices = numSecVertices;
206  fJetTree->SetBranchAddress("Jet_SecVtx_X", secVtx_X);
207  fJetTree->SetBranchAddress("Jet_SecVtx_Y", secVtx_Y);
208  fJetTree->SetBranchAddress("Jet_SecVtx_Z", secVtx_Z);
209  fJetTree->SetBranchAddress("Jet_SecVtx_Mass", secVtx_Mass);
210  fJetTree->SetBranchAddress("Jet_SecVtx_Lxy", secVtx_Lxy);
211  fJetTree->SetBranchAddress("Jet_SecVtx_SigmaLxy", secVtx_SigmaLxy);
212  fJetTree->SetBranchAddress("Jet_SecVtx_Chi2", secVtx_Chi2);
213  fJetTree->SetBranchAddress("Jet_SecVtx_Dispersion", secVtx_Dispersion);
214  }
215 
216  // Set jet shape observables
217  fBuffer_Shape_Mass = jet->M();
218 
219  // Set Monte Carlo information
221  {
222  fBuffer_Jet_MC_MotherParton = motherParton;
223  fBuffer_Jet_MC_MotherHadron = motherHadron;
224  fBuffer_Jet_MC_MotherIC = partonInitialCollision;
225  fBuffer_Jet_MC_MatchedPt = matchedPt;
226  fBuffer_Jet_MC_TruePtFraction = truePtFraction;
227  }
228 
230  {
231  fBuffer_NumTriggerTracks = numTriggerTracks;
232  fJetTree->SetBranchAddress("Jet_TriggerTrack_Pt", triggerTrackPt);
233  fJetTree->SetBranchAddress("Jet_TriggerTrack_dEta", triggerTrackDeltaEta);
234  fJetTree->SetBranchAddress("Jet_TriggerTrack_dPhi", triggerTrackDeltaPhi);
235  }
236 
237  // Add buffered jet to tree
238  fJetTree->Fill();
239 
240  return kTRUE;
241 }
242 
243 
244 //________________________________________________________________________
246 {
247  // Create the tree with active branches
248 
249  void* dummy = 0; // for some branches, do not need a buffer pointer yet
250 
251 
252  if(fInitialized)
253  AliFatal("Tree is already initialized.");
254 
255  // ### Prepare the jet tree
256  fJetTree = new TTree(Form("JetTree_%s", GetName()), "");
257 
258  fJetTree->Branch("Jet_Pt",&fBuffer_JetPt,"Jet_Pt/F");
259  fJetTree->Branch("Jet_Phi",&fBuffer_JetPhi,"Jet_Phi/F");
260  fJetTree->Branch("Jet_Eta",&fBuffer_JetEta,"Jet_Eta/F");
261  fJetTree->Branch("Jet_Area",&fBuffer_JetArea,"Jet_Area/F");
262  fJetTree->Branch("Jet_NumConstituents",&fBuffer_NumConstituents,"Jet_NumConstituents/I");
263 
265  {
266  fJetTree->Branch("Event_BackgroundDensity",&fBuffer_Event_BackgroundDensity,"Event_BackgroundDensity/F");
267  fJetTree->Branch("Event_Vertex_X",&fBuffer_Event_Vertex_X,"Event_Vertex_X/F");
268  fJetTree->Branch("Event_Vertex_Y",&fBuffer_Event_Vertex_Y,"Event_Vertex_Y/F");
269  fJetTree->Branch("Event_Vertex_Z",&fBuffer_Event_Vertex_Z,"Event_Vertex_Z/F");
270  fJetTree->Branch("Event_Centrality",&fBuffer_Event_Centrality,"Event_Centrality/F");
271  fJetTree->Branch("Event_ID",&fBuffer_Event_ID,"Event_ID/L");
272  fJetTree->Branch("Event_MagneticField",&fBuffer_Event_MagneticField,"Event_MagneticField/F");
273 
275  {
276  fJetTree->Branch("Event_PtHard",&fBuffer_Event_PtHard,"Event_PtHard/F");
277  fJetTree->Branch("Event_Weight",&fBuffer_Event_Weight,"Event_Weight/F");
278  fJetTree->Branch("Event_ImpactParameter",&fBuffer_Event_ImpactParameter,"Event_ImpactParameter/F");
279  }
280  }
281 
283  {
284  fJetTree->Branch("Jet_Const_Pt",fBuffer_Const_Pt,"Jet_Const_Pt[Jet_NumConstituents]/F");
285  fJetTree->Branch("Jet_Const_Phi",fBuffer_Const_Phi,"Jet_Const_Phi[Jet_NumConstituents]/F");
286  fJetTree->Branch("Jet_Const_Eta",fBuffer_Const_Eta,"Jet_Const_Eta[Jet_NumConstituents]/F");
287  fJetTree->Branch("Jet_Const_Charge",fBuffer_Const_Charge,"Jet_Const_Charge[Jet_NumConstituents]/F");
289  fJetTree->Branch("Jet_Const_Label",fBuffer_Const_Label,"Jet_Const_Label[Jet_NumConstituents]/I");
290  }
291 
293  {
294  fJetTree->Branch("Jet_Const_IPd",&dummy,"Jet_Const_IPd[Jet_NumConstituents]/F");
295  fJetTree->Branch("Jet_Const_IPz",&dummy,"Jet_Const_IPz[Jet_NumConstituents]/F");
296  fJetTree->Branch("Jet_Const_CovIPd",&dummy,"Jet_Const_CovIPd[Jet_NumConstituents]/F");
297  fJetTree->Branch("Jet_Const_CovIPz",&dummy,"Jet_Const_CovIPz[Jet_NumConstituents]/F");
298 
299  fJetTree->Branch("Jet_Const_ProdVtx_X",fBuffer_Const_ProdVtx_X,"Jet_Const_ProdVtx_X[Jet_NumConstituents]/F");
300  fJetTree->Branch("Jet_Const_ProdVtx_Y",fBuffer_Const_ProdVtx_Y,"Jet_Const_ProdVtx_Y[Jet_NumConstituents]/F");
301  fJetTree->Branch("Jet_Const_ProdVtx_Z",fBuffer_Const_ProdVtx_Z,"Jet_Const_ProdVtx_Z[Jet_NumConstituents]/F");
302  }
303 
305  {
306  fJetTree->Branch("Jet_Const_PID_ITS",&dummy,"Jet_Const_PID_ITS[Jet_NumConstituents]/F");
307  fJetTree->Branch("Jet_Const_PID_TPC",&dummy,"Jet_Const_PID_TPC[Jet_NumConstituents]/F");
308  fJetTree->Branch("Jet_Const_PID_TOF",&dummy,"Jet_Const_PID_TOF[Jet_NumConstituents]/F");
309  fJetTree->Branch("Jet_Const_PID_TRD",&dummy,"Jet_Const_PID_TRD[Jet_NumConstituents]/F");
310 
311  fJetTree->Branch("Jet_Const_PID_Reconstructed",&dummy,"Jet_Const_PID_Reconstructed[Jet_NumConstituents]/S");
313  fJetTree->Branch("Jet_Const_PID_Truth",&dummy,"Jet_Const_PID_Truth[Jet_NumConstituents]/I");
314  }
315 
316  if(fSaveJetShapes)
317  {
318  fJetTree->Branch("Jet_Shape_Mass",&fBuffer_Shape_Mass,"Jet_Shape_Mass/F");
319  }
320 
322  {
323  fJetTree->Branch("Jet_MC_MotherParton",&fBuffer_Jet_MC_MotherParton,"Jet_MC_MotherParton/I");
324  fJetTree->Branch("Jet_MC_MotherHadron",&fBuffer_Jet_MC_MotherHadron,"Jet_MC_MotherHadron/I");
325  fJetTree->Branch("Jet_MC_MotherIC",&fBuffer_Jet_MC_MotherIC,"Jet_MC_MotherIC/I");
326  fJetTree->Branch("Jet_MC_MatchedPt",&fBuffer_Jet_MC_MatchedPt,"Jet_MC_MatchedPt/F");
327  fJetTree->Branch("Jet_MC_TruePtFraction",&fBuffer_Jet_MC_TruePtFraction,"Jet_MC_TruePtFraction/F");
328  }
329 
331  {
332  fJetTree->Branch("Jet_NumSecVertices",&fBuffer_NumSecVertices,"Jet_NumSecVertices/I");
333 
334  fJetTree->Branch("Jet_SecVtx_X",&dummy,"Jet_SecVtx_X[Jet_NumSecVertices]/F");
335  fJetTree->Branch("Jet_SecVtx_Y",&dummy,"Jet_SecVtx_Y[Jet_NumSecVertices]/F");
336  fJetTree->Branch("Jet_SecVtx_Z",&dummy,"Jet_SecVtx_Z[Jet_NumSecVertices]/F");
337  fJetTree->Branch("Jet_SecVtx_Mass",&dummy,"Jet_SecVtx_Mass[Jet_NumSecVertices]/F");
338  fJetTree->Branch("Jet_SecVtx_Lxy",&dummy,"Jet_SecVtx_Lxy[Jet_NumSecVertices]/F");
339  fJetTree->Branch("Jet_SecVtx_SigmaLxy",&dummy,"Jet_SecVtx_SigmaLxy[Jet_NumSecVertices]/F");
340  fJetTree->Branch("Jet_SecVtx_Chi2",&dummy,"Jet_SecVtx_Chi2[Jet_NumSecVertices]/F");
341  fJetTree->Branch("Jet_SecVtx_Dispersion",&dummy,"Jet_SecVtx_Dispersion[Jet_NumSecVertices]/F");
342  }
343 
344  // Trigger track requirement active -> Save trigger track
346  {
347  fJetTree->Branch("Jet_NumTriggerTracks",&fBuffer_NumTriggerTracks,"Jet_NumTriggerTracks/I");
348  fJetTree->Branch("Jet_TriggerTrack_Pt",&dummy,"Jet_TriggerTrack_Pt[Jet_NumTriggerTracks]/F");
349  fJetTree->Branch("Jet_TriggerTrack_dEta",&dummy,"Jet_TriggerTrack_dEta[Jet_NumTriggerTracks]/F");
350  fJetTree->Branch("Jet_TriggerTrack_dPhi",&dummy,"Jet_TriggerTrack_dPhi[Jet_NumTriggerTracks]/F");
351  }
352 
353  fInitialized = kTRUE;
354 }
355 
356 //________________________________________________________________________
358  AliAnalysisTaskEmcalJet("AliAnalysisTaskJetExtractor", kTRUE),
359  fJetTree(0),
360  fVtxTagger(0),
361  fHadronMatchingRadius(0.4),
362  fTrueJetMatchingRadius(0.4),
363  fSecondaryVertexMaxChi2(1e10),
364  fSecondaryVertexMaxDispersion(0.05),
365  fCalculateSecondaryVertices(kTRUE),
366  fVertexerCuts(0),
367  fSetEmcalJetFlavour(0),
368  fEventPercentage(1.0),
369  fTruthMinLabel(0),
370  fTruthMaxLabel(100000),
371  fSaveTrackPDGCode(kTRUE),
372  fEventWeight(0.),
373  fImpactParameter(0.),
374  fEventCut_TriggerTrackMinPt(0),
375  fEventCut_TriggerTrackMaxPt(0),
376  fEventCut_TriggerTrackMinLabel(-9999999),
377  fEventCut_TriggerTrackMaxLabel(+9999999),
378  fJetsCont(0),
379  fTracksCont(0),
380  fTruthParticleArray(0),
381  fTruthJetsArrayName(""),
382  fTruthJetsRhoName(""),
383  fTruthParticleArrayName("mcparticles"),
384  fTriggerTracks_Pt(),
385  fTriggerTracks_Eta(),
386  fTriggerTracks_Phi()
387 {
389  fJetTree = new AliEmcalJetTree(GetName());
392  fJetTree->SetSaveJetShapes(kTRUE);
393  DefineOutput(2, TTree::Class());
394 }
395 
396 //________________________________________________________________________
398  AliAnalysisTaskEmcalJet(name, kTRUE),
399  fJetTree(0),
400  fVtxTagger(0),
406  fVertexerCuts(0),
408  fEventPercentage(1.0),
409  fTruthMinLabel(0),
410  fTruthMaxLabel(100000),
411  fSaveTrackPDGCode(kTRUE),
412  fEventWeight(0.),
413  fImpactParameter(0.),
418  fJetsCont(0),
419  fTracksCont(0),
422  fTruthJetsRhoName(""),
423  fTruthParticleArrayName("mcparticles"),
427 {
429  fJetTree = new AliEmcalJetTree(GetName());
432  fJetTree->SetSaveJetShapes(kTRUE);
433  DefineOutput(2, TTree::Class());
434 }
435 
436 //________________________________________________________________________
438 {
439  // Destructor.
440 }
441 
442 //________________________________________________________________________
444 {
446 
447  // ### Basic container settings
449  if(!fJetsCont)
450  AliFatal("Jet input container not found!");
451  fJetsCont->PrintCuts();
453  if(!fTracksCont)
454  AliFatal("Particle input container not found attached to jets!");
455 
456  // Activate saving of trigger tracks if this is demanded
459 
460  // ### Initialize the jet tree (settings must all be given at this stage)
462  OpenFile(2);
463  PostData(2, fJetTree->GetTreePointer());
464 
465  PrintConfig();
466 
467  // ### Add control histograms (already some created in base task)
468  AddHistogram2D<TH2D>("hTrackCount", "Number of tracks in acceptance vs. centrality", "COLZ", 500, 0., 5000., 100, 0, 100, "N tracks","Centrality", "dN^{Events}/dN^{Tracks}");
469  AddHistogram2D<TH2D>("hBackgroundPt", "Background p_{T} distribution", "", 150, 0., 150., 100, 0, 100, "Background p_{T} (GeV/c)", "Centrality", "dN^{Events}/dp_{T}");
470 
471  AddHistogram2D<TH2D>("hJetPtRaw", "Jets p_{T} distribution (raw)", "COLZ", 300, 0., 300., 100, 0, 100, "p_{T, jet} (GeV/c)", "Centrality", "dN^{Jets}/dp_{T}");
472  AddHistogram2D<TH2D>("hJetPt", "Jets p_{T} distribution (background subtracted)", "COLZ", 400, -100., 300., 100, 0, 100, "p_{T, jet} (GeV/c)", "Centrality", "dN^{Jets}/dp_{T}");
473  AddHistogram2D<TH2D>("hJetPtExtracted", "Extracted jets p_{T} distribution (background subtracted)", "COLZ", 400, -100., 300., 100, 0, 100, "p_{T, jet} (GeV/c)", "Centrality", "dN^{Jets}/dp_{T}");
474  AddHistogram2D<TH2D>("hJetPhiEta", "Jet angular distribution #phi/#eta", "COLZ", 180, 0., 2*TMath::Pi(), 100, -2.5, 2.5, "#phi", "#eta", "dN^{Jets}/d#phi d#eta");
475  AddHistogram2D<TH2D>("hJetArea", "Jet area", "COLZ", 200, 0., 2., 100, 0, 100, "Jet A", "Centrality", "dN^{Jets}/dA");
476 
477  AddHistogram2D<TH2D>("hDeltaPt", "#delta p_{T} distribution", "", 400, -100., 300., 100, 0, 100, "p_{T, cone} (GeV/c)", "Centrality", "dN^{Tracks}/dp_{T}");
478 
479  AddHistogram2D<TH2D>("hConstituentPt", "Jet constituent p_{T} distribution", "COLZ", 400, 0., 300., 100, 0, 100, "p_{T, const} (GeV/c)", "Centrality", "dN^{Const}/dp_{T}");
480  AddHistogram2D<TH2D>("hConstituentPhiEta", "Jet constituent relative #phi/#eta distribution", "COLZ", 120, -0.6, 0.6, 120, -0.6, 0.6, "#Delta#phi", "#Delta#eta", "dN^{Const}/d#phi d#eta");
481 
482  AddHistogram1D<TH1D>("hExtractionPercentage", "Extracted jets p_{T} distribution (background subtracted)", "e1p", 400, -100., 300., "p_{T, jet} (GeV/c)", "Extracted percentage");
483 
484  // Track QA plots
485  AddHistogram2D<TH2D>("hTrackPt", "Tracks p_{T} distribution", "", 300, 0., 300., 100, 0, 100, "p_{T} (GeV/c)", "Centrality", "dN^{Tracks}/dp_{T}");
486  AddHistogram2D<TH2D>("hTrackPhi", "Track angular distribution in #phi", "LEGO2", 180, 0., 2*TMath::Pi(), 100, 0, 100, "#phi", "Centrality", "dN^{Tracks}/(d#phi)");
487  AddHistogram2D<TH2D>("hTrackEta", "Track angular distribution in #eta", "LEGO2", 100, -2.5, 2.5, 100, 0, 100, "#eta", "Centrality", "dN^{Tracks}/(d#eta)");
488  AddHistogram2D<TH2D>("hTrackPhiEta", "Track angular distribution #phi/#eta", "COLZ", 180, 0., 2*TMath::Pi(), 100, -2.5, 2.5, "#phi", "#eta", "dN^{Tracks}/d#phi d#eta");
489 
490  AddHistogram2D<TH2D>("hTrackEtaPt", "Track angular distribution in #eta vs. p_{T}", "LEGO2", 100, -2.5, 2.5, 300, 0., 300., "#eta", "p_{T} (GeV/c)", "dN^{Tracks}/(d#eta dp_{T})");
491  AddHistogram2D<TH2D>("hTrackPhiPt", "Track angular distribution in #phi vs. p_{T}", "LEGO2", 180, 0, 2*TMath::Pi(), 300, 0., 300., "#phi", "p_{T} (GeV/c)", "dN^{Tracks}/(d#phi dp_{T})");
492 }
493 
494 
495 //________________________________________________________________________
497 {
499  if (fTruthParticleArrayName != "")
500  fTruthParticleArray = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fTruthParticleArrayName.Data()));
501 
502  // ### Prepare vertexer
504  {
505  if(!fVertexerCuts)
506  AliFatal("VertexerCuts not given but secondary vertex calculation turned on.");
507  fVtxTagger = new AliHFJetsTaggingVertex();
508  fVtxTagger->SetCuts(fVertexerCuts);
509  }
510 
511  // ### Save tree extraction percentage to histogram
512  std::vector<Float_t> extractionPtBins = fJetTree->GetExtractionPercentagePtBins();
513  std::vector<Float_t> extractionPercentages = fJetTree->GetExtractionPercentages();
514 
515  for(Int_t i=0; i<static_cast<Int_t>(extractionPercentages.size()); i++)
516  {
517  Double_t percentage = extractionPercentages[i];
518  for(Int_t pt=static_cast<Int_t>(extractionPtBins[i*2]); pt<static_cast<Int_t>(extractionPtBins[i*2+1]); pt++)
519  {
520  FillHistogram("hExtractionPercentage", pt, percentage);
521  }
522  }
523 
524  // ### Add HIJING/JEWEL histograms (in case we need it)
525  if(MCEvent())
526  {
527  if(dynamic_cast<AliGenHijingEventHeader*>(MCEvent()->GenEventHeader()))
528  AddHistogram1D<TH1D>("hHIJING_Ncoll", "N_{coll} from HIJING", "e1p", 1000, 0., 5000, "N_{coll}", "dN^{Events}/dN^{N_{coll}}");
529 
530  if(dynamic_cast<AliGenHepMCEventHeader*>(MCEvent()->GenEventHeader()))
531  {
532  AddHistogram1D<TH1D>("hJEWEL_NProduced", "NProduced from JEWEL", "e1p", 1000, 0., 1000, "Nproduced", "dN^{Events}/dN^{Nproduced}");
533  AddHistogram1D<TH1D>("hJEWEL_ImpactParameter", "Impact parameter from JEWEL", "e1p", 1000, 0., 1, "IP", "dN^{Events}/dN^{IP}");
534  TProfile* evweight = new TProfile("hJEWEL_EventWeight", "Event weight from JEWEL", 1, 0,1);
535  fOutput->Add(evweight);
536  }
537  }
538 
539 }
540 
541 //________________________________________________________________________
543 {
544  // ################################### EVENT SELECTION
545 
546  if(!IsTriggerTrackInEvent())
547  return kFALSE;
548  if(gRandom->Rndm() >= fEventPercentage)
549  return kFALSE;
550 
551  // ################################### EVENT PROPERTIES
553 
554  // Load vertex if possible
555  Double_t vtxX = 0;
556  Double_t vtxY = 0;
557  Double_t vtxZ = 0;
558  Long64_t eventID = 0;
559  const AliVVertex* myVertex = InputEvent()->GetPrimaryVertex();
560  if(!myVertex && MCEvent())
561  myVertex = MCEvent()->GetPrimaryVertex();
562  if(myVertex)
563  {
564  vtxX = myVertex->GetX();
565  vtxY = myVertex->GetY();
566  vtxZ = myVertex->GetZ();
567  }
568  // Get event ID from header
569  AliVHeader* eventIDHeader = InputEvent()->GetHeader();
570  if(eventIDHeader)
571  eventID = eventIDHeader->GetEventIdAsLong();
572 
573  // If available, get Ncoll from HIJING
574  if(MCEvent())
575  {
576  if(AliGenHijingEventHeader* mcHeader = dynamic_cast<AliGenHijingEventHeader*>(MCEvent()->GenEventHeader()))
577  {
578  Double_t ncoll = mcHeader->NN() + mcHeader->NNw() + mcHeader->NwN() + mcHeader->NwNw();
579  FillHistogram("hHIJING_Ncoll", ncoll);
580  }
581  if(AliGenHepMCEventHeader* mcHeader = dynamic_cast<AliGenHepMCEventHeader*>(MCEvent()->GenEventHeader()))
582  {
583  fEventWeight = mcHeader->EventWeight();
584  fImpactParameter = mcHeader->impact_parameter();
585 
586  TProfile* tmpHist = static_cast<TProfile*>(fOutput->FindObject("hJEWEL_EventWeight"));
587  tmpHist->Fill(0., fEventWeight);
588  FillHistogram("hJEWEL_NProduced", mcHeader->NProduced());
589  FillHistogram("hJEWEL_ImpactParameter", fImpactParameter);
590  }
591  }
592 
593  // ################################### MAIN JET LOOP
594  fJetsCont->ResetCurrentID();
595  while(AliEmcalJet *jet = fJetsCont->GetNextAcceptJet())
596  {
598 
599  Double_t matchedJetPt = 0;
600  Double_t truePtFraction = 0;
601  Int_t currentJetType_HM = 0;
602  Int_t currentJetType_PM = 0;
603  Int_t currentJetType_IC = 0;
605  {
606  // Get jet type from MC (hadron matching, parton matching definition - for HF jets)
607  GetJetType(jet, currentJetType_HM, currentJetType_PM, currentJetType_IC);
608  // Get true pT estimators
609  GetJetTruePt(jet, matchedJetPt, truePtFraction);
610  }
611 
612  // ### CONSTITUENT LOOP: Retrieve PID values + impact parameters
613  std::vector<Float_t> vecSigITS; std::vector<Float_t> vecSigTPC; std::vector<Float_t> vecSigTRD; std::vector<Float_t> vecSigTOF; std::vector<Short_t> vecRecoPID; std::vector<Int_t> vecTruePID;
614  std::vector<Float_t> vec_d0; std::vector<Float_t> vec_d0cov; std::vector<Float_t> vec_z0; std::vector<Float_t> vec_z0cov;
615 
617  for(Int_t i = 0; i < jet->GetNumberOfTracks(); i++)
618  {
619  AliVParticle* particle = static_cast<AliVParticle*>(jet->TrackAt(i, fTracksCont->GetArray()));
620  if(!particle) continue;
622  {
623  Float_t sigITS = 0; Float_t sigTPC = 0; Float_t sigTOF = 0; Float_t sigTRD = 0; Short_t recoPID = 0; Int_t truePID = 0;
624  AddPIDInformation(particle, sigITS, sigTPC, sigTOF, sigTRD, recoPID, truePID);
625  vecSigITS.push_back(sigITS); vecSigTPC.push_back(sigTPC); vecSigTOF.push_back(sigTOF); vecSigTRD.push_back(sigTRD); vecRecoPID.push_back(recoPID); vecTruePID.push_back(truePID);
626  }
628  {
629  Float_t d0 = 0; Float_t d0cov = 0; Float_t z0 = 0; Float_t z0cov = 0;
630  GetTrackImpactParameters(myVertex, dynamic_cast<AliAODTrack*>(particle), d0, d0cov, z0, z0cov);
631  vec_d0.push_back(d0); vec_d0cov.push_back(d0cov); vec_z0.push_back(z0); vec_z0cov.push_back(z0cov);
632  }
633  }
634 
635  // Reconstruct secondary vertices
636  std::vector<Float_t> secVtx_X; std::vector<Float_t> secVtx_Y; std::vector<Float_t> secVtx_Z; std::vector<Float_t> secVtx_Mass; std::vector<Float_t> secVtx_Lxy; std::vector<Float_t> secVtx_SigmaLxy; std::vector<Float_t> secVtx_Chi2; std::vector<Float_t> secVtx_Dispersion;
638  ReconstructSecondaryVertices(myVertex, jet, secVtx_X, secVtx_Y, secVtx_Z, secVtx_Mass, secVtx_Lxy, secVtx_SigmaLxy, secVtx_Chi2, secVtx_Dispersion);
639 
640  // Now change trigger tracks eta/phi to dEta/dPhi relative to jet
641  std::vector<Float_t> triggerTracks_dEta(fTriggerTracks_Eta);
642  std::vector<Float_t> triggerTracks_dPhi(fTriggerTracks_Phi);
643  for(UInt_t i=0; i<triggerTracks_dEta.size(); i++)
644  {
645  triggerTracks_dEta[i] = jet->Eta() - fTriggerTracks_Eta[i];
646  triggerTracks_dPhi[i] = TMath::Min(TMath::Abs(jet->Phi() - fTriggerTracks_Phi[i]), TMath::TwoPi() - TMath::Abs(jet->Phi() - fTriggerTracks_Phi[i]));
647  if( ((TMath::Abs(jet->Phi() - fTriggerTracks_Phi[i]) <= TMath::Pi()) && (jet->Phi() - fTriggerTracks_Phi[i] > 0))
648  || ((TMath::Abs(jet->Phi() - fTriggerTracks_Phi[i]) > TMath::Pi()) && (jet->Phi() - fTriggerTracks_Phi[i] < 0)) )
649  triggerTracks_dPhi[i] = -triggerTracks_dPhi[i];
650  }
651 
652  // Fill jet to tree
653  Bool_t accepted = fJetTree->AddJetToTree(jet, fJetsCont->GetRhoVal(), vtxX, vtxY, vtxZ, fCent, eventID, InputEvent()->GetMagneticField(), fTracksCont,
654  currentJetType_PM,currentJetType_HM,currentJetType_IC,matchedJetPt,truePtFraction,fPtHard,fEventWeight,fImpactParameter,
655  vecSigITS, vecSigTPC, vecSigTOF, vecSigTRD, vecRecoPID, vecTruePID,
656  fTriggerTracks_Pt, triggerTracks_dEta, triggerTracks_dPhi,
657  vec_d0, vec_z0, vec_d0cov, vec_z0cov,
658  secVtx_X, secVtx_Y, secVtx_Z, secVtx_Mass, secVtx_Lxy, secVtx_SigmaLxy, secVtx_Chi2, secVtx_Dispersion);
659 
660  if(accepted)
661  FillHistogram("hJetPtExtracted", jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), fCent);
662  }
663 
664  // ################################### PARTICLE PROPERTIES
665  fTracksCont->ResetCurrentID();
666  while(AliVTrack *track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle()))
668 
669  return kTRUE;
670 }
671 
672 //________________________________________________________________________
674 {
675  // Cut for trigger track requirement
677  {
678  // Clear vector of trigger tracks
679  fTriggerTracks_Pt.clear();
680  fTriggerTracks_Eta.clear();
681  fTriggerTracks_Phi.clear();
682 
683  // Go through all tracks and check whether trigger tracks can be found
684  fTracksCont->ResetCurrentID();
685  while(AliVTrack *track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle()))
686  {
687  if( (track->GetLabel() >= fEventCut_TriggerTrackMinLabel) && (track->GetLabel() < fEventCut_TriggerTrackMaxLabel) )
688  if( (track->Pt() >= fEventCut_TriggerTrackMinPt) && (track->Pt() < fEventCut_TriggerTrackMaxPt) )
689  {
690  fTriggerTracks_Pt.push_back(track->Pt());
691  fTriggerTracks_Eta.push_back(track->Eta());
692  fTriggerTracks_Phi.push_back(track->Phi());
693  }
694  }
695  // No particle has been found that fulfills requirement -> Do not accept event
696  if(fTriggerTracks_Pt.size() == 0)
697  return kFALSE;
698  }
699  return kTRUE;
700 }
701 
702 
703 //________________________________________________________________________
704 void AliAnalysisTaskJetExtractor::GetJetTruePt(AliEmcalJet* jet, Double_t& matchedJetPt, Double_t& truePtFraction)
705 {
706  // #################################################################################
707  // ##### METHOD 1: If fTruthJetsArrayName is set, a matching jet is searched for
708  Double_t bestMatchDeltaR = 999.;
709  if(fTruthJetsArrayName != "")
710  {
711  // "True" background
712  AliRhoParameter* rho = static_cast<AliRhoParameter*>(InputEvent()->FindListObject(fTruthJetsRhoName.Data()));
713  Double_t trueRho = 0;
714  if(rho)
715  trueRho = rho->GetVal();
716 
717  TClonesArray* truthArray = static_cast<TClonesArray*>(InputEvent()->FindListObject(Form("%s", fTruthJetsArrayName.Data())));
718 
719  // Loop over all true jets to find the best match
720  matchedJetPt = 0;
721  if(truthArray)
722  for(Int_t i=0; i<truthArray->GetEntries(); i++)
723  {
724  AliEmcalJet* truthJet = static_cast<AliEmcalJet*>(truthArray->At(i));
725  if(truthJet->Pt() < 0.15)
726  continue;
727 
728  Double_t deltaEta = (truthJet->Eta()-jet->Eta());
729  Double_t deltaPhi = TMath::Min(TMath::Abs(truthJet->Phi()-jet->Phi()),TMath::TwoPi() - TMath::Abs(truthJet->Phi()-jet->Phi()));
730  Double_t deltaR = TMath::Sqrt(deltaEta*deltaEta + deltaPhi*deltaPhi);
731 
732  // Cut jets too far away
733  if (deltaR > fTrueJetMatchingRadius)
734  continue;
735 
736  // Search for the best match
737  if(deltaR < bestMatchDeltaR)
738  {
739  bestMatchDeltaR = deltaR;
740  matchedJetPt = truthJet->Pt() - truthJet->Area()* trueRho;
741  }
742  }
743  }
744 
745  // #################################################################################
746  // ##### METHOD 2: Calculate fraction of "true" pT -- pT which is not from a toy
747  Double_t pt_nonMC = 0.;
748  Double_t pt_all = 0.;
749  truePtFraction = 0;
750 
751  for(Int_t iConst = 0; iConst < jet->GetNumberOfTracks(); iConst++)
752  {
753  // Loop over all valid jet constituents
754  AliVParticle* particle = static_cast<AliVParticle*>(jet->TrackAt(iConst, fTracksCont->GetArray()));
755  if(!particle) continue;
756  if(particle->Pt() < 1e-6) continue;
757 
758  // Particles marked w/ labels within label range are considered from toy
759  if( (particle->GetLabel() >= fTruthMinLabel) && (particle->GetLabel() < fTruthMaxLabel))
760  pt_nonMC += particle->Pt();
761  pt_all += particle->Pt();
762  }
763  if(pt_all)
764  truePtFraction = (pt_nonMC/pt_all);
765 
766 }
767 
768 
769 //________________________________________________________________________
771 {
773 
775  return;
776 
777  typeHM = 0;
778  typePM = 0;
779 
780  AliAODMCParticle* parton[2];
781  parton[0] = (AliAODMCParticle*) fVtxTagger->IsMCJetParton(fTruthParticleArray, jet, radius); // method 1 (parton matching)
782  parton[1] = (AliAODMCParticle*) fVtxTagger->IsMCJetMeson(fTruthParticleArray, jet, radius); // method 2 (hadron matching)
783 
784  if (parton[0]) {
785  Int_t pdg = TMath::Abs(parton[0]->PdgCode());
786  typePM = pdg;
787  }
788 
789  if (!parton[1])
790  {
791  // No HF jet (according to hadron matching) -- now also separate jets in udg (1) and s-jets (3)
792  if(IsStrangeJet(jet))
793  typeHM = 3;
794  else
795  typeHM = 1;
796  }
797  else {
798  Int_t pdg = TMath::Abs(parton[1]->PdgCode());
799  if(fVtxTagger->IsDMeson(pdg)) typeHM = 4;
800  else if (fVtxTagger->IsBMeson(pdg)) typeHM = 5;
801  }
802 
803  // Set flavour of AliEmcalJet object (set ith bit while i corresponds to type)
805  jet->AddFlavourTag(static_cast<Int_t>(TMath::Power(2, typeHM)));
806 
807 
808  const AliEmcalPythiaInfo* partonsInfo = GetPythiaInfo();
809  typeIC = 0;
810  if (partonsInfo)
811  {
812  // Get primary partons directions
813  Double_t parton1phi = partonsInfo->GetPartonPhi6();
814  Double_t parton1eta = partonsInfo->GetPartonEta6();
815  Double_t parton2phi = partonsInfo->GetPartonPhi7();
816  Double_t parton2eta = partonsInfo->GetPartonEta7();
817 
818 
819  Double_t delta1Eta = (parton1eta-jet->Eta());
820  Double_t delta1Phi = TMath::Min(TMath::Abs(parton1phi-jet->Phi()),TMath::TwoPi() - TMath::Abs(parton1phi-jet->Phi()));
821  Double_t delta1R = TMath::Sqrt(delta1Eta*delta1Eta + delta1Phi*delta1Phi);
822  Double_t delta2Eta = (parton2eta-jet->Eta());
823  Double_t delta2Phi = TMath::Min(TMath::Abs(parton2phi-jet->Phi()),TMath::TwoPi() - TMath::Abs(parton2phi-jet->Phi()));
824  Double_t delta2R = TMath::Sqrt(delta2Eta*delta2Eta + delta2Phi*delta2Phi);
825 
826  // Check if one of the partons if closer than matching criterion
827  Bool_t matched = (delta1R < fJetsCont->GetJetRadius()/2.) || (delta2R < fJetsCont->GetJetRadius()/2.);
828 
829  // Matching criterion fulfilled -> Set flag to closest
830  if(matched)
831  {
832  if(delta1R < delta2R)
833  typeIC = partonsInfo->GetPartonFlag6();
834  else
835  typeIC = partonsInfo->GetPartonFlag7();
836  }
837  }
838 }
839 
840 //________________________________________________________________________
841 void AliAnalysisTaskJetExtractor::GetTrackImpactParameters(const AliVVertex* vtx, AliAODTrack* track, Float_t& d0, Float_t& d0cov, Float_t& z0, Float_t& z0cov)
842 {
843  if (track)
844  {
845  Double_t d0rphiz[2],covd0[3];
846  Bool_t isDCA=track->PropagateToDCA(vtx,InputEvent()->GetMagneticField(),3.0,d0rphiz,covd0);
847  if(isDCA)
848  {
849  d0 = d0rphiz[0];
850  z0 = d0rphiz[1];
851  d0cov = covd0[0];
852  z0cov = covd0[2];
853  }
854  }
855 }
856 
857 //________________________________________________________________________
858 void AliAnalysisTaskJetExtractor::ReconstructSecondaryVertices(const AliVVertex* primVtx, const AliEmcalJet* jet, std::vector<Float_t>& secVtx_X, std::vector<Float_t>& secVtx_Y, std::vector<Float_t>& secVtx_Z, std::vector<Float_t>& secVtx_Mass, std::vector<Float_t>& secVtx_Lxy, std::vector<Float_t>& secVtx_SigmaLxy, std::vector<Float_t>& secVtx_Chi2, std::vector<Float_t>& secVtx_Dispersion)
859 {
860  if(!primVtx)
861  return;
862 
863  // Create ESD vertex from the existing AliVVertex
864  Double_t vtxPos[3] = {primVtx->GetX(), primVtx->GetY(), primVtx->GetZ()};
865  Double_t covMatrix[6] = {0};
866  primVtx->GetCovarianceMatrix(covMatrix);
867  AliESDVertex* esdVtx = new AliESDVertex(vtxPos, covMatrix, primVtx->GetChi2(), primVtx->GetNContributors());
868 
869  TClonesArray* secVertexArr = 0;
870  vctr_pair_dbl_int arrDispersion;
871  arrDispersion.reserve(5);
873  {
874  //###########################################################################
875  // ********* Calculate secondary vertices
876  // Derived from AliAnalysisTaskEmcalJetBtagSV
877  secVertexArr = new TClonesArray("AliAODVertex");
878  Int_t nDauRejCount = 0;
879  Int_t nVtx = fVtxTagger->FindVertices(jet,
880  fTracksCont->GetArray(),
881  (AliAODEvent*)InputEvent(),
882  esdVtx,
883  InputEvent()->GetMagneticField(),
884  secVertexArr,
885  0,
886  arrDispersion,
887  nDauRejCount);
888 
889 
890  if(nVtx < 0)
891  {
892  secVertexArr->Clear();
893  delete secVertexArr;
894  return;
895  }
896  //###########################################################################
897  }
898  else // Load HF vertex branch from AOD event, if possible
899  {
900  secVertexArr = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("VerticesHF"));
901  if(!secVertexArr)
902  return;
903  }
904 
905  // Loop over all potential secondary vertices
906  for(Int_t i=0; i<secVertexArr->GetEntriesFast(); i++)
907  {
908  AliAODVertex* secVtx = (AliAODVertex*)(secVertexArr->UncheckedAt(i));
910  if((strcmp(secVtx->GetParent()->ClassName(), "AliAODRecoDecayHF3Prong")))
911  continue;
912 
913  // Calculate vtx distance
914  Double_t effX = secVtx->GetX() - esdVtx->GetX();
915  Double_t effY = secVtx->GetY() - esdVtx->GetY();
916  //Double_t effZ = secVtx->GetZ() - esdVtx->GetZ();
917 
918  // ##### Vertex properties
919  // vertex dispersion
920  Double_t dispersion = arrDispersion[i].first;
921 
922  // invariant mass
923  Double_t mass = fVtxTagger->GetVertexInvariantMass(secVtx);
924 
925  // signed length
926  Double_t Lxy = TMath::Sqrt(effX*effX + effY*effY);
927  Double_t jetP[3]; jet->PxPyPz(jetP);
928  Double_t signLxy = effX * jetP[0] + effY * jetP[1];
929  if (signLxy < 0.) Lxy *= -1.;
930 
931  Double_t sigmaLxy = 0;
932  AliAODVertex* aodVtx = (AliAODVertex*)(primVtx);
933  if (aodVtx)
934  sigmaLxy = aodVtx->ErrorDistanceXYToVertex(secVtx);
935 
936  // Add secondary vertices if they fulfill the conditions
937  if( (dispersion > fSecondaryVertexMaxDispersion) || (TMath::Abs(secVtx->GetChi2perNDF()) > fSecondaryVertexMaxChi2) )
938  continue;
939 
940  secVtx_X.push_back(secVtx->GetX()); secVtx_Y.push_back(secVtx->GetY()); secVtx_Z.push_back(secVtx->GetZ()); secVtx_Chi2.push_back(secVtx->GetChi2perNDF());
941  secVtx_Dispersion.push_back(dispersion); secVtx_Mass.push_back(mass); secVtx_Lxy.push_back(Lxy); secVtx_SigmaLxy.push_back(sigmaLxy);
942  }
943 
945  {
946  secVertexArr->Clear();
947  delete secVertexArr;
948  }
949  delete esdVtx;
950 }
951 
952 //________________________________________________________________________
954 {
955  // Do hadron matching jet type tagging using mcparticles
956  // ... if not explicitly deactivated
958  {
959  for(Int_t i=0; i<fTruthParticleArray->GetEntries();i++)
960  {
961  AliAODMCParticle* part = (AliAODMCParticle*)fTruthParticleArray->At(i);
962  if(!part) continue;
963 
964  // Check if the particle has strangeness
965  Int_t absPDG = TMath::Abs(part->PdgCode());
966  if ((absPDG > 300 && absPDG < 400) || (absPDG > 3000 && absPDG < 4000))
967  {
968  // Check if particle is in a radius around the jet
969  Double_t rsquared = (part->Eta() - jet->Eta())*(part->Eta() - jet->Eta()) + (part->Phi() - jet->Phi())*(part->Phi() - jet->Phi());
971  continue;
972  else
973  return kTRUE;
974  }
975  }
976  }
977  // As fallback, the MC stack will be tried
978  else if(MCEvent() && (MCEvent()->Stack()))
979  {
980  AliStack* stack = MCEvent()->Stack();
981  // Go through the whole particle stack
982  for(Int_t i=0; i<stack->GetNtrack(); i++)
983  {
984  TParticle *part = stack->Particle(i);
985  if(!part) continue;
986 
987  // Check if the particle has strangeness
988  Int_t absPDG = TMath::Abs(part->GetPdgCode());
989  if ((absPDG > 300 && absPDG < 400) || (absPDG > 3000 && absPDG < 4000))
990  {
991  // Check if particle is in a radius around the jet
992  Double_t rsquared = (part->Eta() - jet->Eta())*(part->Eta() - jet->Eta()) + (part->Phi() - jet->Phi())*(part->Phi() - jet->Phi());
994  continue;
995  else
996  return kTRUE;
997  }
998  }
999  }
1000  return kFALSE;
1001 
1002 }
1003 //________________________________________________________________________
1005 {
1006  // ### Event control plots
1008  FillHistogram("hBackgroundPt", fJetsCont->GetRhoVal(), fCent);
1009 }
1010 
1011 //________________________________________________________________________
1013 {
1014  // ### Jet control plots
1015  FillHistogram("hJetPtRaw", jet->Pt(), fCent);
1016  FillHistogram("hJetPt", jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), fCent);
1017  FillHistogram("hJetPhiEta", jet->Phi(), jet->Eta());
1018  FillHistogram("hJetArea", jet->Area(), fCent);
1019 
1020  // ### Jet constituent plots
1021  for(Int_t i = 0; i < jet->GetNumberOfTracks(); i++)
1022  {
1023  AliVParticle* jetConst = static_cast<AliVParticle*>(jet->TrackAt(i, fTracksCont->GetArray()));
1024  if(!jetConst) continue;
1025 
1026  // Constituent eta/phi (relative to jet)
1027  Double_t deltaEta = jet->Eta() - jetConst->Eta();
1028  Double_t deltaPhi = TMath::Min(TMath::Abs(jet->Phi() - jetConst->Phi()), TMath::TwoPi() - TMath::Abs(jet->Phi() - jetConst->Phi()));
1029  if(!((jet->Phi() - jetConst->Phi() < 0) && (jet->Phi() - jetConst->Phi() <= TMath::Pi())))
1030  deltaPhi = -deltaPhi;
1031 
1032  FillHistogram("hConstituentPt", jetConst->Pt(), fCent);
1033  FillHistogram("hConstituentPhiEta", deltaPhi, deltaEta);
1034  }
1035 
1036  // ### Random cone / delta pT plots
1037  const Int_t kNumRandomConesPerEvent = 4;
1038  for(Int_t iCone=0; iCone<kNumRandomConesPerEvent; iCone++)
1039  {
1040  // Throw random cone
1041  Double_t tmpRandConeEta = fJetsCont->GetJetEtaMin() + gRandom->Rndm()*TMath::Abs(fJetsCont->GetJetEtaMax()-fJetsCont->GetJetEtaMin());
1042  Double_t tmpRandConePhi = gRandom->Rndm()*TMath::TwoPi();
1043  Double_t tmpRandConePt = 0;
1044  // Fill pT that is in cone
1045  fTracksCont->ResetCurrentID();
1046  while(AliVTrack *track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle()))
1047  if(IsTrackInCone(track, tmpRandConeEta, tmpRandConePhi, fJetsCont->GetJetRadius()))
1048  tmpRandConePt += track->Pt();
1049 
1050  // Fill histograms
1051  FillHistogram("hDeltaPt", tmpRandConePt - fJetsCont->GetRhoVal()*fJetsCont->GetJetRadius()*fJetsCont->GetJetRadius()*TMath::Pi(), fCent);
1052  }
1053 }
1054 
1055 //________________________________________________________________________
1057 {
1058  FillHistogram("hTrackPt", track->Pt(), fCent);
1059  FillHistogram("hTrackPhi", track->Phi(), fCent);
1060  FillHistogram("hTrackEta", track->Eta(), fCent);
1061  FillHistogram("hTrackEtaPt", track->Eta(), track->Pt());
1062  FillHistogram("hTrackPhiPt", track->Phi(), track->Pt());
1063  FillHistogram("hTrackPhiEta", track->Phi(), track->Eta());
1064 }
1065 
1066 //________________________________________________________________________
1067 void AliAnalysisTaskJetExtractor::AddPIDInformation(AliVParticle* particle, Float_t& sigITS, Float_t& sigTPC, Float_t& sigTOF, Float_t& sigTRD, Short_t& recoPID, Int_t& truePID)
1068 {
1069  truePID = 9;
1070  if(!particle) return;
1071 
1072  // If we have AODs, retrieve particle PID signals
1073  AliAODTrack* aodtrack = dynamic_cast<AliAODTrack*>(particle);
1074 
1075  if(aodtrack)
1076  {
1077  // Get AOD value from reco
1078  recoPID = aodtrack->GetMostProbablePID();
1079  AliAODPid* pidObj = aodtrack->GetDetPid();
1080  if(!pidObj)
1081  return;
1082 
1083  sigITS = pidObj->GetITSsignal();
1084  sigTPC = pidObj->GetTPCsignal();
1085  sigTOF = pidObj->GetTOFsignal();
1086  sigTRD = pidObj->GetTRDsignal();
1087  }
1088 
1089  // Get truth values if we are on MC
1091  {
1092  for(Int_t i=0; i<fTruthParticleArray->GetEntries();i++)
1093  {
1094  AliAODMCParticle* mcParticle = dynamic_cast<AliAODMCParticle*>(fTruthParticleArray->At(i));
1095  if(!mcParticle) continue;
1096 
1097  if (mcParticle->GetLabel() == particle->GetLabel())
1098  {
1099  if(fSaveTrackPDGCode)
1100  {
1101  truePID = mcParticle->PdgCode();
1102  }
1103  else // Use same convention as for PID in AODs
1104  {
1105  if(TMath::Abs(mcParticle->PdgCode()) == 2212) // proton
1106  truePID = 4;
1107  else if (TMath::Abs(mcParticle->PdgCode()) == 211) // pion
1108  truePID = 2;
1109  else if (TMath::Abs(mcParticle->PdgCode()) == 321) // kaon
1110  truePID = 3;
1111  else if (TMath::Abs(mcParticle->PdgCode()) == 11) // electron
1112  truePID = 0;
1113  else if (TMath::Abs(mcParticle->PdgCode()) == 13) // muon
1114  truePID = 1;
1115  else if (TMath::Abs(mcParticle->PdgCode()) == 700201) // deuteron
1116  truePID = 5;
1117  else if (TMath::Abs(mcParticle->PdgCode()) == 700301) // triton
1118  truePID = 6;
1119  else if (TMath::Abs(mcParticle->PdgCode()) == 700302) // He3
1120  truePID = 7;
1121  else if (TMath::Abs(mcParticle->PdgCode()) == 700202) // alpha
1122  truePID = 8;
1123  else
1124  truePID = 9;
1125  }
1126 
1127  break;
1128  }
1129  }
1130  }
1131 }
1132 
1133 //________________________________________________________________________
1135 {
1136  // Print properties for extraction
1137  // to be implemented later
1138 }
1139 
1140 
1141 //########################################################################
1142 // HELPERS
1143 //########################################################################
1144 
1145 //________________________________________________________________________
1146 inline Bool_t AliAnalysisTaskJetExtractor::IsTrackInCone(AliVParticle* track, Double_t eta, Double_t phi, Double_t radius)
1147 {
1148  // This is to use a full cone in phi even at the edges of phi (2pi -> 0) (0 -> 2pi)
1149  Double_t trackPhi = 0.0;
1150  if (track->Phi() > (TMath::TwoPi() - (radius-phi)))
1151  trackPhi = track->Phi() - TMath::TwoPi();
1152  else if (track->Phi() < (phi+radius - TMath::TwoPi()))
1153  trackPhi = track->Phi() + TMath::TwoPi();
1154  else
1155  trackPhi = track->Phi();
1156 
1157  if ( TMath::Abs(trackPhi-phi)*TMath::Abs(trackPhi-phi) + TMath::Abs(track->Eta()-eta)*TMath::Abs(track->Eta()-eta) <= radius*radius)
1158  return kTRUE;
1159 
1160  return kFALSE;
1161 }
1162 
1163 //________________________________________________________________________
1165 {
1166  TH1* tmpHist = static_cast<TH1*>(fOutput->FindObject(key));
1167  if(!tmpHist)
1168  {
1169  AliError(Form("Cannot find histogram <%s> ",key)) ;
1170  return;
1171  }
1172 
1173  tmpHist->Fill(x);
1174 }
1175 
1176 //________________________________________________________________________
1178 {
1179  TH1* tmpHist = static_cast<TH1*>(fOutput->FindObject(key));
1180  if(!tmpHist)
1181  {
1182  AliError(Form("Cannot find histogram <%s> ",key));
1183  return;
1184  }
1185 
1186  if (tmpHist->IsA()->GetBaseClass("TH1"))
1187  static_cast<TH1*>(tmpHist)->Fill(x,y); // Fill x with y
1188  else if (tmpHist->IsA()->GetBaseClass("TH2"))
1189  static_cast<TH2*>(tmpHist)->Fill(x,y); // Fill x,y with 1
1190 }
1191 
1192 //________________________________________________________________________
1194 {
1195  TH2* tmpHist = static_cast<TH2*>(fOutput->FindObject(key));
1196  if(!tmpHist)
1197  {
1198  AliError(Form("Cannot find histogram <%s> ",key));
1199  return;
1200  }
1201 
1202  tmpHist->Fill(x,y,add);
1203 }
1204 
1205 //________________________________________________________________________
1207 {
1208  TH3* tmpHist = static_cast<TH3*>(fOutput->FindObject(key));
1209  if(!tmpHist)
1210  {
1211  AliError(Form("Cannot find histogram <%s> ",key));
1212  return;
1213  }
1214 
1215  if(add)
1216  tmpHist->Fill(x,y,z,add);
1217  else
1218  tmpHist->Fill(x,y,z);
1219 }
1220 
1221 
1222 //________________________________________________________________________
1223 template <class T> T* AliAnalysisTaskJetExtractor::AddHistogram1D(const char* name, const char* title, const char* options, Int_t xBins, Double_t xMin, Double_t xMax, const char* xTitle, const char* yTitle)
1224 {
1225  T* tmpHist = new T(name, title, xBins, xMin, xMax);
1226 
1227  tmpHist->GetXaxis()->SetTitle(xTitle);
1228  tmpHist->GetYaxis()->SetTitle(yTitle);
1229  tmpHist->SetOption(options);
1230  tmpHist->SetMarkerStyle(kFullCircle);
1231  tmpHist->Sumw2();
1232 
1233  fOutput->Add(tmpHist);
1234 
1235  return tmpHist;
1236 }
1237 
1238 //________________________________________________________________________
1239 template <class T> T* AliAnalysisTaskJetExtractor::AddHistogram2D(const char* name, const char* title, const char* options, Int_t xBins, Double_t xMin, Double_t xMax, Int_t yBins, Double_t yMin, Double_t yMax, const char* xTitle, const char* yTitle, const char* zTitle)
1240 {
1241  T* tmpHist = new T(name, title, xBins, xMin, xMax, yBins, yMin, yMax);
1242  tmpHist->GetXaxis()->SetTitle(xTitle);
1243  tmpHist->GetYaxis()->SetTitle(yTitle);
1244  tmpHist->GetZaxis()->SetTitle(zTitle);
1245  tmpHist->SetOption(options);
1246  tmpHist->SetMarkerStyle(kFullCircle);
1247  tmpHist->Sumw2();
1248 
1249  fOutput->Add(tmpHist);
1250 
1251  return tmpHist;
1252 }
1253 
1254 //________________________________________________________________________
1255 template <class T> T* AliAnalysisTaskJetExtractor::AddHistogram3D(const char* name, const char* title, const char* options, Int_t xBins, Double_t xMin, Double_t xMax, Int_t yBins, Double_t yMin, Double_t yMax, Int_t zBins, Double_t zMin, Double_t zMax, const char* xTitle, const char* yTitle, const char* zTitle)
1256 {
1257  T* tmpHist = new T(name, title, xBins, xMin, xMax, yBins, yMin, yMax, zBins, zMin, zMax);
1258  tmpHist->GetXaxis()->SetTitle(xTitle);
1259  tmpHist->GetYaxis()->SetTitle(yTitle);
1260  tmpHist->GetZaxis()->SetTitle(zTitle);
1261  tmpHist->SetOption(options);
1262  tmpHist->SetMarkerStyle(kFullCircle);
1263  tmpHist->Sumw2();
1264 
1265  fOutput->Add(tmpHist);
1266 
1267  return tmpHist;
1268 }
1269 
1270 //________________________________________________________________________
1272 {
1273  // Called once at the end of the analysis.
1274 }
1275 
1276 // ### ADDTASK MACRO
1277 //________________________________________________________________________
1278 AliAnalysisTaskJetExtractor* AliAnalysisTaskJetExtractor::AddTaskJetExtractor(const char* configFile, AliRDHFJetsCutsVertex* vertexerCuts, const char* taskNameSuffix)
1279 {
1280  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
1281  TString trackArray = "tracks";
1282  TString jetArray = "Jet_AKTChargedR040_tracks_pT0150_E_scheme";
1283  TString rhoObject = "";
1284  Double_t jetRadius = 0.4;
1285  Double_t minJetEta = 0.5;
1286  Double_t minJetPt = 0.15;
1287  Double_t minTrackPt = 0.15;
1288  Double_t minJetAreaPerc = 0.557;
1289  TString suffix = "";
1290  if(taskNameSuffix != 0)
1291  suffix = taskNameSuffix;
1292 
1293  // ###### Load configuration from YAML files
1295  fYAMLConfig.AddConfiguration("alien:///alice/cern.ch/user/r/rhaake/ConfigFiles/JetExtractor_BaseConfig.yaml", "baseConfig");
1296  fYAMLConfig.AddConfiguration(configFile, "customConfig");
1297  fYAMLConfig.Initialize();
1298 
1299  fYAMLConfig.GetProperty("TrackArray", trackArray);
1300  fYAMLConfig.GetProperty("JetArray", jetArray);
1301  fYAMLConfig.GetProperty("RhoName", rhoObject, kFALSE);
1302  fYAMLConfig.GetProperty("JetRadius", jetRadius);
1303  fYAMLConfig.GetProperty("MinJetEta", minJetEta);
1304  fYAMLConfig.GetProperty("MinJetPt", minJetPt);
1305  fYAMLConfig.GetProperty("MinTrackPt", minTrackPt);
1306  fYAMLConfig.GetProperty("MinJetAreaPerc", minJetAreaPerc);
1307  fYAMLConfig.Print();
1308 
1309  // ###### Task name
1310  TString name("AliAnalysisTaskJetExtractor");
1311  if (jetArray != "") {
1312  name += "_";
1313  name += jetArray;
1314  }
1315  if (rhoObject != "") {
1316  name += "_";
1317  name += rhoObject;
1318  }
1319  if (suffix != "") {
1320  name += "_";
1321  name += suffix;
1322  }
1323 
1324  // ###### Setup task with default settings
1326  myTask->SetNeedEmcalGeom(kFALSE);
1327  myTask->SetVzRange(-10.,10.);
1328 
1329  // Particle container and track pt cut
1330  AliParticleContainer* trackCont = 0;
1331  if(trackArray == "mcparticles")
1332  trackCont = myTask->AddMCParticleContainer(trackArray);
1333  else if(trackArray =="mctracks")
1334  trackCont = myTask->AddParticleContainer(trackArray);
1335  else
1336  trackCont = myTask->AddTrackContainer(trackArray);
1337  trackCont->SetParticlePtCut(minTrackPt);
1338 
1339  // Secondary vertex cuts (default settings from PWGHF)
1340  // (can be overwritten by using myTask->SetVertexerCuts(cuts) from outside macro)
1341  AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts("AliESDtrackCuts", "default");
1342  esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
1343  esdTrackCuts->SetMinNClustersTPC(90);
1344  esdTrackCuts->SetMaxChi2PerClusterTPC(4);
1345  esdTrackCuts->SetRequireTPCRefit(kTRUE);
1346  esdTrackCuts->SetRequireITSRefit(kTRUE);
1347  esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD, AliESDtrackCuts::kAny);
1348  esdTrackCuts->SetMinDCAToVertexXY(0.);
1349  esdTrackCuts->SetEtaRange(-0.8, 0.8);
1350  esdTrackCuts->SetPtRange(1.0, 1.e10);
1351  vertexerCuts->AddTrackCuts(esdTrackCuts);
1352  vertexerCuts->SetNprongs(3);
1353  vertexerCuts->SetMinPtHardestTrack(1.0);//default 0.3
1354  vertexerCuts->SetSecVtxWithKF(kFALSE);//default with StrLinMinDist
1355  vertexerCuts->SetImpParCut(0.);//default 0
1356  vertexerCuts->SetDistPrimSec(0.);//default 0
1357  vertexerCuts->SetCospCut(-1);//default -1
1358  myTask->SetVertexerCuts(vertexerCuts);
1359 
1360  // Jet container
1361  AliJetContainer *jetCont = myTask->AddJetContainer(jetArray,6,jetRadius);
1362  if (jetCont) {
1363  jetCont->SetRhoName(rhoObject);
1364  jetCont->SetPercAreaCut(minJetAreaPerc);
1365  jetCont->SetJetPtCut(minJetPt);
1366  jetCont->SetLeadingHadronType(0);
1367  jetCont->SetPtBiasJetTrack(minTrackPt);
1368  jetCont->SetJetEtaLimits(-minJetEta, +minJetEta);
1369  jetCont->ConnectParticleContainer(trackCont);
1370  jetCont->SetMaxTrackPt(1000);
1371  }
1372 
1373  mgr->AddTask(myTask);
1374 
1375  // ###### Connect inputs/outputs
1376  mgr->ConnectInput (myTask, 0, mgr->GetCommonInputContainer() );
1377  mgr->ConnectOutput (myTask, 1, mgr->CreateContainer(Form("%s_histos", name.Data()), AliEmcalList::Class(), AliAnalysisManager::kOutputContainer, Form("%s:ChargedJetsHadronCF", mgr->GetCommonFileName())) );
1378  mgr->ConnectOutput (myTask, 2, mgr->CreateContainer(Form("%s_tree", name.Data()), TTree::Class(), AliAnalysisManager::kOutputContainer, mgr->GetCommonFileName()) );
1379 
1380  return myTask;
1381 }
Bool_t fSaveConstituentsIP
save arrays of constituent impact parameters
void ExecOnce()
Perform steps needed to initialize the analysis.
Int_t pdg
Bool_t fSaveSecondaryVertices
save reconstructed sec. vertex properties
static AliAnalysisTaskJetExtractor * AddTaskJetExtractor(const char *configFile, AliRDHFJetsCutsVertex *vertexerCuts, const char *taskNameSuffix=0)
void AddFlavourTag(Int_t tag)
Definition: AliEmcalJet.h:351
Double_t Area() const
Definition: AliEmcalJet.h:130
void SetParticlePtCut(Double_t cut)
virtual AliVParticle * GetNextAcceptParticle()
void SetSaveConstituents(Bool_t val)
Double_t GetRhoVal() const
double Double_t
Definition: External.C:58
TString fTruthParticleArrayName
Array name of MC particles in event (mcparticles)
Double_t GetJetEtaMin() const
const char * title
Definition: MakeQAPdf.C:27
void GetJetType(AliEmcalJet *jet, Int_t &typeHM, Int_t &typePM, Int_t &typeIC)
Float_t GetPartonEta6() const
void GetTrackImpactParameters(const AliVVertex *vtx, AliAODTrack *track, Float_t &d0, Float_t &d0cov, Float_t &z0, Float_t &z0cov)
Bool_t fSetEmcalJetFlavour
if set, the flavour property of the AliEmcalJets will be set
AliJetContainer * GetJetContainer(Int_t i=0) const
Definition: External.C:244
Float_t GetPartonEta7() const
TTree * fJetTree
! tree structure
Double_t Eta() const
Definition: AliEmcalJet.h:121
long long Long64_t
Definition: External.C:43
T * AddHistogram3D(const char *name="CustomHistogram", const char *title="NO_TITLE", const char *options="", Int_t xBins=100, Double_t xMin=0.0, Double_t xMax=20.0, Int_t yBins=100, Double_t yMin=0.0, Double_t yMax=20.0, Int_t zBins=100, Double_t zMin=0.0, Double_t zMax=20.0, const char *xTitle="x axis", const char *yTitle="y axis", const char *zTitle="z axis")
Float_t fBuffer_Event_BackgroundDensity
! array buffer
void SetLeadingHadronType(Int_t t)
Double_t Phi() const
Definition: AliEmcalJet.h:117
std::vector< Float_t > GetExtractionPercentages()
Bool_t fSaveMCInformation
save MC information
Double_t mass
void SetPtBiasJetTrack(Float_t b)
centrality
Int_t fTruthMinLabel
min track label to consider it as true particle
Double_t fSecondaryVertexMaxChi2
Max chi2 of secondary vertex (others will be discarded)
Double_t GetJetEtaMax() const
Float_t GetPartonPhi7() const
Double_t fSecondaryVertexMaxDispersion
Max dispersion of secondary vertex (others will be discarded)
AliJetContainer * AddJetContainer(const char *n, TString defaultCutType, Float_t jetRadius=0.4)
Double_t fEventCut_TriggerTrackMinPt
Event requirement, trigger track min pT.
std::vector< Float_t > GetExtractionPercentagePtBins()
Bool_t fSaveEventProperties
save general event properties (bgrd. etc.)
std::vector< Float_t > fTriggerTracks_Pt
If trigger track requi. is used -> save pT.
bool GetProperty(std::vector< std::string > propertyPath, const std::string &propertyName, T &property, const bool requiredProperty) const
Bool_t fSaveJetShapes
save jet shapes
void SetPercAreaCut(Float_t p)
Bool_t fSaveConstituentPID
save arrays of constituent PID parameters
Double_t fEventWeight
save event weight for each event (only works for JEWEL)
AliStack * stack
void SetVzRange(Double_t min, Double_t max)
const Int_t kMaxNumConstituents
Bool_t IsTrackInCone(AliVParticle *track, Double_t eta, Double_t phi, Double_t radius)
TRandom * gRandom
T * AddHistogram2D(const char *name="CustomHistogram", const char *title="NO_TITLE", const char *options="", Int_t xBins=100, Double_t xMin=0.0, Double_t xMax=20.0, Int_t yBins=100, Double_t yMin=0.0, Double_t yMax=20.0, const char *xTitle="x axis", const char *yTitle="y axis", const char *zTitle="z axis")
TString fTruthJetsRhoName
Array name for particle-level rho.
Int_t GetPartonFlag7() const
std::vector< Int_t > fExtractionJetTypes_HM
Jet-types that will be extracted with this tree (hadron matching)
Container for particles within the EMCAL framework.
Float_t GetPartonPhi6() const
Double_t fTrueJetMatchingRadius
Matching radius to true jet.
Int_t TrackAt(Int_t idx) const
Definition: AliEmcalJet.h:160
std::vector< Float_t > fTriggerTracks_Phi
If trigger track requi. is used -> save phi.
UShort_t GetNumberOfTracks() const
Definition: AliEmcalJet.h:139
UShort_t T(UShort_t m, UShort_t t)
Definition: RingBits.C:60
Int_t fEventCut_TriggerTrackMaxLabel
Event requirement, trigger track max label (can be used to selected embedded particles) ...
Float_t * fBuffer_Const_Phi
! array buffer
Float_t fBuffer_Event_ImpactParameter
! array buffer
Int_t GetPartonFlag6() const
void SetRhoName(const char *n)
Float_t fBuffer_Jet_MC_MatchedPt
! array buffer
Bool_t fSaveTriggerTracks
save event trigger track
AliParticleContainer * GetParticleContainer() const
AliParticleContainer * fTracksCont
! Tracks
std::vector< Float_t > fExtractionPercentages
Percentages which will be extracted for a given pT bin.
TClonesArray * fTruthParticleArray
! Array of MC particles in event (mcparticles)
Double_t fEventPercentage
percentage (0, 1] which will be extracted
void FillHistogram3D(const char *key, Double_t x, Double_t y, Double_t z, Double_t add=0)
int Int_t
Definition: External.C:63
Int_t fBuffer_Jet_MC_MotherParton
! array buffer
Float_t fBuffer_Event_PtHard
! array buffer
void SetJetPtCut(Float_t cut)
void SetSaveJetShapes(Bool_t val)
unsigned int UInt_t
Definition: External.C:33
Float_t * fBuffer_Const_Charge
! array buffer
float Float_t
Definition: External.C:68
Float_t fBuffer_Event_Weight
! array buffer
T * AddHistogram1D(const char *name="CustomHistogram", const char *title="NO_TITLE", const char *options="", Int_t xBins=100, Double_t xMin=0.0, Double_t xMax=20.0, const char *xTitle="x axis", const char *yTitle="y axis")
Int_t * fBuffer_Const_Label
! array buffer
AliParticleContainer * AddParticleContainer(const char *n)
Create new particle container and attach it to the task.
Double_t fImpactParameter
save IP for each event (only works for JEWEL)
Float_t fBuffer_JetPhi
! array buffer
Int_t fTruthMaxLabel
max track label to consider it as true particle
Float_t fBuffer_JetArea
! array buffer
std::vector< Int_t > fExtractionJetTypes_PM
Jet-types that will be extracted with this tree (parton matching)
std::vector< Float_t > fExtractionPercentagePtBins
pT-bins associated with fExtractionPercentages
Double_t fCent
!event centrality
Long64_t fBuffer_Event_ID
! array buffer
Float_t fBuffer_Event_Centrality
! array buffer
Int_t fBuffer_Jet_MC_MotherIC
! array buffer
std::ostream & Print(std::ostream &in, const int index=-1) const
Float_t * fBuffer_Const_ProdVtx_X
! array buffer
Bool_t fInitialized
init state of tree
Bool_t Run()
Run function. This is the core function of the analysis and contains the user code. Therefore users have to implement this function.
Analysis task that implements AliEmcalJetTree to extract jets to a tree.
void AddPIDInformation(AliVParticle *particle, Float_t &sigITS, Float_t &sigTPC, Float_t &sigTOF, Float_t &sigTRD, Short_t &recoPID, Int_t &truePID)
AliMCParticleContainer * AddMCParticleContainer(const char *n)
Create new container for MC particles and attach it to the task.
AliEmcalJet * GetNextAcceptJet()
AliHFJetsTaggingVertex * fVtxTagger
! class for sec. vertexing
short Short_t
Definition: External.C:23
TString fTruthJetsArrayName
Array name for particle-level jets.
void ConnectParticleContainer(AliParticleContainer *c)
Double_t Pt() const
Definition: AliEmcalJet.h:109
AliRDHFJetsCutsVertex * fVertexerCuts
Cuts used for the vertexer (given in add task macro)
Int_t fEventCut_TriggerTrackMinLabel
Event requirement, trigger track min label (can be used to selected embedded particles) ...
void SetSaveEventProperties(Bool_t val)
Float_t fPtHard
!event -hard
void FillHistogram(const char *key, Double_t x)
void FillTrackControlHistograms(AliVTrack *track)
Float_t GetJetRadius() const
AliEmcalList * fOutput
!output list
Float_t fBuffer_JetPt
! array buffer
Definition: External.C:220
void ReconstructSecondaryVertices(const AliVVertex *primVtx, const AliEmcalJet *jet, std::vector< Float_t > &secVtx_X, std::vector< Float_t > &secVtx_Y, std::vector< Float_t > &secVtx_Z, std::vector< Float_t > &secVtx_Mass, std::vector< Float_t > &secVtx_Lxy, std::vector< Float_t > &secVtx_SigmaLxy, std::vector< Float_t > &secVtx_Chi2, std::vector< Float_t > &secVtx_Dispersion)
Float_t fBuffer_Event_Vertex_Y
! array buffer
AliTrackContainer * AddTrackContainer(const char *n)
Create new track container and attach it to the task.
Bool_t fCalculateSecondaryVertices
Calculate the secondary vertices (instead of loading)
Float_t * fBuffer_Const_ProdVtx_Y
! array buffer
void GetJetTruePt(AliEmcalJet *jet, Double_t &matchedJetPt, Double_t &truePtFraction)
Bool_t fSaveConstituents
save arrays of constituent basic properties
YAML configuration class for AliPhysics.
void SetVertexerCuts(AliRDHFJetsCutsVertex *val)
Store some informaion about a Pythia eventThis class is used to store some information about a Pythia...
void SetMakeGeneralHistograms(Bool_t g)
void SetNeedEmcalGeom(Bool_t n)
Double_t fHadronMatchingRadius
Matching radius to search for beauty/charm hadrons around jet.
Base task in the EMCAL jet framework.
Float_t * fBuffer_Const_Pt
! array buffer
Represent a jet reconstructed using the EMCal jet framework.
Definition: AliEmcalJet.h:51
const AliEmcalPythiaInfo * GetPythiaInfo() const
Bool_t PxPyPz(Double_t p[3]) const
Definition: AliEmcalJet.h:111
Float_t * fBuffer_Const_Eta
! array buffer
Float_t fBuffer_JetEta
! array buffer
int AddConfiguration(std::string configurationFilename, std::string configurationName="")
Bool_t fSaveTrackPDGCode
save PDG code instead of code defined for AOD pid
Declaration of class AliEmcalPythiaInfo.
const char Option_t
Definition: External.C:48
Int_t fBuffer_NumConstituents
! array buffer
void UserCreateOutputObjects()
Main initialization function on the worker.
Float_t fBuffer_Event_Vertex_X
! array buffer
bool Bool_t
Definition: External.C:53
Double_t yMin
Int_t GetNAcceptedParticles() const
void SetMaxTrackPt(Float_t b)
void SetJetEtaLimits(Float_t min, Float_t max)
Int_t fBuffer_Jet_MC_MotherHadron
! array buffer
Double_t fEventCut_TriggerTrackMaxPt
Event requirement, trigger track max pT.
std::vector< Float_t > fTriggerTracks_Eta
If trigger track requi. is used -> save eta.
Float_t fBuffer_Shape_Mass
! array buffer
Double_t M() const
Definition: AliEmcalJet.h:120
Float_t * fBuffer_Const_ProdVtx_Z
! array buffer
Container for jet within the EMCAL jet framework.
Float_t fBuffer_Jet_MC_TruePtFraction
! array buffer
Class managing creation of a tree containing jets.
Bool_t AddJetToTree(AliEmcalJet *jet, Float_t bgrdDensity=0, Float_t vertexX=0, Float_t vertexY=0, Float_t vertexZ=0, Float_t centrality=0, Long64_t eventID=0, Float_t magField=0, AliParticleContainer *trackCont=0, Int_t motherParton=0, Int_t motherHadron=0, Int_t partonInitialCollision=0, Float_t matchedPt=0, Float_t truePtFraction=0, Float_t ptHard=0, Float_t eventWeight=0, Float_t impactParameter=0, std::vector< Float_t > &trackPID_ITS=DEFAULT_VECTOR_FLOAT, std::vector< Float_t > &trackPID_TPC=DEFAULT_VECTOR_FLOAT, std::vector< Float_t > &trackPID_TOF=DEFAULT_VECTOR_FLOAT, std::vector< Float_t > &trackPID_TRD=DEFAULT_VECTOR_FLOAT, std::vector< Short_t > &trackPID_Reco=DEFAULT_VECTOR_SHORT, std::vector< Int_t > &trackPID_Truth=DEFAULT_VECTOR_INT, std::vector< Float_t > &triggerTrackPt=DEFAULT_VECTOR_FLOAT, std::vector< Float_t > &triggerTrackDeltaEta=DEFAULT_VECTOR_FLOAT, std::vector< Float_t > &triggerTrackDeltaPhi=DEFAULT_VECTOR_FLOAT, std::vector< Float_t > &trackIP_d0=DEFAULT_VECTOR_FLOAT, std::vector< Float_t > &trackIP_z0=DEFAULT_VECTOR_FLOAT, std::vector< Float_t > &trackIP_d0cov=DEFAULT_VECTOR_FLOAT, std::vector< Float_t > &trackIP_z0cov=DEFAULT_VECTOR_FLOAT, std::vector< Float_t > &secVtx_X=DEFAULT_VECTOR_FLOAT, std::vector< Float_t > &secVtx_Y=DEFAULT_VECTOR_FLOAT, std::vector< Float_t > &secVtx_Z=DEFAULT_VECTOR_FLOAT, std::vector< Float_t > &secVtx_Mass=DEFAULT_VECTOR_FLOAT, std::vector< Float_t > &secVtx_Lxy=DEFAULT_VECTOR_FLOAT, std::vector< Float_t > &secVtx_SigmaLxy=DEFAULT_VECTOR_FLOAT, std::vector< Float_t > &secVtx_Chi2=DEFAULT_VECTOR_FLOAT, std::vector< Float_t > &secVtx_Dispersion=DEFAULT_VECTOR_FLOAT)
Definition: External.C:196
Double_t yMax
TList * OpenFile(const char *fname)
Definition: DrawAnaELoss.C:65
Float_t fBuffer_Event_Vertex_Z
! array buffer
void SetSaveTriggerTracks(Bool_t val)
Float_t fBuffer_Event_MagneticField
! array buffer