AliPhysics  f1bf8b7 (f1bf8b7)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliAnalysisTaskChargedJetsHadronCF.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-2016, 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 <algorithm>
17 #include <vector>
18 #include <TClonesArray.h>
19 #include <TF1.h>
20 #include <TH1F.h>
21 #include <TH2F.h>
22 #include <TH3F.h>
23 #include <THn.h>
24 #include <TTree.h>
25 #include <TList.h>
26 #include <TLorentzVector.h>
27 
28 #include "AliEmcalPythiaInfo.h"
29 #include "AliMCEvent.h"
30 #include "AliPythia.h"
31 #include "AliStack.h"
32 
33 #include "AliVTrack.h"
34 #include "AliVHeader.h"
35 #include "AliEmcalJet.h"
36 #include "AliRhoParameter.h"
37 #include "AliLog.h"
38 #include "AliJetContainer.h"
39 #include "AliTrackContainer.h"
40 #include "AliAODTrack.h"
41 #include "AliAODMCParticle.h"
42 #include "AliAODPid.h"
43 #include "AliPicoTrack.h"
44 #include "AliVParticle.h"
45 #include "TRandom3.h"
47 #include "AliBasicParticle.h"
48 
50 
51 //________________________________________________________________________
53 {
54 // dummy destructor
55 }
56 
60 //________________________________________________________________________
63  fJetsCont(0),
64  fTracksCont(0),
65  fEventExtractionPercentage(0),
66  fEventExtractionMinJetPt(0),
67  fEventExtractionMaxJetPt(0),
68  fConstPtFilterBit(1024),
69  fNumberOfCentralityBins(10),
70  fJetsOutput(),
71  fTracksOutput(),
72  fJetParticleArrayName("JetsDPhiBasicParticles"),
73  fTrackParticleArrayName(""),
74  fJetEmbeddingArray(),
75  fJetEmbeddingArrayName(""),
76  fJetEmbeddingTrackArrayName(""),
77  fJetEmbeddingMaxDistance(0.3),
78  fJetEmbeddingNumMatchedJets(2),
79  fJetEmbeddingUsePerTrackMCPercentage(kTRUE),
80  fJetEmbeddingUseBgrdForMCPercentage(kFALSE),
81  fJetEmbeddingCreatePtPlotPerCut(kFALSE),
82  fJetEmbeddingCuts(),
83  fJetVetoArray(),
84  fJetVetoArrayName(""),
85  fJetVetoJetByJet(1),
86  fMatchedJets(),
87  fRandom(0),
88  fTracksTree(0),
89  fTreeBufferTrack(0),
90  fTreeBufferPID(0),
91  fTreeBufferPDG(0),
92  fTrackExtractionPercentagePower(0),
93  fNumRandomConesPerEvent(10),
94  fJetOutputMode(0),
95  fLeadingJet(0),
96  fSubleadingJet(0),
97  fInitialPartonMatchedJet1(0),
98  fInitialPartonMatchedJet2(0),
99  fAcceptedJets(0),
100  fAcceptedTracks(0)
101 {
102  // Default constructor.
103  SetMakeGeneralHistograms(kTRUE);
104  fRandom = new TRandom3(0);
105 }
106 
107 
108 //________________________________________________________________________
110  AliAnalysisTaskEmcalJet(name, kTRUE),
111  fJetsCont(0),
112  fTracksCont(0),
113  fEventExtractionPercentage(0),
114  fEventExtractionMinJetPt(0),
115  fEventExtractionMaxJetPt(0),
116  fConstPtFilterBit(1024),
117  fNumberOfCentralityBins(10),
118  fJetsOutput(),
119  fTracksOutput(),
120  fJetParticleArrayName("JetsDPhiBasicParticles"),
121  fTrackParticleArrayName(""),
122  fJetEmbeddingArray(),
123  fJetEmbeddingArrayName(""),
124  fJetEmbeddingTrackArrayName(""),
125  fJetEmbeddingMaxDistance(0.3),
126  fJetEmbeddingNumMatchedJets(2),
127  fJetEmbeddingUsePerTrackMCPercentage(kTRUE),
128  fJetEmbeddingUseBgrdForMCPercentage(kFALSE),
129  fJetEmbeddingCreatePtPlotPerCut(kFALSE),
130  fJetEmbeddingCuts(),
131  fJetVetoArray(),
132  fJetVetoArrayName(""),
133  fJetVetoJetByJet(1),
134  fMatchedJets(),
135  fRandom(0),
136  fTracksTree(0),
137  fTreeBufferTrack(0),
138  fTreeBufferPID(0),
139  fTreeBufferPDG(0),
140  fTrackExtractionPercentagePower(0),
141  fNumRandomConesPerEvent(10),
142  fJetOutputMode(0),
143  fLeadingJet(0),
144  fSubleadingJet(0),
145  fInitialPartonMatchedJet1(0),
146  fInitialPartonMatchedJet2(0),
147  fAcceptedJets(0),
148  fAcceptedTracks(0)
149 {
150  // Constructor
152  fRandom = new TRandom3(0);
153 }
154 
155 //________________________________________________________________________
157 {
158  // Destructor.
159 }
160 
161 //________________________________________________________________________
163 {
165 
166  // ### Basic container settings
168  if(fJetsCont) { //get particles connected to jets
169  fJetsCont->PrintCuts();
171  } else { //no jets, just analysis tracks
173  }
174  if(fTracksCont) fTracksCont->SetClassName("AliVTrack");
175 
176  // ### Create all histograms
177  fHistEventRejection->GetXaxis()->SetBinLabel(15,"JetVeto");
178 
179  // Track QA plots
180  AddHistogram2D<TH2D>("hTrackPt", "Tracks p_{T} distribution", "", 300, 0., 300., fNumberOfCentralityBins, 0, 100, "p_{T} (GeV/c)", "Centrality", "dN^{Tracks}/dp_{T}");
181  AddHistogram2D<TH2D>("hTrackPhi", "Track angular distribution in #phi", "LEGO2", 180, 0., 2*TMath::Pi(), fNumberOfCentralityBins, 0, 100, "#phi", "Centrality", "dN^{Tracks}/(d#phi)");
182  AddHistogram2D<TH2D>("hTrackEta", "Track angular distribution in #eta", "LEGO2", 100, -2.5, 2.5, fNumberOfCentralityBins, 0, 100, "#eta", "Centrality", "dN^{Tracks}/(d#eta)");
183  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");
184 
185  AddHistogram2D<TH2D>("hLeadingTrackPt", "Leading tracks p_{T} distribution", "", 300, 0., 300., fNumberOfCentralityBins, 0, 100, "p_{T} (GeV/c)", "Centrality", "dN^{Tracks}/dp_{T}");
186  AddHistogram2D<TH2D>("hLeadingTrackPhi", "Leading tracks angular distribution in #phi", "LEGO2", 180, 0., 2*TMath::Pi(), fNumberOfCentralityBins, 0, 100, "#phi", "Centrality", "dN^{Tracks}/(d#phi)");
187  AddHistogram2D<TH2D>("hLeadingTrackEta", "Leading tracks angular distribution in #eta", "LEGO2", 100, -2.5, 2.5, fNumberOfCentralityBins, 0, 100, "#eta", "Centrality", "dN^{Tracks}/(d#eta)");
188  AddHistogram2D<TH2D>("hLeadingTrackPhiEta", "Track angular distribution #phi/#eta", "COLZ", 180, 0., 2*TMath::Pi(), 100, -2.5, 2.5, "#phi", "#eta", "dN^{Tracks}/d#phi d#eta");
189 
190  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})");
191  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})");
192 
193  // Create plots for each embedding cut
194  for(Int_t i = -1; i<static_cast<Int_t>(fJetEmbeddingCuts.size()); i++)
195  {
196  const char* appendix = "";
197  if(i>-1)
198  {
200  appendix = Form("_%s", currentCut.fCutName.Data());
201 
202  // Don't double-add cuts
203  if( static_cast<TH1*>(fOutput->FindObject(Form("hJetPtRaw%s", appendix))) )
204  continue;
205  }
206  // Jet QA plots
207  AddHistogram2D<TH2D>(Form("hJetPtRaw%s", appendix), "Jets p_{T} distribution (no bgrd. corr.)", "", 300, 0., 300., fNumberOfCentralityBins, 0, 100, "p_{T, jet} (GeV/c)", "Centrality", "dN^{Jets}/dp_{T}");
208  AddHistogram2D<TH2D>(Form("hJetPt%s", appendix), "Jets p_{T} distribution (background subtracted)", "", 400, -100., 300., fNumberOfCentralityBins, 0, 100, "p_{T, jet} (GeV/c)", "Centrality", "dN^{Jets}/dp_{T}");
209  AddHistogram2D<TH2D>(Form("hJetPhi%s", appendix), "Jet angular distribution #phi", "LEGO2", 180, 0., 2*TMath::Pi(), fNumberOfCentralityBins, 0, 100, "#phi", "Centrality", "dN^{Jets}/d#phi");
210  AddHistogram2D<TH2D>(Form("hJetEta%s", appendix), "Jet angular distribution #eta", "LEGO2", 100, -2.5, 2.5, fNumberOfCentralityBins, 0, 100, "#eta","Centrality","dN^{Jets}/d#eta");
211  AddHistogram2D<TH2D>(Form("hJetPhiPt%s", appendix), "Jet angular distribution #phi vs. p_{T}", "LEGO2", 180, 0., 2*TMath::Pi(), 400, -100., 300., "#phi", "p_{T, jet} (GeV/c)", "dN^{Jets}/d#phi dp_{T}");
212  AddHistogram2D<TH2D>(Form("hJetEtaPt%s", appendix), "Jet angular distribution #eta vs. p_{T}", "LEGO2", 100, -2.5, 2.5, 400, -100., 300., "#eta","p_{T, jet} (GeV/c)","dN^{Jets}/d#eta dp_{T}");
213  AddHistogram2D<TH2D>(Form("hJetPhiEta%s", appendix), "Jet angular distribution #phi/#eta", "COLZ", 180, 0., 2*TMath::Pi(), 100, -2.5, 2.5, "#phi", "#eta", "dN^{Jets}/d#phi d#eta");
214  AddHistogram2D<TH2D>(Form("hJetArea%s", appendix), "Jet area", "LEGO2", 200, 0., 2., fNumberOfCentralityBins, 0, 100, "Jet A", "Centrality", "dN^{Jets}/dA");
215  AddHistogram2D<TH2D>(Form("hJetAreaPt%s", appendix), "Jet area vs. p_{T}", "LEGO2", 200, 0., 2., 400, -100., 300., "Jet A", "p_{T, jet} (GeV/c)", "dN^{Jets}/dA dp_{T}");
216  AddHistogram2D<TH2D>(Form("hJetPtLeadingHadron%s", appendix), "Jet leading hadron p_{T} distribution vs. jet p_{T}", "", 300, 0., 300., 300, 0., 300., "p_{T, jet} (GeV/c)", "p_{T,lead had} (GeV/c)", "dN^{Jets}/dp_{T}dp_{T,had}");
217 
218  // Leading/subleading ...
219  AddHistogram2D<TH2D>(Form("hLeadingJetPtRaw%s", appendix), "Jets p_{T} distribution (no bgrd. corr.)", "", 300, 0., 300., fNumberOfCentralityBins, 0, 100, "p_{T, jet} (GeV/c)", "Centrality", "dN^{Jets}/dp_{T}");
220  AddHistogram2D<TH2D>(Form("hLeadingJetPt%s", appendix), "Jets p_{T} distribution (background subtracted)", "", 400, -100., 300., fNumberOfCentralityBins, 0, 100, "p_{T, jet} (GeV/c)", "Centrality", "dN^{Jets}/dp_{T}");
221  AddHistogram2D<TH2D>(Form("hLeadingJetPhi%s", appendix), "Jet angular distribution #phi", "LEGO2", 180, 0., 2*TMath::Pi(), fNumberOfCentralityBins, 0, 100, "#phi", "Centrality", "dN^{Jets}/d#phi");
222  AddHistogram2D<TH2D>(Form("hLeadingJetEta%s", appendix), "Jet angular distribution #eta", "LEGO2", 100, -2.5, 2.5, fNumberOfCentralityBins, 0, 100, "#eta","Centrality","dN^{Jets}/d#eta");
223  AddHistogram2D<TH2D>(Form("hLeadingJetPhiPt%s", appendix), "Jet angular distribution #phi vs. p_{T}", "LEGO2", 180, 0., 2*TMath::Pi(), 400, -100., 300., "#phi", "p_{T, jet} (GeV/c)", "dN^{Jets}/d#phi dp_{T}");
224  AddHistogram2D<TH2D>(Form("hLeadingJetEtaPt%s", appendix), "Jet angular distribution #eta vs. p_{T}", "LEGO2", 100, -2.5, 2.5, 400, -100., 300., "#eta","p_{T, jet} (GeV/c)","dN^{Jets}/d#eta dp_{T}");
225  AddHistogram2D<TH2D>(Form("hLeadingJetPhiEta%s", appendix), "Jet angular distribution #phi/#eta", "COLZ", 180, 0., 2*TMath::Pi(), 100, -2.5, 2.5, "#phi", "#eta", "dN^{Jets}/d#phi d#eta");
226  AddHistogram2D<TH2D>(Form("hLeadingJetArea%s", appendix), "Jet area", "LEGO2", 200, 0., 2., fNumberOfCentralityBins, 0, 100, "Jet A", "Centrality", "dN^{Jets}/dA");
227  AddHistogram2D<TH2D>(Form("hLeadingJetAreaPt%s", appendix), "Jet area vs. p_{T}", "LEGO2", 200, 0., 2., 400, -100., 300., "Jet A", "p_{T, jet} (GeV/c)", "dN^{Jets}/dA dp_{T}");
228  AddHistogram2D<TH2D>(Form("hLeadingJetPtLeadingHadron%s", appendix), "Jet leading hadron p_{T} distribution vs. jet p_{T}", "", 300, 0., 300., 300, 0., 300., "p_{T, jet} (GeV/c)", "p_{T,lead had} (GeV/c)", "dN^{Jets}/dp_{T}dp_{T,had}");
229 
230  AddHistogram2D<TH2D>(Form("hSubleadingJetPtRaw%s", appendix), "Jets p_{T} distribution (no bgrd. corr.)", "", 300, 0., 300., fNumberOfCentralityBins, 0, 100, "p_{T, jet} (GeV/c)", "Centrality", "dN^{Jets}/dp_{T}");
231  AddHistogram2D<TH2D>(Form("hSubleadingJetPt%s", appendix), "Jets p_{T} distribution (background subtracted)", "", 400, -100., 300., fNumberOfCentralityBins, 0, 100, "p_{T, jet} (GeV/c)", "Centrality", "dN^{Jets}/dp_{T}");
232  AddHistogram2D<TH2D>(Form("hSubleadingJetPhi%s", appendix), "Jet angular distribution #phi", "LEGO2", 180, 0., 2*TMath::Pi(), fNumberOfCentralityBins, 0, 100, "#phi", "Centrality", "dN^{Jets}/d#phi");
233  AddHistogram2D<TH2D>(Form("hSubleadingJetEta%s", appendix), "Jet angular distribution #eta", "LEGO2", 100, -2.5, 2.5, fNumberOfCentralityBins, 0, 100, "#eta","Centrality","dN^{Jets}/d#eta");
234  AddHistogram2D<TH2D>(Form("hSubleadingJetPhiPt%s", appendix), "Jet angular distribution #phi vs. p_{T}", "LEGO2", 180, 0., 2*TMath::Pi(), 400, -100., 300., "#phi", "p_{T, jet} (GeV/c)", "dN^{Jets}/d#phi dp_{T}");
235  AddHistogram2D<TH2D>(Form("hSubleadingJetEtaPt%s", appendix), "Jet angular distribution #eta vs. p_{T}", "LEGO2", 100, -2.5, 2.5, 400, -100., 300., "#eta","p_{T, jet} (GeV/c)","dN^{Jets}/d#eta dp_{T}");
236  AddHistogram2D<TH2D>(Form("hSubleadingJetPhiEta%s", appendix), "Jet angular distribution #phi/#eta", "COLZ", 180, 0., 2*TMath::Pi(), 100, -2.5, 2.5, "#phi", "#eta", "dN^{Jets}/d#phi d#eta");
237  AddHistogram2D<TH2D>(Form("hSubleadingJetArea%s", appendix), "Jet area", "LEGO2", 200, 0., 2., fNumberOfCentralityBins, 0, 100, "Jet A", "Centrality", "dN^{Jets}/dA");
238  AddHistogram2D<TH2D>(Form("hSubleadingJetAreaPt%s", appendix), "Jet area vs. p_{T}", "LEGO2", 200, 0., 2., 400, -100., 300., "Jet A", "p_{T, jet} (GeV/c)", "dN^{Jets}/dA dp_{T}");
239  AddHistogram2D<TH2D>(Form("hSubleadingJetPtLeadingHadron%s", appendix), "Jet leading hadron p_{T} distribution vs. jet p_{T}", "", 300, 0., 300., 300, 0., 300., "p_{T, jet} (GeV/c)", "p_{T,lead had} (GeV/c)", "dN^{Jets}/dp_{T}dp_{T,had}");
240 
241  AddHistogram2D<TH2D>(Form("hJetConstituentPt_Cent0_100%s", appendix), "Jet constituent p_{T} distribution vs. jet p_T (background subtracted)", "", 400, -100., 300., 300, 0., 300., "p_{T, jet} (GeV/c)", "p_{T, track} (GeV/c)", "dN^{Tracks}/d^{2}p_{T}");
242  AddHistogram2D<TH2D>(Form("hJetConstituentPt_Cent0_10%s", appendix), "Jet constituent p_{T} distribution vs. jet p_T (background subtracted), 0-10 centrality", "", 400, -100., 300., 300, 0., 300., "p_{T, jet} (GeV/c)", "p_{T, track} (GeV/c)", "dN^{Tracks}/d^{2}p_{T}");
243  AddHistogram2D<TH2D>(Form("hJetConstituentPt_Cent10_30%s", appendix), "Jet constituent p_{T} distribution vs. jet p_T (background subtracted), 10-30 centrality", "", 400, -100., 300., 300, 0., 300., "p_{T, jet} (GeV/c)", "p_{T, track} (GeV/c)", "dN^{Tracks}/d^{2}p_{T}");
244  AddHistogram2D<TH2D>(Form("hJetConstituentPt_Cent30_50%s", appendix), "Jet constituent p_{T} distribution vs. jet p_T (background subtracted), 30-50 centrality", "", 400, -100., 300., 300, 0., 300., "p_{T, jet} (GeV/c)", "p_{T, track} (GeV/c)", "dN^{Tracks}/d^{2}p_{T}");
245  AddHistogram2D<TH2D>(Form("hJetConstituentPt_Cent50_90%s", appendix), "Jet constituent p_{T} distribution vs. jet p_T (background subtracted), 50-90 centrality", "", 400, -100., 300., 300, 0., 300., "p_{T, jet} (GeV/c)", "p_{T, track} (GeV/c)", "dN^{Tracks}/d^{2}p_{T}");
246 
247  AddHistogram2D<TH2D>(Form("hJetConstituentPt_Cent0_100_FilterBit%i%s", fConstPtFilterBit, appendix), "Jet constituent p_{T} distribution vs. jet p_T (background subtracted)", "", 400, -100., 300., 300, 0., 300., "p_{T, jet} (GeV/c)", "p_{T, track} (GeV/c)", "dN^{Tracks}/d^{2}p_{T}");
248  AddHistogram2D<TH2D>(Form("hJetConstituentPt_Cent0_10_FilterBit%i%s", fConstPtFilterBit, appendix), "Jet constituent p_{T} distribution vs. jet p_T (background subtracted), 0-10 centrality", "", 400, -100., 300., 300, 0., 300., "p_{T, jet} (GeV/c)", "p_{T, track} (GeV/c)", "dN^{Tracks}/d^{2}p_{T}");
249  AddHistogram2D<TH2D>(Form("hJetConstituentPt_Cent10_30_FilterBit%i%s", fConstPtFilterBit, appendix), "Jet constituent p_{T} distribution vs. jet p_T (background subtracted), 10-30 centrality", "", 400, -100., 300., 300, 0., 300., "p_{T, jet} (GeV/c)", "p_{T, track} (GeV/c)", "dN^{Tracks}/d^{2}p_{T}");
250  AddHistogram2D<TH2D>(Form("hJetConstituentPt_Cent30_50_FilterBit%i%s", fConstPtFilterBit, appendix), "Jet constituent p_{T} distribution vs. jet p_T (background subtracted), 30-50 centrality", "", 400, -100., 300., 300, 0., 300., "p_{T, jet} (GeV/c)", "p_{T, track} (GeV/c)", "dN^{Tracks}/d^{2}p_{T}");
251  AddHistogram2D<TH2D>(Form("hJetConstituentPt_Cent50_90_FilterBit%i%s", fConstPtFilterBit, appendix), "Jet constituent p_{T} distribution vs. jet p_T (background subtracted), 50-90 centrality", "", 400, -100., 300., 300, 0., 300., "p_{T, jet} (GeV/c)", "p_{T, track} (GeV/c)", "dN^{Tracks}/d^{2}p_{T}");
252 
253  AddHistogram2D<TH2D>(Form("hJetConstituentCount_Cent0_100%s", appendix), "Jet constituent count vs. jet p_T (background subtracted)", "", 400, -100., 300., 200, 0., 200., "p_{T, jet} (GeV/c)", "Count", "dN^{Jets}/dNdp_{T}");
254  AddHistogram2D<TH2D>(Form("hJetConstituentCount_Cent0_10%s", appendix), "Jet constituent count vs. jet p_T (background subtracted), 0-10 centrality", "", 400, -100., 300., 200, 0., 200., "p_{T, jet} (GeV/c)", "Count", "dN^{Jets}/dNdp_{T}");
255  }
256 
257  // Embedding plots
258  if(fJetOutputMode == 4 || fJetOutputMode == 5)
259  {
260  Double_t maxRatio = 1.;
262  maxRatio = 2.;
263 
264  for(Int_t i = -1; i<static_cast<Int_t>(fJetEmbeddingCuts.size()); i++)
265  {
266  const char* appendix = "";
267  if(i>-1)
268  {
270  appendix = Form("_%s", currentCut.fCutName.Data());
271 
272  // Don't double-add cuts
273  if( static_cast<TH1*>(fOutput->FindObject(Form("hEmbeddingDeltaR%s", appendix))) )
274  continue;
275  }
276  AddHistogram2D<TH2D>(Form("hEmbeddingDeltaR%s", appendix), "Matched jet #Delta R distribution", "", 200, -50., 150., 100, 0, 1.0, "p_{T, jet} (GeV/c)", "#Delta R", "dN^{Matched}/dp_{T}dR");
277  AddHistogram2D<TH2D>(Form("hEmbeddingDeltaEta%s", appendix), "Matched jet #Delta #eta distribution", "", 200, -50., 150., 100, -1.0, 1.0, "p_{T, jet} (GeV/c)", "#Delta #eta", "dN^{Matched}/dp_{T}d#eta");
278  AddHistogram2D<TH2D>(Form("hEmbeddingDeltaPhi%s", appendix), "Matched jet #Delta #phi distribution", "", 200, -50., 150., 100, -1.0, 1.0, "p_{T, jet} (GeV/c)", "#Delta #phi", "dN^{Matched}/dp_{T}d#phi");
279  AddHistogram1D<TH1D>(Form("hEmbeddingJetPt%s", appendix), "Embedded jets p_{T} distribution", "", 200, -50., 150., "p_{T, jet} (GeV/c)", "dN/dp_{T}");
280  AddHistogram2D<TH2D>(Form("hEmbeddingJetPhiEta%s", appendix), "Embedded jet angular distribution #phi/#eta", "COLZ", 180, 0., 2*TMath::Pi(), 100, -2.5, 2.5, "#phi", "#eta", "dN^{Jets}/d#phi d#eta");
281 
283  {
284  AddHistogram3D<TH3D>(Form("hEmbeddingPtCorr010%s", appendix), "Matched jet p_{T} distributions (0-10% centrality)", "", 150, 0., 150., 150, 0., 150., 100, 0., maxRatio, "p_{T, MC jet} (GeV/c)", "p_{T, emb} (GeV/c)", "% MC");
285  AddHistogram3D<TH3D>(Form("hEmbeddingPtCorr1030%s", appendix), "Matched jet p_{T} distributions (10-30% centrality)", "", 150, 0., 150., 150, 0., 150., 100, 0., maxRatio, "p_{T, MC jet} (GeV/c)", "p_{T, emb} (GeV/c)", "% MC");
286  AddHistogram3D<TH3D>(Form("hEmbeddingPtCorr3050%s", appendix), "Matched jet p_{T} distributions (30-50% centrality)", "", 150, 0., 150., 150, 0., 150., 100, 0., maxRatio, "p_{T, MC jet} (GeV/c)", "p_{T, emb} (GeV/c)", "% MC");
287  AddHistogram3D<TH3D>(Form("hEmbeddingPtCorr5090%s", appendix), "Matched jet p_{T} distributions (50-90% centrality)", "", 150, 0., 150., 150, 0., 150., 100, 0., maxRatio, "p_{T, MC jet} (GeV/c)", "p_{T, emb} (GeV/c)", "% MC");
288  }
289  else
290  {
291  AddHistogram2D<TH2D>(Form("hEmbeddingPtCorr010_Above20%s", appendix), "Matched jet p_{T} distributions, MC ratio > 20% (0-10% centrality)", "", 150, 0., 150., 150, 0., 150.,"p_{T, MC jet} (GeV/c)", "p_{T, emb} (GeV/c)", "% MC");
292  AddHistogram2D<TH2D>(Form("hEmbeddingPtCorr1030_Above20%s", appendix),"Matched jet p_{T} distributions, MC ratio > 20% (10-30% centrality)","", 150, 0., 150., 150, 0., 150.,"p_{T, MC jet} (GeV/c)", "p_{T, emb} (GeV/c)", "% MC");
293  AddHistogram2D<TH2D>(Form("hEmbeddingPtCorr3050_Above20%s", appendix),"Matched jet p_{T} distributions, MC ratio > 20% (30-50% centrality)","", 150, 0., 150., 150, 0., 150.,"p_{T, MC jet} (GeV/c)", "p_{T, emb} (GeV/c)", "% MC");
294  AddHistogram2D<TH2D>(Form("hEmbeddingPtCorr5090_Above20%s", appendix),"Matched jet p_{T} distributions, MC ratio > 20% (50-90% centrality)","", 150, 0., 150., 150, 0., 150.,"p_{T, MC jet} (GeV/c)", "p_{T, emb} (GeV/c)", "% MC");
295  }
296  }
297 
299  {
300  AddHistogram3D<TH3D>("hEmbeddingPtCorr010", "Matched jet p_{T} distributions (0-10% centrality)", "", 150, 0., 150., 150, 0., 150., 100, 0., maxRatio, "p_{T, MC jet} (GeV/c)", "p_{T, emb} (GeV/c)", "% MC");
301  AddHistogram3D<TH3D>("hEmbeddingPtCorr1030", "Matched jet p_{T} distributions (10-30% centrality)", "", 150, 0., 150., 150, 0., 150., 100, 0., maxRatio, "p_{T, MC jet} (GeV/c)", "p_{T, emb} (GeV/c)", "% MC");
302  AddHistogram3D<TH3D>("hEmbeddingPtCorr3050", "Matched jet p_{T} distributions (30-50% centrality)", "", 150, 0., 150., 150, 0., 150., 100, 0., maxRatio, "p_{T, MC jet} (GeV/c)", "p_{T, emb} (GeV/c)", "% MC");
303  AddHistogram3D<TH3D>("hEmbeddingPtCorr5090", "Matched jet p_{T} distributions (50-90% centrality)", "", 150, 0., 150., 150, 0., 150., 100, 0., maxRatio, "p_{T, MC jet} (GeV/c)", "p_{T, emb} (GeV/c)", "% MC");
304  }
305  }
306 
307  // Random cone plots, background, ...
308  AddHistogram2D<TH2D>("hRandomConePt", "Random cone p_{T} distribution", "", 400, -100., 300., fNumberOfCentralityBins, 0, 100, "p_{T, cone} (GeV/c)", "Centrality", "dN^{Tracks}/dp_{T}");
309  AddHistogram2D<TH2D>("hRandomConePtCut3GeV", "Random cone p_{T} distribution, cut p_{T} > 3 GeV/c", "", 400, -100., 300., fNumberOfCentralityBins, 0, 100, "p_{T, cone} (GeV/c)", "Centrality", "dN^{Tracks}/dp_{T}");
310  AddHistogram2D<TH2D>("hRandomConeRawPt", "Random cone p_{T} distribution (no bgrd. correction)", "", 300, 0., 300., fNumberOfCentralityBins, 0, 100, "p_{T, cone} (GeV/c)", "Centrality", "dN^{Tracks}/dp_{T}");
311  AddHistogram2D<TH2D>("hRandomConeRawPtCut3GeV", "Random cone p_{T} distribution (no bgrd. correction), cut p_{T} > 3 GeV/c", "", 300, 0., 300., fNumberOfCentralityBins, 0, 100, "p_{T, cone} (GeV/c)", "Centrality", "dN^{Tracks}/dp_{T}");
312 
313  AddHistogram2D<TH2D>("hTrackCount", "Number of tracks in acceptance vs. centrality", "LEGO2", 500, 0., 5000., fNumberOfCentralityBins, 0, 100, "N tracks","Centrality", "dN^{Events}/dN^{Tracks}");
314  AddHistogram2D<TH2D>("hJetCount", "Number of jets in acceptance vs. centrality", "LEGO2", 100, 0., 100., fNumberOfCentralityBins, 0, 100, "N Jets","Centrality", "dN^{Events}/dN^{Jets}");
315  AddHistogram2D<TH2D>("hBackgroundPt", "Background p_{T} distribution", "", 150, 0., 150., fNumberOfCentralityBins, 0, 100, "Background p_{T} (GeV/c)", "Centrality", "dN^{Events}/dp_{T}");
316 
317 
318  for(Int_t i = -1; i<static_cast<Int_t>(fJetEmbeddingCuts.size()); i++)
319  {
320  const char* appendix = "";
321  if(i>-1)
322  {
324  appendix = Form("_%s", currentCut.fCutName.Data());
325 
326  // Don't double-add cuts
327  if( static_cast<TH1*>(fOutput->FindObject(Form("hBackgroundPtJetPt_Cent0_100%s", appendix))) )
328  continue;
329  }
330  AddHistogram2D<TH2D>(Form("hBackgroundPtJetPt_Cent0_100%s", appendix), "Background p_{T} distribution vs. jet p_{T}", "", 150, 0., 150., 400, -100., 300., "Background p_{T} (GeV/c)", "Jet p_{T} (GeV/c)", "dN^{Events}/dp_{T}");
331  AddHistogram2D<TH2D>(Form("hBackgroundPtJetPt_Cent0_10%s", appendix), "Background p_{T} distribution vs. jet p_{T}", "", 150, 0., 150., 400, -100., 300., "Background p_{T} (GeV/c)", "Jet p_{T} (GeV/c)", "dN^{Events}/dp_{T}");
332  AddHistogram2D<TH2D>(Form("hBackgroundPtConstCount_Cent0_100%s", appendix), "Background p_{T} distribution vs. const. count", "", 150, 0., 150., 200, 0., 200., "Background p_{T} (GeV/c)", "Count", "dN^{Events}/dp_{T}");
333  AddHistogram2D<TH2D>(Form("hBackgroundPtConstCount_Cent0_10%s", appendix), "Background p_{T} distribution vs. const. count", "", 150, 0., 150., 200, 0., 200., "Background p_{T} (GeV/c)", "Count", "dN^{Events}/dp_{T}");
334  }
335 
336  PostData(1, fOutput); // Post data for ALL output slots > 0 here.
337 }
338 
339 
340 //________________________________________________________________________
342 
344 
345 
346  // ### Add the jets as basic correlation particles to the event
347  // This output object carries all accepted jets
348  if (!(fInputEvent->FindListObject(Form("%s", fJetParticleArrayName.Data()))))
349  {
350  fJetsOutput.push_back(new TClonesArray("AliPicoTrack"));
351  fJetsOutput.at(0)->SetName(fJetParticleArrayName.Data());
352  fInputEvent->AddObject(fJetsOutput.at(0));
353  }
354  else
355  AliFatal(Form("%s: Object with name %s already in event!", GetName(), Form("%s", fJetParticleArrayName.Data())));
356 
357  // These output objects carry all jets that pass certain cuts
358  if( (fJetOutputMode==4 || fJetOutputMode==5) )
359  {
360  // before, check if all given output names are OK
361  for(Int_t i = 0; i<fJetEmbeddingCuts.size(); i++)
362  if (fInputEvent->FindListObject(Form("%s", fJetEmbeddingCuts.at(i).fOutputName.Data())))
363  AliFatal(Form("%s: Object with name %s already in event!", GetName(), Form("%s", fJetEmbeddingCuts.at(i).fOutputName.Data())));
364 
365  for(Int_t i = 0; i<fJetEmbeddingCuts.size(); i++)
366  {
367  // If the cut demands a new output stream, add it
368  if (!fInputEvent->FindListObject(Form("%s", fJetEmbeddingCuts.at(i).fOutputName.Data())))
369  {
370  fJetsOutput.push_back(new TClonesArray("AliPicoTrack"));
371  fJetsOutput.at(fJetsOutput.size()-1)->SetName(fJetEmbeddingCuts.at(i).fOutputName.Data());
372  fInputEvent->AddObject(fJetsOutput.at(fJetsOutput.size()-1));
373  fJetEmbeddingCuts.at(i).fArrayIndex = fJetsOutput.size()-1;
374 
375  // Set the array indices for all cuts that use this output stream
376  for(Int_t j = 0; j<fJetEmbeddingCuts.size(); j++)
377  {
378  if(fJetEmbeddingCuts.at(j).fArrayIndex != -1)
379  continue;
380  if(fJetEmbeddingCuts.at(j).fOutputName != fJetEmbeddingCuts.at(i).fOutputName)
381  continue;
382  fJetEmbeddingCuts.at(j).fArrayIndex = fJetEmbeddingCuts.at(i).fArrayIndex;
383  }
384  }
385  }
386  }
387 
388  // ##############################################
389  // ##############################################
390 
391  // ### Prepare the track tree
393  {
394  fTracksTree = new TTree("ExtractedTracks", "ExtractedTracks");
395  fTracksTree->Branch("Kinematics", "AliBasicParticle", &fTreeBufferTrack, 1000);
396  fTracksTree->Branch("PID", "AliAODPid", &fTreeBufferPID, 1000);
397  fTracksTree->Branch("PDG",&fTreeBufferPDG,"a/I");
398  fOutput->Add(fTracksTree);
399  }
400 
401  // ### Add the tracks as basic correlation particles to the event (optional)
402  if(fTrackParticleArrayName != "")
403  {
404  if (!(fInputEvent->FindListObject(Form("%s", fTrackParticleArrayName.Data()))))
405  {
406  fTracksOutput = new TClonesArray("AliPicoTrack");
407  fTracksOutput->SetName(fTrackParticleArrayName.Data());
408  fInputEvent->AddObject(fTracksOutput);
409  }
410  else
411  AliFatal(Form("%s: Object with name %s already in event!", GetName(), Form("%s", fTrackParticleArrayName.Data())));
412  }
413 
414  // ### Import jets for embedding (optional)
415  if(fJetEmbeddingArrayName != "")
416  {
417  fJetEmbeddingArray = static_cast<TClonesArray*>(InputEvent()->FindListObject(Form("%s", fJetEmbeddingArrayName.Data())));
418  if(!fJetEmbeddingArray)
419  AliFatal(Form("Importing jets for embedding failed! Array '%s' not found!", fJetEmbeddingArrayName.Data()));
420  }
421  else if(fJetOutputMode==4 || fJetOutputMode==5)
422  AliFatal(Form("fJetEmbeddingArrayName must be set in jet output mode 4 or 5."));
423 
424  // ### Import veto jets for matching (optional)
425  if(fJetVetoArrayName != "")
426  {
427  fJetVetoArray = static_cast<TClonesArray*>(InputEvent()->FindListObject(Form("%s", fJetVetoArrayName.Data())));
428  if(!fJetVetoArray)
429  AliFatal(Form("Importing jets for veto failed! Array '%s' not found!", fJetVetoArrayName.Data()));
430  }
431 }
432 
433 //________________________________________________________________________
435 {
436  if(fJetOutputMode==3) // output leading&subleading jet
437  {
438  if((jet!=fLeadingJet) && (jet!=fSubleadingJet))
439  return kFALSE;
440  }
441  else if(fJetOutputMode==1) // output the leading jet
442  {
443  if(jet!=fLeadingJet)
444  return kFALSE;
445  }
446  else if(fJetOutputMode==2) // output the subleading jet
447  {
448  if(jet!=fSubleadingJet)
449  return kFALSE;
450  }
451  else if(fJetOutputMode==6)
452  {
454  return kFALSE;
455  }
456 
457  if(fJetOutputMode==4) // matching jets only
458  return (std::find(fMatchedJets.begin(), fMatchedJets.end(), jet) != fMatchedJets.end());
459  else if(fJetOutputMode==5) // non-matching jets only
460  return (std::find(fMatchedJets.begin(), fMatchedJets.end(), jet) == fMatchedJets.end());
461 
462  return kTRUE;
463 }
464 
465 //________________________________________________________________________
467 {
468  TString appendix("");
469  if(cutName)
470  appendix = Form("_%s", cutName);
471 
472  // All jets
473  FillHistogram(Form("hJetPtRaw%s", appendix.Data()), jet->Pt(), fCent);
474  FillHistogram(Form("hJetPt%s", appendix.Data()), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), fCent);
475  FillHistogram(Form("hJetPhi%s", appendix.Data()), jet->Phi(), fCent);
476  FillHistogram(Form("hJetEta%s", appendix.Data()), jet->Eta(), fCent);
477  FillHistogram(Form("hJetEtaPt%s", appendix.Data()), jet->Eta(), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area());
478  FillHistogram(Form("hJetPhiPt%s", appendix.Data()), jet->Phi(), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area());
479  FillHistogram(Form("hJetPhiEta%s", appendix.Data()), jet->Phi(), jet->Eta());
480  FillHistogram(Form("hJetArea%s", appendix.Data()), jet->Area(), fCent);
481  FillHistogram(Form("hJetAreaPt%s", appendix.Data()), jet->Area(), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area());
482  FillHistogram(Form("hJetPtLeadingHadron%s", appendix.Data()), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), fJetsCont->GetLeadingHadronPt(jet));
483 
484  FillHistogram(Form("hBackgroundPtJetPt_Cent0_100%s", appendix.Data()), fJetsCont->GetRhoVal(), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area());
485  if( (fCent >= 0) && (fCent < 10) )
486  FillHistogram(Form("hBackgroundPtJetPt_Cent0_10%s", appendix.Data()), fJetsCont->GetRhoVal(), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area());
487 
488  // Leading jet plots
489  if(jet==fLeadingJet)
490  {
491  FillHistogram(Form("hLeadingJetPtRaw%s", appendix.Data()), jet->Pt(), fCent);
492  FillHistogram(Form("hLeadingJetPt%s", appendix.Data()), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), fCent);
493  FillHistogram(Form("hLeadingJetPhi%s", appendix.Data()), jet->Phi(), fCent);
494  FillHistogram(Form("hLeadingJetEta%s", appendix.Data()), jet->Eta(), fCent);
495  FillHistogram(Form("hLeadingJetEtaPt%s", appendix.Data()), jet->Eta(), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area());
496  FillHistogram(Form("hLeadingJetPhiPt%s", appendix.Data()), jet->Phi(), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area());
497  FillHistogram(Form("hLeadingJetPhiEta%s", appendix.Data()), jet->Phi(), jet->Eta());
498  FillHistogram(Form("hLeadingJetArea%s", appendix.Data()), jet->Area(), fCent);
499  FillHistogram(Form("hLeadingJetAreaPt%s", appendix.Data()), jet->Area(), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area());
500  FillHistogram(Form("hLeadingJetPtLeadingHadron%s", appendix.Data()), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), fJetsCont->GetLeadingHadronPt(jet));
501  }
502 
503  // Subleading jet plot
504  else if(jet==fSubleadingJet)
505  {
506  FillHistogram(Form("hSubleadingJetPtRaw%s", appendix.Data()), jet->Pt(), fCent);
507  FillHistogram(Form("hSubleadingJetPt%s", appendix.Data()), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), fCent);
508  FillHistogram(Form("hSubleadingJetPhi%s", appendix.Data()), jet->Phi(), fCent);
509  FillHistogram(Form("hSubleadingJetEta%s", appendix.Data()), jet->Eta(), fCent);
510  FillHistogram(Form("hSubleadingJetEtaPt%s", appendix.Data()), jet->Eta(), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area());
511  FillHistogram(Form("hSubleadingJetPhiPt%s", appendix.Data()), jet->Phi(), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area());
512  FillHistogram(Form("hSubleadingJetPhiEta%s", appendix.Data()), jet->Phi(), jet->Eta());
513  FillHistogram(Form("hSubleadingJetArea%s", appendix.Data()), jet->Area(), fCent);
514  FillHistogram(Form("hSubleadingJetAreaPt%s", appendix.Data()), jet->Area(), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area());
515  FillHistogram(Form("hSubleadingJetPtLeadingHadron%s", appendix.Data()), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), fJetsCont->GetLeadingHadronPt(jet));
516  }
517 
518  // ####### Jet constituent plots
519  for(Int_t i = 0; i < jet->GetNumberOfTracks(); i++)
520  {
521  AliVParticle* constituent = static_cast<AliVParticle*>(jet->TrackAt(i, fTracksCont->GetArray()));
522  if(!constituent)
523  continue;
524 
525  Bool_t filterConditionFulfilled = kFALSE;
526  AliAODTrack* aodTrack = static_cast<AliAODTrack*>(constituent);
527  if (aodTrack)
528  filterConditionFulfilled = aodTrack->TestFilterBit(fConstPtFilterBit);
529 
530  // Fill jet constituent plots
531  FillHistogram(Form("hJetConstituentPt_Cent0_100%s", appendix.Data()), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), constituent->Pt());
532  if(filterConditionFulfilled)
533  FillHistogram(Form("hJetConstituentPt_Cent0_100_FilterBit%i%s", fConstPtFilterBit, appendix.Data()), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), constituent->Pt());
534  if( (fCent >= 0) && (fCent < 10) )
535  {
536  FillHistogram(Form("hJetConstituentPt_Cent0_10%s", appendix.Data()), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), constituent->Pt());
537  if(filterConditionFulfilled)
538  FillHistogram(Form("hJetConstituentPt_Cent0_10_FilterBit%i%s", fConstPtFilterBit, appendix.Data()), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), constituent->Pt());
539  }
540  else if( (fCent >= 10) && (fCent < 30) )
541  {
542  FillHistogram(Form("hJetConstituentPt_Cent10_30%s", appendix.Data()), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), constituent->Pt());
543  if(filterConditionFulfilled)
544  FillHistogram(Form("hJetConstituentPt_Cent10_30_FilterBit%i%s", fConstPtFilterBit, appendix.Data()), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), constituent->Pt());
545  }
546  else if( (fCent >= 30) && (fCent < 50) )
547  {
548  FillHistogram(Form("hJetConstituentPt_Cent30_50%s", appendix.Data()), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), constituent->Pt());
549  if(filterConditionFulfilled)
550  FillHistogram(Form("hJetConstituentPt_Cent30_50_FilterBit%i%s", fConstPtFilterBit, appendix.Data()), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), constituent->Pt());
551  }
552  else if( (fCent >= 50) && (fCent < 90) )
553  {
554  FillHistogram(Form("hJetConstituentPt_Cent50_90%s", appendix.Data()), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), constituent->Pt());
555  if(filterConditionFulfilled)
556  FillHistogram(Form("hJetConstituentPt_Cent50_90_FilterBit%i%s", fConstPtFilterBit, appendix.Data()), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), constituent->Pt());
557  }
558  }
559 
560  FillHistogram(Form("hJetConstituentCount_Cent0_100%s", appendix.Data()), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), jet->GetNumberOfTracks());
561  if( (fCent >= 0) && (fCent < 10) )
562  FillHistogram(Form("hJetConstituentCount_Cent0_10%s", appendix.Data()), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), jet->GetNumberOfTracks());
563 
564  FillHistogram(Form("hBackgroundPtConstCount_Cent0_100%s", appendix.Data()), fJetsCont->GetRhoVal(), jet->GetNumberOfTracks());
565  if( (fCent >= 0) && (fCent < 10) )
566  FillHistogram(Form("hBackgroundPtConstCount_Cent0_10%s", appendix.Data()), fJetsCont->GetRhoVal(), jet->GetNumberOfTracks());
567 
568  // ####### Embedding plots
569  if( (fJetOutputMode == 4) || (fJetOutputMode == 5))
570  {
571  AliEmcalJet* refJet = GetReferenceJet(jet);
572  Double_t deltaEta = (jet->Eta()-refJet->Eta());
573  Double_t deltaPhi = TMath::Min(TMath::Abs(jet->Phi()-refJet->Phi()),TMath::TwoPi() - TMath::Abs(jet->Phi()-refJet->Phi()));
574  if(jet->Phi() < refJet->Phi())
575  deltaPhi = -deltaPhi;
576 
577  Double_t deltaR = TMath::Sqrt(deltaEta*deltaEta + deltaPhi*deltaPhi);
578  FillHistogram(Form("hEmbeddingDeltaR%s", appendix.Data()), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), deltaR);
579  FillHistogram(Form("hEmbeddingDeltaEta%s", appendix.Data()), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), deltaPhi);
580  FillHistogram(Form("hEmbeddingDeltaPhi%s", appendix.Data()), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), deltaEta);
581  FillHistogram(Form("hEmbeddingJetPt%s", appendix.Data()), refJet->Pt());
582  FillHistogram(Form("hEmbeddingJetPhiEta%s", appendix.Data()), refJet->Phi(), refJet->Eta());
583 
584  // Only create 3D plots for each cut on demand
585  Double_t trackRatio = 0.;
586  Double_t ptRatio = 0.;
587  GetTrackMCRatios(jet, refJet, trackRatio, ptRatio);
588 
589  if(fCent >= 0 && fCent < 10)
590  {
591  if((appendix == "") || fJetEmbeddingCreatePtPlotPerCut)
592  FillHistogram3D(Form("hEmbeddingPtCorr010%s", appendix.Data()), refJet->Pt(), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), ptRatio);
593  if(ptRatio >= 0.2)
594  FillHistogram(Form("hEmbeddingPtCorr010_Above20%s", appendix.Data()), refJet->Pt(), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area());
595  }
596  else if (fCent >= 10 && fCent < 30)
597  {
598  if((appendix == "") || fJetEmbeddingCreatePtPlotPerCut)
599  FillHistogram3D(Form("hEmbeddingPtCorr1030%s", appendix.Data()), refJet->Pt(), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), ptRatio);
600  if(ptRatio >= 0.2)
601  FillHistogram(Form("hEmbeddingPtCorr1030_Above20%s", appendix.Data()), refJet->Pt(), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area());
602  }
603  else if (fCent >= 30 && fCent < 50)
604  {
605  if((appendix == "") || fJetEmbeddingCreatePtPlotPerCut)
606  FillHistogram3D(Form("hEmbeddingPtCorr3050%s", appendix.Data()), refJet->Pt(), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), ptRatio);
607  if(ptRatio >= 0.2)
608  FillHistogram(Form("hEmbeddingPtCorr3050_Above20%s", appendix.Data()), refJet->Pt(), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area());
609  }
610  else if (fCent >= 50 && fCent < 90)
611  {
612  if((appendix == "") || fJetEmbeddingCreatePtPlotPerCut)
613  FillHistogram3D(Form("hEmbeddingPtCorr5090%s", appendix.Data()), refJet->Pt(), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area(), ptRatio);
614  if(ptRatio >= 0.2)
615  FillHistogram(Form("hEmbeddingPtCorr5090_Above20%s", appendix.Data()), refJet->Pt(), jet->Pt() - fJetsCont->GetRhoVal()*jet->Area());
616  }
617  }
618 }
619 
620 //________________________________________________________________________
622 {
623  FillHistogram("hTrackPt", track->Pt(), fCent);
624  FillHistogram("hTrackPhi", track->Phi(), fCent);
625  FillHistogram("hTrackEta", track->Eta(), fCent);
626  FillHistogram("hTrackEtaPt", track->Eta(), track->Pt());
627  FillHistogram("hTrackPhiPt", track->Phi(), track->Pt());
628  FillHistogram("hTrackPhiEta", track->Phi(), track->Eta());
629 }
630 
631 //________________________________________________________________________
633 {
634  Double_t tmpPt = jet->Pt() - fJetsCont->GetRhoVal()*jet->Area();
635  new ((*fJetsOutput.at(arrayIndex))[jetsAlreadyInArray]) AliPicoTrack(tmpPt, jet->Eta(), jet->Phi(), jet->Charge(), 0, 0);
636  jetsAlreadyInArray++;
637 }
638 
639 //________________________________________________________________________
641 {
642  if(fTrackParticleArrayName != "")
643  {
644  new ((*fTracksOutput)[fAcceptedTracks]) AliPicoTrack(track->Pt(), track->Eta(), track->Phi(), track->Charge(), 0, 0); // only Pt,Eta,Phi are interesting for correlations;
645  fAcceptedTracks++;
646  }
647 }
648 
649 //________________________________________________________________________
651 {
652  // Only allow when we have aod tracks
653  AliAODTrack* aodtrack = dynamic_cast<AliAODTrack*>(track);
654  if(!aodtrack)
655  return;
656 
657  // Discard tracks according to their pT (below 20 GeV)
658  if(track->Pt() < 20.)
659  if(fRandom->Rndm() >= TMath::Power((1/20. * track->Pt()), fTrackExtractionPercentagePower))
660  return;
661 
662  // Create basic particle from track and extract pid object
663  fTreeBufferPID = aodtrack->GetDetPid();
664 
665  Int_t truthPID = 0;
666 
667  // Get truth values if we are on MC
668  TClonesArray* fTruthParticleArray = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("mcparticles"));
669  if(fTruthParticleArray)
670  {
671  for(Int_t i=0; i<fTruthParticleArray->GetEntries();i++)
672  {
673  AliAODMCParticle* mcParticle = dynamic_cast<AliAODMCParticle*>(fTruthParticleArray->At(i));
674  if(!mcParticle) continue;
675 
676  if (mcParticle->GetLabel() == aodtrack->GetLabel())
677  {
678  truthPID = mcParticle->PdgCode();
679  break;
680  }
681  }
682  }
683 
684  fTreeBufferPDG = truthPID;
685 
686  AliBasicParticle basicParticle(aodtrack->Eta(), aodtrack->Phi(), aodtrack->Pt(), aodtrack->Charge());
687  fTreeBufferTrack = &basicParticle;
688 
689  fTracksTree->Fill();
690 }
691 
692 //________________________________________________________________________
694 {
695  // Check jet pT threshold
697  return;
698 
699  // Discard jets statistically
700  if(fRandom->Rndm() >= fEventExtractionPercentage)
701  return;
702 
703  static Int_t numSavedEvents = 0;
704  numSavedEvents++;
705 
706 
707  AddHistogram2D<TH2D>(Form("Event%i", numSavedEvents), "Event display", "COLZ", 180, 0., 2*TMath::Pi(), 100, -2.5, 2.5, "#phi", "#eta", "dN^{Tracks}/d#phi d#eta");
708  fTracksCont->ResetCurrentID();
709  while(AliVTrack *track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle()))
710  FillHistogram(Form("Event%i", numSavedEvents), track->Phi(), track->Eta(), track->Pt());
711 
712 }
713 
714 //________________________________________________________________________
716 {
718 
719  // Jet veto:
720  // 1. jet-by-jet mode: veto is active if jets overlaps with a suitable jet
721  // 2. other mode: veto is active if suitable jet is in sample
722  AliEmcalJet* vetoJet = 0;
723  if(!fJetVetoJetByJet)
724  vetoJet = GetLeadingVetoJet();
725 
726  // ####### Jet loop
727  fAcceptedJets = 0;
728  for(Int_t i = 0; i<fJetEmbeddingCuts.size(); i++)
729  fJetEmbeddingCuts.at(i).fAcceptedJets = 0;
730  fJetsCont->ResetCurrentID();
731  while(AliEmcalJet *jet = fJetsCont->GetNextAcceptJet())
732  {
733  if(!IsJetSelected(jet))
734  continue;
735 
736  // Plots + output jets regardless of embedding cut
737  FillHistogramsJets(jet, 0);
739 
740  // Plots + output jets for each embedding cut
741  if(fJetVetoJetByJet)
742  vetoJet = GetVetoJet(jet);
743  for(Int_t i = 0; i<fJetEmbeddingCuts.size(); i++)
744  {
746  AliEmcalJet* refJet = GetReferenceJet(jet);
747  Double_t trackRatio = 0.;
748  Double_t ptRatio = 0.;
749  GetTrackMCRatios(jet, refJet, trackRatio, ptRatio);
750 
751  Double_t vetoJetPt = 0.;
752  if(vetoJet)
753  vetoJetPt = vetoJet->Pt() - vetoJet->Area()*fJetsCont->GetRhoVal();
754 
755  if(!currentCut.IsCutFulfilled(jet->Pt(), refJet->Pt(), fCent, ptRatio, vetoJetPt))
756  continue;
757 
758  FillHistogramsJets(jet, currentCut.fCutName.Data());
759  AddJetToOutputArray(jet, currentCut.fArrayIndex, currentCut.fAcceptedJets);
760  }
761  }
762 
763  // ####### Particle loop
764  fAcceptedTracks = 0;
765  fTracksCont->ResetCurrentID();
766  Int_t trackcount = 0;
767  while(AliVTrack *track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle()))
768  {
769  // Track plots
770  FillHistogramsTracks(track);
771  // Add track to output array
772  trackcount++;
773  AddTrackToOutputArray(track);
775  AddTrackToTree(track);
776  }
777 
778  // Add event to output tree
780  AddEventToTree();
781 
782  // ######### Random cone sampling
783  for(Int_t iCone=0; iCone<fNumRandomConesPerEvent; iCone++)
784  {
785  // Throw random cone
786  Double_t tmpRandConeEta = fJetsCont->GetJetEtaMin() + fRandom->Rndm()*TMath::Abs(fJetsCont->GetJetEtaMax()-fJetsCont->GetJetEtaMin());
787  Double_t tmpRandConePhi = fRandom->Rndm()*TMath::TwoPi();
788  Double_t tmpRandConePt = 0;
789  Double_t tmpRandConePt3GeV = 0;
790  // Fill pT that is in cone
791  fTracksCont->ResetCurrentID();
792  while(AliVTrack *track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle()))
793  if(IsTrackInCone(track, tmpRandConeEta, tmpRandConePhi, fJetsCont->GetJetRadius()))
794  {
795  tmpRandConePt += track->Pt();
796  if (track->Pt() > 3.0)
797  tmpRandConePt3GeV += track->Pt();
798  }
799  // Fill histograms
800  FillHistogram("hRandomConePt", tmpRandConePt - fJetsCont->GetRhoVal()*fJetsCont->GetJetRadius()*fJetsCont->GetJetRadius()*TMath::Pi(), fCent);
801  FillHistogram("hRandomConePtCut3GeV", tmpRandConePt3GeV - fJetsCont->GetRhoVal()*fJetsCont->GetJetRadius()*fJetsCont->GetJetRadius()*TMath::Pi(), fCent);
802  FillHistogram("hRandomConeRawPt", tmpRandConePt, fCent);
803  FillHistogram("hRandomConeRawPtCut3GeV", tmpRandConePt3GeV, fCent);
804  }
805 
806  // ####### Event properties
807  FillHistogram("hBackgroundPt", fJetsCont->GetRhoVal(), fCent);
808  FillHistogram("hJetCount", fAcceptedJets, fCent);
809  FillHistogram("hTrackCount", trackcount, fCent);
810  // NOTE: It is possible to use fTracksCont->GetLeadingParticle() since we do not apply additional track cuts
811  AliVTrack* leadTrack = static_cast<AliVTrack*>(fTracksCont->GetLeadingParticle());
812  if(leadTrack)
813  {
814  FillHistogram("hLeadingTrackPt", leadTrack->Pt(), fCent);
815  FillHistogram("hLeadingTrackPhi", leadTrack->Phi(), fCent);
816  FillHistogram("hLeadingTrackEta", leadTrack->Eta(), fCent);
817  FillHistogram("hLeadingTrackPhiEta", leadTrack->Phi(), leadTrack->Eta());
818  }
819 
820  return kTRUE;
821 }
822 
823 //########################################################################
824 // HELPERS
825 //########################################################################
826 
827 //________________________________________________________________________
829 {
830  if(!fPythiaInfo)
831  AliError("fPythiaInfo object not available. Is it activated with SetGeneratePythiaInfoObject()?");
832 
833  Double_t bestMatchDeltaR1 = 999.;
834  Double_t bestMatchDeltaR2 = 999.;
835 
836  fJetsCont->ResetCurrentID();
837  while(AliEmcalJet *jet = fJetsCont->GetNextAcceptJet())
838  {
839  // Check via geometrical matching if jet is connected to the initial collision
840  Double_t deltaEta1 = TMath::Abs(jet->Eta()-fPythiaInfo->GetPartonEta6());
841  Double_t deltaEta2 = TMath::Abs(jet->Eta()-fPythiaInfo->GetPartonEta7());
842  Double_t deltaPhi1 = TMath::Min(TMath::Abs(jet->Phi()-fPythiaInfo->GetPartonPhi6()),TMath::TwoPi() - TMath::Abs(jet->Phi()-fPythiaInfo->GetPartonPhi6()));
843  Double_t deltaPhi2 = TMath::Min(TMath::Abs(jet->Phi()-fPythiaInfo->GetPartonPhi7()),TMath::TwoPi() - TMath::Abs(jet->Phi()-fPythiaInfo->GetPartonPhi7()));
844 
845  Double_t deltaR1 = TMath::Sqrt(deltaEta1*deltaEta1 + deltaPhi1*deltaPhi1);
846  Double_t deltaR2 = TMath::Sqrt(deltaEta2*deltaEta2 + deltaPhi2*deltaPhi2);
847 
848  if(deltaR1 < bestMatchDeltaR1)
849  {
850  bestMatchDeltaR1 = deltaR1;
852  }
853  if(deltaR2 < bestMatchDeltaR2)
854  {
855  bestMatchDeltaR2 = deltaR2;
857  }
858  }
859 
860  if(bestMatchDeltaR1 > 0.3)
862  if(bestMatchDeltaR2 > 0.3)
864 }
865 
866 //________________________________________________________________________
868 {
869  fMatchedJets.clear();
870  fMatchedJetsReference.clear();
871 
872  AliEmcalJet* jetLeading = 0;
873  AliEmcalJet* jetSubLeading = 0;
874  GetLeadingJetsInArray(fJetEmbeddingArray, "rho", jetLeading, jetSubLeading);
875 
876  // ############ Search for all matches with leading/subleading jets
877  for(Int_t i=0; i<fJetEmbeddingNumMatchedJets; i++)
878  {
879  AliEmcalJet* probeJet = 0;
880  if(i==0)
881  probeJet = jetLeading;
882  else if(i==1)
883  probeJet = jetSubLeading;
884 
885  if(!probeJet)
886  continue;
887 
888  AliEmcalJet* matchedJet = 0;
889  AliEmcalJet* matchedJetReference = 0;
890  Double_t bestMatchDeltaR = 999.;
891  fJetsCont->ResetCurrentID();
892  // Loop over all embedded jets to find the best match
893  while(AliEmcalJet* embeddedJet = fJetsCont->GetNextAcceptJet())
894  {
895  Double_t deltaEta = (embeddedJet->Eta()-probeJet->Eta());
896  Double_t deltaPhi = TMath::Min(TMath::Abs(embeddedJet->Phi()-probeJet->Phi()),TMath::TwoPi() - TMath::Abs(embeddedJet->Phi()-probeJet->Phi()));
897  Double_t deltaR = TMath::Sqrt(deltaEta*deltaEta + deltaPhi*deltaPhi);
898 
899  // Cut jets too far away
900  if (deltaR > fJetEmbeddingMaxDistance)
901  continue;
902 
903  // Search for the best match
904  if(deltaR < bestMatchDeltaR)
905  {
906  bestMatchDeltaR = deltaR;
907  matchedJet = embeddedJet;
908  matchedJetReference = probeJet;
909  }
910  }
911  // Put matched jet to a list
912  if(matchedJet && matchedJetReference)
913  {
914  fMatchedJets.push_back(matchedJet);
915  fMatchedJetsReference.push_back(matchedJetReference);
916  }
917  }
918 }
919 
920 
921 //________________________________________________________________________
923 {
924  Int_t tracksFromMC = 0;
925  Int_t tracksTotal = 0;
926  Double_t ptFromMC = 0;
927  Double_t ptTotal = 0;
928 
929  if(fJetEmbeddingUsePerTrackMCPercentage) // Calculate MC track percentage from tracks of the matched jet
930  {
931  for(Int_t j = 0; j < jet->GetNumberOfTracks(); j++)
932  {
933  AliVParticle* constituent = static_cast<AliVParticle*>(jet->TrackAt(j, fTracksCont->GetArray()));
934  if(!constituent)
935  continue;
936 
937  TClonesArray* mcArray = static_cast<TClonesArray*>(InputEvent()->FindListObject(Form("%s", fJetEmbeddingTrackArrayName.Data())));
938  Bool_t foundInMC = kFALSE;
939  for(Int_t k=0; k<mcJet->GetNumberOfTracks(); k++)
940  {
941  AliVParticle* mcConstituent = static_cast<AliVParticle*>(mcJet->TrackAt(k, mcArray));
942  if(!mcConstituent)
943  continue;
944  if(mcConstituent->GetLabel() == constituent->GetLabel())
945  {
946  foundInMC = kTRUE;
947  break;
948  }
949  }
950 
951  if(foundInMC)
952  {
953  tracksFromMC++;
954  ptFromMC += constituent->Pt();
955  }
956  tracksTotal++;
957  ptTotal += constituent->Pt();
958 
959  }
960  }
961  else // Calculate MC track percentage from all tracks in MC
962  {
963  for(Int_t j = 0; j < jet->GetNumberOfTracks(); j++)
964  {
965  AliVParticle* constituent = static_cast<AliVParticle*>(jet->TrackAt(j, fTracksCont->GetArray()));
966  if(!constituent)
967  continue;
968 
969  // Plots on MC percentage in jet
970  if(constituent->GetLabel() > 10000)
971  {
972  tracksFromMC++;
973  ptFromMC += constituent->Pt();
974  }
975  tracksTotal++;
976  ptTotal += constituent->Pt();
977  }
978  }
979 
980 
981  trackRatio = 0.;
982  if(tracksTotal)
983  trackRatio = tracksFromMC/((Double_t)tracksTotal);
984 
986  ptTotal = jet->Pt() - fJetsCont->GetRhoVal()*jet->Area();
987 
988  ptRatio = 0.;
989  if(ptTotal)
990  ptRatio = ptFromMC/ptTotal;
991 
992 }
993 
994 //________________________________________________________________________
996 {
997  Double_t leadingVetoJetPt = -999.;
998  AliEmcalJet* leadingVetoJet = 0;
999  // Search for the 'leading' overlapping jet in the veto sample
1000  if(fJetVetoArray && jet)
1001  {
1002  for(Int_t j=0; j<fJetVetoArray->GetEntries(); j++)
1003  {
1004  UInt_t dummy = 0;
1005  AliEmcalJet* vetoJet = static_cast<AliEmcalJet*>(fJetVetoArray->At(j));
1006  // Check if veto jet is accepted
1007  if(!fJetsCont->AcceptJet(vetoJet , dummy))
1008  continue;
1009 
1010  // Check matching distance
1011  Double_t vetoPt = vetoJet->Pt() - vetoJet->Area()*fJetsCont->GetRhoVal();
1012  Double_t deltaEta = (jet->Eta()-vetoJet->Eta());
1013  Double_t deltaPhi = TMath::Min(TMath::Abs(jet->Phi()-vetoJet->Phi()),TMath::TwoPi() - TMath::Abs(jet->Phi()-vetoJet->Phi()));
1014  Double_t deltaR = TMath::Sqrt(deltaEta*deltaEta + deltaPhi*deltaPhi);
1015 
1016  if ((vetoPt > leadingVetoJetPt) && (deltaR <= fJetEmbeddingMaxDistance))
1017  {
1018  leadingVetoJetPt = vetoPt;
1019  leadingVetoJet = vetoJet;
1020  }
1021  }
1022  }
1023 
1024  return leadingVetoJet;
1025 }
1026 
1027 //________________________________________________________________________
1029 {
1030  Double_t leadingVetoJetPt = -999.;
1031  AliEmcalJet* leadingVetoJet = 0;
1032  // Search for the 'leading' jet in the veto sample
1033  if(fJetVetoArray)
1034  {
1035  for(Int_t j=0; j<fJetVetoArray->GetEntries(); j++)
1036  {
1037  UInt_t dummy = 0;
1038  AliEmcalJet* vetoJet = static_cast<AliEmcalJet*>(fJetVetoArray->At(j));
1039  // Check if veto jet is accepted
1040  if(!fJetsCont->AcceptJet(vetoJet , dummy))
1041  continue;
1042 
1043  // Check matching distance
1044  Double_t vetoPt = vetoJet->Pt() - vetoJet->Area()*fJetsCont->GetRhoVal();
1045  if (vetoPt > leadingVetoJetPt)
1046  {
1047  leadingVetoJetPt = vetoPt;
1048  leadingVetoJet = vetoJet;
1049  }
1050  }
1051  }
1052  return leadingVetoJet;
1053 }
1054 
1055 //________________________________________________________________________
1057 {
1058  // This is to use a full cone in phi even at the edges of phi (2pi -> 0) (0 -> 2pi)
1059  Double_t trackPhi = 0.0;
1060  if (track->Phi() > (TMath::TwoPi() - (radius-phi)))
1061  trackPhi = track->Phi() - TMath::TwoPi();
1062  else if (track->Phi() < (phi+radius - TMath::TwoPi()))
1063  trackPhi = track->Phi() + TMath::TwoPi();
1064  else
1065  trackPhi = track->Phi();
1066 
1067  if ( TMath::Abs(trackPhi-phi)*TMath::Abs(trackPhi-phi) + TMath::Abs(track->Eta()-eta)*TMath::Abs(track->Eta()-eta) <= radius*radius)
1068  return kTRUE;
1069 
1070  return kFALSE;
1071 }
1072 
1073 
1074 //________________________________________________________________________
1076 {
1077  // Calculate leading + subleading jet
1079  if(fJetOutputMode==6)
1081  else if(fJetOutputMode==4 || fJetOutputMode==5)
1082  GetMatchingJets();
1083 }
1084 
1085 //________________________________________________________________________
1086 void AliAnalysisTaskChargedJetsHadronCF::GetLeadingJets(const char* opt, AliEmcalJet*& jetLeading, AliEmcalJet*& jetSubLeading)
1087 {
1088  // Customized from AliJetContainer::GetLeadingJet()
1089  // Get the leading+subleading jet; if opt contains "rho" the sorting is according to pt-A*rho
1090 
1091  TString option(opt);
1092  option.ToLower();
1093 
1094  jetLeading = 0;
1095  jetSubLeading = 0;
1096 
1097  fJetsCont->ResetCurrentID();
1098  Double_t tmpLeadingPt = 0;
1099  Double_t tmpSubleadingPt = 0;
1100 
1101  if (option.Contains("rho")) {
1102  while (AliEmcalJet* jet = fJetsCont->GetNextAcceptJet()) {
1103  if ( (jet->Pt()-jet->Area()*fJetsCont->GetRhoVal()) > tmpLeadingPt )
1104  {
1105  jetSubLeading = jetLeading;
1106  jetLeading = jet;
1107  tmpSubleadingPt = tmpLeadingPt;
1108  tmpLeadingPt = jet->Pt()-jet->Area()*fJetsCont->GetRhoVal();
1109  }
1110  else if ( (jet->Pt()-jet->Area()*fJetsCont->GetRhoVal()) > tmpSubleadingPt )
1111  {
1112  jetSubLeading = jet;
1113  tmpSubleadingPt = jet->Pt()-jet->Area()*fJetsCont->GetRhoVal();
1114  }
1115  }
1116  }
1117  else {
1118  while (AliEmcalJet* jet = fJetsCont->GetNextAcceptJet()) {
1119  if ( (jet->Pt()) > tmpLeadingPt )
1120  {
1121  jetSubLeading = jetLeading;
1122  jetLeading = jet;
1123  tmpSubleadingPt = tmpLeadingPt;
1124  tmpLeadingPt = jet->Pt();
1125  }
1126  else if ( (jet->Pt()) > tmpSubleadingPt )
1127  {
1128  jetSubLeading = jet;
1129  tmpSubleadingPt = jet->Pt();
1130  }
1131  }
1132  }
1133 }
1134 
1135 //________________________________________________________________________
1136 void AliAnalysisTaskChargedJetsHadronCF::GetLeadingJetsInArray(TClonesArray* arr, const char* opt, AliEmcalJet*& jetLeading, AliEmcalJet*& jetSubLeading)
1137 {
1138  // Customized from AliJetContainer::GetLeadingJet()
1139  // Get the leading+subleading jet; if opt contains "rho" the sorting is according to pt-A*rho
1140 
1141  TString option(opt);
1142  option.ToLower();
1143 
1144  jetLeading = 0;
1145  jetSubLeading = 0;
1146  Double_t tmpLeadingPt = -999.;
1147  Double_t tmpSubleadingPt = -999.;
1148 
1149  for(Int_t i=0; i<arr->GetEntries(); i++)
1150  {
1151  AliEmcalJet* jet = static_cast<AliEmcalJet*>(arr->At(i));
1152  UInt_t dummy = 0;
1153  if(!fJetsCont->AcceptJet(jet , dummy))
1154  continue;
1155 
1156  Double_t jetPt = jet->Pt();
1157  if (option.Contains("rho"))
1158  jetPt -= jet->Area()*fJetsCont->GetRhoVal();
1159 
1160  if ( jetPt > tmpLeadingPt )
1161  {
1162  jetSubLeading = jetLeading;
1163  jetLeading = jet;
1164  tmpSubleadingPt = tmpLeadingPt;
1165  tmpLeadingPt = jetPt;
1166  }
1167  else if ( jetPt > tmpSubleadingPt )
1168  {
1169  jetSubLeading = jet;
1170  tmpSubleadingPt = jetPt;
1171  }
1172 
1173  }
1174 }
1175 
1176 
1177 //________________________________________________________________________
1179 {
1180  // Method for the correct logarithmic binning of histograms
1181  TAxis *axis = h->GetAxis(axisNumber);
1182  int bins = axis->GetNbins();
1183 
1184  Double_t from = axis->GetXmin();
1185  Double_t to = axis->GetXmax();
1186  Double_t *newBins = new Double_t[bins + 1];
1187 
1188  newBins[0] = from;
1189  Double_t factor = pow(to/from, 1./bins);
1190 
1191  for (int i = 1; i <= bins; i++) {
1192  newBins[i] = factor * newBins[i-1];
1193  }
1194  axis->Set(bins, newBins);
1195  delete [] newBins;
1196 }
1197 
1198 //________________________________________________________________________
1200 {
1201  std::vector<AliEmcalJet*>::iterator matchedJetFindResult = std::find(fMatchedJets.begin(), fMatchedJets.end(), jet);
1202  if(matchedJetFindResult == fMatchedJets.end())
1203  {
1204  AliError("Checked for a reference jet but it was not found. Check code.");
1205  return 0;
1206  }
1207 
1208  Int_t matchedJetIndex = (matchedJetFindResult - fMatchedJets.begin());
1209  AliEmcalJet* refJet = fMatchedJetsReference[matchedJetIndex];
1210  return refJet;
1211 }
1212 
1213 //________________________________________________________________________
1215 {
1216  TH1* tmpHist = static_cast<TH1*>(fOutput->FindObject(key));
1217  if(!tmpHist)
1218  {
1219  AliError(Form("Cannot find histogram <%s> ",key)) ;
1220  return;
1221  }
1222 
1223  tmpHist->Fill(x);
1224 }
1225 
1226 //________________________________________________________________________
1228 {
1229  TH1* tmpHist = static_cast<TH1*>(fOutput->FindObject(key));
1230  if(!tmpHist)
1231  {
1232  AliError(Form("Cannot find histogram <%s> ",key));
1233  return;
1234  }
1235 
1236  if (tmpHist->IsA()->GetBaseClass("TH1"))
1237  static_cast<TH1*>(tmpHist)->Fill(x,y); // Fill x with y
1238  else if (tmpHist->IsA()->GetBaseClass("TH2"))
1239  static_cast<TH2*>(tmpHist)->Fill(x,y); // Fill x,y with 1
1240 }
1241 
1242 //________________________________________________________________________
1244 {
1245  TH2* tmpHist = static_cast<TH2*>(fOutput->FindObject(key));
1246  if(!tmpHist)
1247  {
1248  AliError(Form("Cannot find histogram <%s> ",key));
1249  return;
1250  }
1251 
1252  tmpHist->Fill(x,y,add);
1253 }
1254 
1255 //________________________________________________________________________
1257 {
1258  TH3* tmpHist = static_cast<TH3*>(fOutput->FindObject(key));
1259  if(!tmpHist)
1260  {
1261  AliError(Form("Cannot find histogram <%s> ",key));
1262  return;
1263  }
1264 
1265  if(add)
1266  tmpHist->Fill(x,y,z,add);
1267  else
1268  tmpHist->Fill(x,y,z);
1269 }
1270 
1271 
1272 //________________________________________________________________________
1273 template <class T> T* AliAnalysisTaskChargedJetsHadronCF::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)
1274 {
1275  T* tmpHist = new T(name, title, xBins, xMin, xMax);
1276 
1277  tmpHist->GetXaxis()->SetTitle(xTitle);
1278  tmpHist->GetYaxis()->SetTitle(yTitle);
1279  tmpHist->SetOption(options);
1280  tmpHist->SetMarkerStyle(kFullCircle);
1281  tmpHist->Sumw2();
1282 
1283  fOutput->Add(tmpHist);
1284 
1285  return tmpHist;
1286 }
1287 
1288 //________________________________________________________________________
1289 template <class T> T* AliAnalysisTaskChargedJetsHadronCF::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)
1290 {
1291  T* tmpHist = new T(name, title, xBins, xMin, xMax, yBins, yMin, yMax);
1292  tmpHist->GetXaxis()->SetTitle(xTitle);
1293  tmpHist->GetYaxis()->SetTitle(yTitle);
1294  tmpHist->GetZaxis()->SetTitle(zTitle);
1295  tmpHist->SetOption(options);
1296  tmpHist->SetMarkerStyle(kFullCircle);
1297  tmpHist->Sumw2();
1298 
1299  fOutput->Add(tmpHist);
1300 
1301  return tmpHist;
1302 }
1303 
1304 //________________________________________________________________________
1305 template <class T> T* AliAnalysisTaskChargedJetsHadronCF::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)
1306 {
1307  T* tmpHist = new T(name, title, xBins, xMin, xMax, yBins, yMin, yMax, zBins, zMin, zMax);
1308  tmpHist->GetXaxis()->SetTitle(xTitle);
1309  tmpHist->GetYaxis()->SetTitle(yTitle);
1310  tmpHist->GetZaxis()->SetTitle(zTitle);
1311  tmpHist->SetOption(options);
1312  tmpHist->SetMarkerStyle(kFullCircle);
1313  tmpHist->Sumw2();
1314 
1315  fOutput->Add(tmpHist);
1316 
1317  return tmpHist;
1318 }
1319 
1320 //________________________________________________________________________
1322 {
1323  // Called once at the end of the analysis.
1324 }
1325 
Int_t fArrayIndex
array index that holds the output array index
Double_t fEventExtractionMaxJetPt
maximum jet pt of recorded events
Short_t Charge() const
Definition: AliEmcalJet.h:116
Double_t Area() const
Definition: AliEmcalJet.h:123
void ExecOnce()
Perform steps needed to initialize the analysis.
Double_t GetRhoVal() const
double Double_t
Definition: External.C:58
Bool_t fJetVetoJetByJet
If true, the jet veto will be applied on a jet-by-jet basis.
Double_t GetJetEtaMin() const
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")
const char * title
Definition: MakeQAPdf.C:26
void GetLeadingJets(const char *opt, AliEmcalJet *&jetLeading, AliEmcalJet *&jetSubLeading)
Float_t GetPartonEta6() const
AliEmcalPythiaInfo * fPythiaInfo
!event parton info
AliJetContainer * GetJetContainer(Int_t i=0) const
Definition: External.C:244
Float_t GetPartonEta7() const
std::vector< AliChargedJetsHadronCFCuts > fJetEmbeddingCuts
Cuts used in jet embedding.
Double_t Eta() const
Definition: AliEmcalJet.h:114
Double_t Phi() const
Definition: AliEmcalJet.h:110
Container with name, TClonesArray and cuts for particles.
Support task for (charged) jet-hadron correlations.
AliBasicParticle * fTreeBufferTrack
! Tree of extracted jets (buffer)
Double_t GetJetEtaMax() const
Float_t GetPartonPhi7() const
virtual Bool_t AcceptJet(Int_t i, UInt_t &rejectionReason) const
Int_t fJetOutputMode
mode which jets are written to array (0: all accepted, 1: leading, 2: subleading, 3: leading+subleadi...
Double_t fEventExtractionMinJetPt
minimum jet pt of recorded events
std::vector< AliEmcalJet * > fMatchedJetsReference
Jets matched in an event (reference)
Double_t fEventExtractionPercentage
percentage of events that is recorded
std::vector< TClonesArray * > fJetsOutput
! vector of arrays of basic correlation particles attached to the event (jets)
AliEmcalJet * fSubleadingJet
! subleading jet (calculated event-by-event)
void FillHistogram3D(const char *key, Double_t x, Double_t y, Double_t z, Double_t add=0)
Float_t GetPartonPhi6() const
Int_t TrackAt(Int_t idx) const
Definition: AliEmcalJet.h:153
Container class of cuts for AliAnalysisTaskChargedJetsHadronCF.
UShort_t GetNumberOfTracks() const
Definition: AliEmcalJet.h:132
UShort_t T(UShort_t m, UShort_t t)
Definition: RingBits.C:60
void GetLeadingJetsInArray(TClonesArray *arr, const char *opt, AliEmcalJet *&jetLeading, AliEmcalJet *&jetSubLeading)
TClonesArray * fJetEmbeddingArray
! Array of generated jets imported into task (for embedding)
AliParticleContainer * GetParticleContainer(Int_t i=0) const
Get particle container attached to this task.
void GetTrackMCRatios(AliEmcalJet *jet, AliEmcalJet *mcJet, Double_t &trackRatio, Double_t &ptRatio)
AliParticleContainer * GetParticleContainer() const
Double_t fTrackExtractionPercentagePower
Extraction percentage for tracks.
Bool_t fJetEmbeddingCreatePtPlotPerCut
create TH3 per cut or only once
int Int_t
Definition: External.C:63
AliAODPid * fTreeBufferPID
! Tree of extracted jets (buffer)
unsigned int UInt_t
Definition: External.C:33
TString fJetParticleArrayName
Name of fJetsOutput array (if one uses only one)
Double_t GetLeadingHadronPt(const AliEmcalJet *jet) const
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")
TString fJetEmbeddingArrayName
Name of array used to match jets.
Int_t fNumberOfCentralityBins
Number of centrality bins.
Int_t fTreeBufferPDG
! Tree of extracted jets (buffer)
virtual AliVParticle * GetLeadingParticle(const char *opt="")
std::vector< AliEmcalJet * > fMatchedJets
Jets matched in an event (embedded)
TString fJetEmbeddingTrackArrayName
Name of array used to match tracks of jets.
Double_t fCent
!event centrality
TString fTrackParticleArrayName
Name of fTracksOutput array.
TClonesArray * fTracksOutput
! Array of basic correlation particles attached to the event (tracks)
AliEmcalJet * fLeadingJet
! leading jet (calculated event-by-event)
AliEmcalJet * GetNextAcceptJet()
void AddJetToOutputArray(AliEmcalJet *jet, Int_t arrayIndex, Int_t &jetsAlreadyInArray)
Double_t Pt() const
Definition: AliEmcalJet.h:102
Bool_t IsTrackInCone(AliVParticle *track, Double_t eta, Double_t phi, Double_t radius)
Bool_t fJetEmbeddingUsePerTrackMCPercentage
When cutting on MC percentage, calculate it per track and not for all MC tracks.
Int_t fAcceptedJets
temporary var that holds how many jets passed
TString fJetVetoArrayName
Name of array used for veto jets.
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.
Float_t GetJetRadius() const
AliEmcalList * fOutput
!output list
Int_t fJetEmbeddingNumMatchedJets
Number of matched leading jets that will be used.
AliEmcalJet * fInitialPartonMatchedJet1
! On PYTHIA data and fJetOutputMode=6, this holds the PDG code of the initial collisions that was mat...
AliEmcalJet * fInitialPartonMatchedJet2
! On PYTHIA data and fJetOutputMode=6, this holds the PDG code of the initial collisions that was mat...
Definition: External.C:220
ClassImp(AliAnalysisTaskCRC) AliAnalysisTaskCRC
TH1 * fHistEventRejection
!book keep reasons for rejecting event
void SetMakeGeneralHistograms(Bool_t g)
Base task in the EMCAL jet framework.
Bool_t IsCutFulfilled(Double_t pt, Double_t mcPt, Double_t cent, Double_t ptRatio, Double_t vetoPt)
Represent a jet reconstructed using the EMCal jet framework.
Definition: AliEmcalJet.h:44
Bool_t fJetEmbeddingUseBgrdForMCPercentage
When cutting on MC percentage, use bgrd. corr to calculate MC percentage.
Declaration of class AliEmcalPythiaInfo.
const char Option_t
Definition: External.C:48
void UserCreateOutputObjects()
Main initialization function on the worker.
virtual AliVParticle * GetNextAcceptParticle()
bool Bool_t
Definition: External.C:53
Double_t fJetEmbeddingMaxDistance
Max distance allowed to accept an embedded jet.
void FillHistogramsJets(AliEmcalJet *jet, const char *cutName)
Int_t fAcceptedJets
! number accepted jets (calculated event-by-event)
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")
Definition: External.C:196
TClonesArray * fJetVetoArray
! Array of jets imported into task used for veto a matching/embedding
Int_t fAcceptedTracks
! number accepted tracks (calculated event-by-event)
Int_t fNumRandomConesPerEvent
Number of random cones thrown in one event.
Int_t fConstPtFilterBit
For const pt plot, filter bit.