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