git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Add the eta and phi direction of the leading particle.
[u/mrichter/AliRoot.git] / PWG4 / PartCorrDep / AliAnaParticleJetLeadingConeCorrelation.cxx
CommitLineData
1c5acb87 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15/* $Id: $ */
16
17//_________________________________________________________________________
18// Class that contains the algorithm for the reconstruction of jet, cone around leading particle
19// The seed is a backward particle (direct photon)
20// 1)Take the trigger particle stored in AliAODPWG4ParticleCorrelation,
21// 2) Search for the highest pt leading particle opposite to the photon within a phi, pt window
22// 3) Take all particles around leading in a cone R with pt larger than threshold and construct the jet
23//
24// Class created from old AliPHOSGammaJet
25// (see AliRoot versions previous Release 4-09)
26//
27//*-- Author: Gustavo Conesa (LNF-INFN)
28//////////////////////////////////////////////////////////////////////////////
29
30
31// --- ROOT system ---
32#include "TH2F.h"
477d6cee 33#include "TClonesArray.h"
1c5acb87 34
35//---- Analysis system ----
36#include "AliAODTrack.h"
37#include "AliAODCaloCluster.h"
38#include "AliCaloTrackReader.h"
39#include "AliNeutralMesonSelection.h"
1c5acb87 40#include "AliAnaParticleJetLeadingConeCorrelation.h"
41#include "AliCaloPID.h"
42#include "AliAODPWG4ParticleCorrelation.h"
477d6cee 43#include "AliFidutialCut.h"
44#include "AliAODTrack.h"
45#include "AliAODCaloCluster.h"
1c5acb87 46
47ClassImp(AliAnaParticleJetLeadingConeCorrelation)
48
49
50//____________________________________________________________________________
51 AliAnaParticleJetLeadingConeCorrelation::AliAnaParticleJetLeadingConeCorrelation() :
52 AliAnaPartCorrBaseClass(), fJetsOnlyInCTS(kFALSE), fPbPb(kFALSE),
53 fSeveralConeAndPtCuts(0), fReMakeJet(0),
54 fDeltaPhiMaxCut(0.), fDeltaPhiMinCut(0.),
55 fLeadingRatioMaxCut(0.), fLeadingRatioMinCut(0.),
56 fJetCTSRatioMaxCut(0.), fJetCTSRatioMinCut(0.),
57 fJetRatioMaxCut(0.), fJetRatioMinCut(0.),
58 fJetNCone(0),fJetNPt(0), fJetCone(0),
59 fJetPtThreshold(0),fJetPtThresPbPb(0),
60 fPtTriggerSelectionCut(0.0), fSelect(0),
477d6cee 61 fSelectIsolated(0),
1c5acb87 62 //Histograms
63 fOutCont(0x0),
64 fhChargedLeadingPt(0),fhChargedLeadingPhi(0),fhChargedLeadingEta(0),
65 fhChargedLeadingDeltaPt(0),fhChargedLeadingDeltaPhi(0),fhChargedLeadingDeltaEta(0),
66 fhChargedLeadingRatioPt(0),
67 fhNeutralLeadingPt(0),fhNeutralLeadingPhi(0),fhNeutralLeadingEta(0),
68 fhNeutralLeadingDeltaPt(0),fhNeutralLeadingDeltaPhi(0),fhNeutralLeadingDeltaEta(0),
69 fhNeutralLeadingRatioPt(0),
70 fhJetPt(0),fhJetRatioPt(0),fhJetDeltaPhi(0), fhJetDeltaEta(0),
71 fhJetLeadingRatioPt(0),fhJetLeadingDeltaPhi(0),fhJetLeadingDeltaEta(0),
72 fhJetFFz(0),fhJetFFxi(0),fhJetFFpt(0),fhJetNTracksInCone(0),
73 fhBkgPt(0),fhBkgRatioPt(0),fhBkgDeltaPhi(0), fhBkgDeltaEta(0),
74 fhBkgLeadingRatioPt(0),fhBkgLeadingDeltaPhi(0),fhBkgLeadingDeltaEta(0),
75 fhBkgFFz(0),fhBkgFFxi(0),fhBkgFFpt(0),fhBkgNTracksInCone(0),
76 //Several cones and thres histograms
77 fhJetPts(),fhJetRatioPts(),fhJetDeltaPhis(), fhJetDeltaEtas(),
78 fhJetLeadingRatioPts(),fhJetLeadingDeltaPhis(),fhJetLeadingDeltaEtas(),
79 fhJetFFzs(),fhJetFFxis(),fhJetFFpts(),fhJetNTracksInCones(),
80 fhBkgPts(),fhBkgRatioPts(),fhBkgDeltaPhis(), fhBkgDeltaEtas(),
81 fhBkgLeadingRatioPts(),fhBkgLeadingDeltaPhis(),fhBkgLeadingDeltaEtas(),
82 fhBkgFFzs(),fhBkgFFxis(),fhBkgFFpts(),fhBkgNTracksInCones()
83{
84 //Default Ctor
477d6cee 85
1c5acb87 86 //Initialize parameters
87
88 for(Int_t i = 0; i<6; i++){
89 fJetXMin1[i] = 0.0 ;
90 fJetXMin2[i] = 0.0 ;
91 fJetXMax1[i] = 0.0 ;
92 fJetXMax2[i] = 0.0 ;
93 fBkgMean[i] = 0.0 ;
94 fBkgRMS[i] = 0.0 ;
95 if( i < 2 ){
96 fJetE1[i] = 0.0 ;
97 fJetE2[i] = 0.0 ;
98 fJetSigma1[i] = 0.0 ;
99 fJetSigma2[i] = 0.0 ;
100 }
101 }
477d6cee 102
1c5acb87 103 //Several cones and thres histograms
104 for(Int_t i = 0; i<5; i++){
105 fJetCones[i] = 0.0 ;
106 fJetNameCones[i] = "" ;
107 fJetPtThres[i] = 0.0 ;
108 fJetNamePtThres[i] = "" ;
109 for(Int_t j = 0; j<5; j++){
110 fhJetPts[i][j]=0 ;
111 fhJetRatioPts[i][j]=0 ;
112 fhJetDeltaPhis[i][j]=0 ;
113 fhJetDeltaEtas[i][j]=0 ;
114 fhJetLeadingRatioPts[i][j]=0 ;
115 fhJetLeadingDeltaPhis[i][j]=0 ;
116 fhJetLeadingDeltaEtas[i][j]=0 ;
117 fhJetFFzs[i][j]=0 ;
118 fhJetFFxis[i][j]=0 ;
119 fhJetFFpts[i][j]=0 ;
120 fhJetNTracksInCones[i][j]=0 ;
121 fhBkgPts[i][j]=0 ;
122 fhBkgRatioPts[i][j]=0 ;
123 fhBkgDeltaPhis[i][j]=0 ;
124 fhBkgDeltaEtas[i][j]=0 ;
125 fhBkgLeadingRatioPts[i][j]=0 ;
126 fhBkgLeadingDeltaPhis[i][j]=0 ;
127 fhBkgLeadingDeltaEtas[i][j]=0 ;
128 fhBkgFFzs[i][j]=0 ;
129 fhBkgFFxis[i][j]=0 ;
130 fhBkgFFpts[i][j]=0 ;
131 fhBkgNTracksInCones[i][j]=0 ;
132 }
133 }
134
135 InitParameters();
136
137}
138
139//____________________________________________________________________________
140AliAnaParticleJetLeadingConeCorrelation::AliAnaParticleJetLeadingConeCorrelation(const AliAnaParticleJetLeadingConeCorrelation & jetlc) :
141 AliAnaPartCorrBaseClass(jetlc), fJetsOnlyInCTS(jetlc.fJetsOnlyInCTS), fPbPb(jetlc.fPbPb),
142 fSeveralConeAndPtCuts(jetlc.fSeveralConeAndPtCuts), fReMakeJet(jetlc. fReMakeJet),
143 fDeltaPhiMaxCut(jetlc. fDeltaPhiMaxCut), fDeltaPhiMinCut(jetlc.fDeltaPhiMinCut),
144 fLeadingRatioMaxCut(jetlc.fLeadingRatioMaxCut), fLeadingRatioMinCut(jetlc.fLeadingRatioMinCut),
145 fJetCTSRatioMaxCut(jetlc.fJetCTSRatioMaxCut),
146 fJetCTSRatioMinCut(jetlc.fJetCTSRatioMinCut), fJetRatioMaxCut(jetlc.fJetRatioMaxCut),
147 fJetRatioMinCut(jetlc.fJetRatioMinCut), fJetNCone(jetlc.fJetNCone),
148 fJetNPt(jetlc.fJetNPt), fJetCone(jetlc.fJetCone),
149 fJetPtThreshold(jetlc.fJetPtThreshold),fJetPtThresPbPb(jetlc.fJetPtThresPbPb),
477d6cee 150 fPtTriggerSelectionCut(jetlc.fPtTriggerSelectionCut), fSelect(jetlc.fSelect),
151 fSelectIsolated(jetlc.fSelectIsolated),
1c5acb87 152 //Histograms
153 fOutCont(jetlc. fOutCont),
154 fhChargedLeadingPt(jetlc.fhChargedLeadingPt), fhChargedLeadingPhi(jetlc.fhChargedLeadingPhi),
155 fhChargedLeadingEta(jetlc.fhChargedLeadingEta), fhChargedLeadingDeltaPt(jetlc.fhChargedLeadingDeltaPt),
156 fhChargedLeadingDeltaPhi(jetlc.fhChargedLeadingDeltaPhi),fhChargedLeadingDeltaEta(jetlc.fhChargedLeadingDeltaEta),
157 fhChargedLeadingRatioPt(jetlc.fhChargedLeadingRatioPt),
158 fhNeutralLeadingPt(jetlc.fhNeutralLeadingPt),fhNeutralLeadingPhi(jetlc.fhNeutralLeadingPhi),
159 fhNeutralLeadingEta(jetlc.fhNeutralLeadingEta), fhNeutralLeadingDeltaPt(jetlc.fhNeutralLeadingDeltaPt),
160 fhNeutralLeadingDeltaPhi(jetlc.fhNeutralLeadingDeltaPhi),fhNeutralLeadingDeltaEta(jetlc.fhNeutralLeadingDeltaEta),
161 fhNeutralLeadingRatioPt(jetlc.fhNeutralLeadingRatioPt),
162 fhJetPt(jetlc.fhJetPt),fhJetRatioPt(jetlc.fhJetRatioPt),fhJetDeltaPhi(jetlc.fhJetDeltaPhi),
163 fhJetDeltaEta(jetlc.fhJetDeltaEta), fhJetLeadingRatioPt(jetlc.fhJetLeadingRatioPt),
164 fhJetLeadingDeltaPhi(jetlc.fhJetLeadingDeltaPhi),fhJetLeadingDeltaEta(jetlc.fhJetLeadingDeltaEta),
165 fhJetFFz(jetlc.fhJetFFz),fhJetFFxi(jetlc.fhJetFFxi),fhJetFFpt(jetlc.fhJetFFpt),
166 fhJetNTracksInCone(jetlc.fhJetNTracksInCone),
167 fhBkgPt(jetlc.fhBkgPt),fhBkgRatioPt(jetlc.fhBkgRatioPt),fhBkgDeltaPhi(jetlc.fhBkgDeltaPhi),
168 fhBkgDeltaEta(jetlc.fhBkgDeltaEta), fhBkgLeadingRatioPt(jetlc.fhBkgLeadingRatioPt),
169 fhBkgLeadingDeltaPhi(jetlc.fhBkgLeadingDeltaPhi),fhBkgLeadingDeltaEta(jetlc.fhBkgLeadingDeltaEta),
170 fhBkgFFz(jetlc.fhBkgFFz),fhBkgFFxi(jetlc.fhBkgFFxi),fhBkgFFpt(jetlc.fhBkgFFpt),
171 fhBkgNTracksInCone(jetlc.fhBkgNTracksInCone),
172 //Several cones and thres histograms
173 fhJetPts(),fhJetRatioPts(),fhJetDeltaPhis(), fhJetDeltaEtas(),
174 fhJetLeadingRatioPts(),fhJetLeadingDeltaPhis(),fhJetLeadingDeltaEtas(),
175 fhJetFFzs(),fhJetFFxis(),fhJetFFpts(),fhJetNTracksInCones(),
176 fhBkgPts(),fhBkgRatioPts(),fhBkgDeltaPhis(), fhBkgDeltaEtas(),
177 fhBkgLeadingRatioPts(),fhBkgLeadingDeltaPhis(),fhBkgLeadingDeltaEtas(),
178 fhBkgFFzs(),fhBkgFFxis(),fhBkgFFpts(),fhBkgNTracksInCones()
179{
180 // cpy ctor
181
182 for(Int_t i = 0; i<6; i++){
183 fJetXMin1[i] = jetlc.fJetXMin1[i] ;
184 fJetXMin2[i] = jetlc.fJetXMin2[i] ;
185 fJetXMax1[i] = jetlc.fJetXMax1[i] ;
186 fJetXMax2[i] = jetlc.fJetXMax2[i] ;
187 fBkgMean[i] = jetlc.fBkgMean[i] ;
188 fBkgRMS[i] = jetlc.fBkgRMS[i] ;
189 if( i < 2 ){
190 fJetE1[i] = jetlc.fJetE1[i] ;
191 fJetE2[i] = jetlc.fJetE2[i] ;
192 fJetSigma1[i] = jetlc.fJetSigma1[i] ;
193 fJetSigma2[i] = jetlc.fJetSigma2[i] ;
194 }
195 }
477d6cee 196
1c5acb87 197 //Several cones and thres histograms
198 for(Int_t i = 0; i<5; i++){
199 fJetCones[i] = jetlc.fJetCones[i] ;
200 fJetNameCones[i] = jetlc.fJetNameCones[i] ;
201 fJetPtThres[i] = jetlc.fJetPtThres[i] ;
202 fJetNamePtThres[i] = jetlc.fJetNamePtThres[i] ;
203 for(Int_t j = 0; j<5; j++){
204 fhJetPts[i][j] = jetlc.fhJetPts[i][j] ;
205 fhJetRatioPts[i][j] = jetlc.fhJetRatioPts[i][j] ;
206 fhJetDeltaPhis[i][j] = jetlc.fhJetDeltaPhis[i][j] ;
207 fhJetDeltaEtas[i][j] = jetlc.fhJetDeltaEtas[i][j] ;
208 fhJetLeadingRatioPts[i][j] = jetlc.fhJetLeadingRatioPts[i][j] ;
209 fhJetLeadingDeltaPhis[i][j] = jetlc.fhJetLeadingDeltaPhis[i][j] ;
210 fhJetLeadingDeltaEtas[i][j] = jetlc.fhJetLeadingDeltaEtas[i][j] ;
211 fhJetFFzs[i][j] = jetlc.fhJetFFzs[i][j] ;
212 fhJetFFxis[i][j] = jetlc.fhJetFFxis[i][j] ;
213 fhJetFFpts[i][j] = jetlc.fhJetFFpts[i][j] ;
214 fhJetNTracksInCones[i][j] = fhJetNTracksInCones[i][j] ;
215 fhBkgPts[i][j] = jetlc.fhBkgPts[i][j] ;
216 fhBkgRatioPts[i][j] = jetlc.fhBkgRatioPts[i][j] ;
217 fhBkgDeltaPhis[i][j] = jetlc.fhBkgDeltaPhis[i][j] ;
218 fhBkgDeltaEtas[i][j] = jetlc.fhBkgDeltaEtas[i][j] ;
219 fhBkgLeadingRatioPts[i][j] = jetlc.fhBkgLeadingRatioPts[i][j] ;
220 fhBkgLeadingDeltaPhis[i][j] = jetlc.fhBkgLeadingDeltaPhis[i][j] ;
221 fhBkgLeadingDeltaEtas[i][j] = jetlc.fhBkgLeadingDeltaEtas[i][j] ;
222 fhBkgFFzs[i][j] = jetlc.fhBkgFFzs[i][j] ;
223 fhBkgFFxis[i][j] = jetlc.fhBkgFFxis[i][j] ;
224 fhBkgFFpts[i][j] = jetlc.fhBkgFFpts[i][j] ;
225 fhBkgNTracksInCones[i][j] = jetlc.fhBkgNTracksInCones[i][j] ;
226 }
227 }
228}
229
230//_________________________________________________________________________
231AliAnaParticleJetLeadingConeCorrelation & AliAnaParticleJetLeadingConeCorrelation::operator = (const AliAnaParticleJetLeadingConeCorrelation & jetlc)
232{
233 // assignment operator
477d6cee 234
1c5acb87 235 if(this == &jetlc)return *this;
236 ((AliAnaPartCorrBaseClass *)this)->operator=(jetlc);
1c5acb87 237
477d6cee 238 fSeveralConeAndPtCuts = jetlc.fSeveralConeAndPtCuts ;
239 fPbPb = jetlc.fPbPb ;
240 fReMakeJet = jetlc.fReMakeJet ;
241 fJetsOnlyInCTS = jetlc.fJetsOnlyInCTS;
242
243 fDeltaPhiMaxCut = jetlc.fDeltaPhiMaxCut ;
244 fDeltaPhiMinCut = jetlc.fDeltaPhiMinCut ;
245 fLeadingRatioMaxCut = jetlc.fLeadingRatioMaxCut ;
246 fLeadingRatioMinCut = jetlc.fLeadingRatioMinCut ;
247
248 fJetCTSRatioMaxCut = jetlc.fJetCTSRatioMaxCut ;
249 fJetCTSRatioMinCut = jetlc.fJetCTSRatioMinCut ;
250 fJetRatioMaxCut = jetlc.fJetRatioMaxCut ;
251 fJetRatioMinCut = jetlc.fJetRatioMinCut ;
1c5acb87 252
477d6cee 253 fJetNCone = jetlc.fJetNCone ;
254 fJetNPt = jetlc.fJetNPt ; fJetCone = jetlc.fJetCone ;
255 fJetPtThreshold = jetlc.fJetPtThreshold ;
256 fJetPtThresPbPb = jetlc.fJetPtThresPbPb ;
1c5acb87 257 fPtTriggerSelectionCut = jetlc.fPtTriggerSelectionCut ;
477d6cee 258 fSelect = jetlc.fSelect ;
259 fSelectIsolated = jetlc.fSelectIsolated ;
1c5acb87 260
261 for(Int_t i = 0; i<6; i++){
262 fJetXMin1[i] = jetlc.fJetXMin1[i] ;
263 fJetXMin2[i] = jetlc.fJetXMin2[i] ;
264 fJetXMax1[i] = jetlc.fJetXMax1[i] ;
265 fJetXMax2[i] = jetlc.fJetXMax2[i] ;
266 fBkgMean[i] = jetlc.fBkgMean[i] ;
267 fBkgRMS[i] = jetlc.fBkgRMS[i] ;
268 if( i < 2 ){
269 fJetE1[i] = jetlc.fJetE1[i] ;
270 fJetE2[i] = jetlc.fJetE2[i] ;
271 fJetSigma1[i] = jetlc.fJetSigma1[i] ;
272 fJetSigma2[i] = jetlc.fJetSigma2[i] ;
273 }
274 }
275
276 //Histograms
277 fOutCont = jetlc. fOutCont ;
278 fhChargedLeadingPt = jetlc.fhChargedLeadingPt; fhChargedLeadingPhi = jetlc.fhChargedLeadingPhi;
279 fhChargedLeadingEta = jetlc.fhChargedLeadingEta; fhChargedLeadingDeltaPt = jetlc.fhChargedLeadingDeltaPt;
280 fhChargedLeadingDeltaPhi = jetlc.fhChargedLeadingDeltaPhi;fhChargedLeadingDeltaEta = jetlc.fhChargedLeadingDeltaEta;
281 fhChargedLeadingRatioPt = jetlc.fhChargedLeadingRatioPt;
282 fhNeutralLeadingPt = jetlc.fhNeutralLeadingPt;fhNeutralLeadingPhi = jetlc.fhNeutralLeadingPhi;
283 fhNeutralLeadingEta = jetlc.fhNeutralLeadingEta; fhNeutralLeadingDeltaPt = jetlc.fhNeutralLeadingDeltaPt;
284 fhNeutralLeadingDeltaPhi = jetlc.fhNeutralLeadingDeltaPhi;fhNeutralLeadingDeltaEta = jetlc.fhNeutralLeadingDeltaEta;
285 fhNeutralLeadingRatioPt = jetlc.fhNeutralLeadingRatioPt;
286 fhJetPt = jetlc.fhJetPt;fhJetRatioPt = jetlc.fhJetRatioPt;fhJetDeltaPhi = jetlc.fhJetDeltaPhi;
287 fhJetDeltaEta = jetlc.fhJetDeltaEta; fhJetLeadingRatioPt = jetlc.fhJetLeadingRatioPt;
288 fhJetLeadingDeltaPhi = jetlc.fhJetLeadingDeltaPhi;fhJetLeadingDeltaEta = jetlc.fhJetLeadingDeltaEta;
289 fhJetFFz = jetlc.fhJetFFz;fhJetFFxi = jetlc.fhJetFFxi;fhJetFFpt = jetlc.fhJetFFpt;
290 fhJetNTracksInCone = jetlc.fhJetNTracksInCone;
291 fhBkgPt = jetlc.fhBkgPt;fhBkgRatioPt = jetlc.fhBkgRatioPt;fhBkgDeltaPhi = jetlc.fhBkgDeltaPhi;
292 fhBkgDeltaEta = jetlc.fhBkgDeltaEta; fhBkgLeadingRatioPt = jetlc.fhBkgLeadingRatioPt;
293 fhBkgLeadingDeltaPhi = jetlc.fhBkgLeadingDeltaPhi;fhBkgLeadingDeltaEta = jetlc.fhBkgLeadingDeltaEta;
294 fhBkgFFz = jetlc.fhBkgFFz;fhBkgFFxi = jetlc.fhBkgFFxi;fhBkgFFpt = jetlc.fhBkgFFpt;
295 fhBkgNTracksInCone = jetlc.fhBkgNTracksInCone;
296
297
298 //Several cones and thres histograms
299 for(Int_t i = 0; i<5; i++){
300 fJetCones[i] = jetlc.fJetCones[i] ;
301 fJetNameCones[i] = jetlc.fJetNameCones[i] ;
302 fJetPtThres[i] = jetlc.fJetPtThres[i] ;
303 fJetNamePtThres[i] = jetlc.fJetNamePtThres[i] ;
304
305 for(Int_t j = 0; j<5; j++){
306 fhJetPts[i][j] = jetlc.fhJetPts[i][j] ;
307 fhJetRatioPts[i][j] = jetlc.fhJetRatioPts[i][j] ;
308 fhJetDeltaPhis[i][j] = jetlc.fhJetDeltaPhis[i][j] ;
309 fhJetDeltaEtas[i][j] = jetlc.fhJetDeltaEtas[i][j] ;
310 fhJetLeadingRatioPts[i][j] = jetlc.fhJetLeadingRatioPts[i][j] ;
311 fhJetLeadingDeltaPhis[i][j] = jetlc.fhJetLeadingDeltaPhis[i][j] ;
312 fhJetLeadingDeltaEtas[i][j] = jetlc.fhJetLeadingDeltaEtas[i][j] ;
313 fhJetFFzs[i][j] = jetlc.fhJetFFzs[i][j] ;
314 fhJetFFxis[i][j] = jetlc.fhJetFFxis[i][j] ;
315 fhJetFFpts[i][j] = jetlc.fhJetFFpts[i][j] ;
316 fhJetNTracksInCones[i][j] = fhJetNTracksInCones[i][j] ;
317 fhBkgPts[i][j] = jetlc.fhBkgPts[i][j] ;
318 fhBkgRatioPts[i][j] = jetlc.fhBkgRatioPts[i][j] ;
319 fhBkgDeltaPhis[i][j] = jetlc.fhBkgDeltaPhis[i][j] ;
320 fhBkgDeltaEtas[i][j] = jetlc.fhBkgDeltaEtas[i][j] ;
321 fhBkgLeadingRatioPts[i][j] = jetlc.fhBkgLeadingRatioPts[i][j] ;
322 fhBkgLeadingDeltaPhis[i][j] = jetlc.fhBkgLeadingDeltaPhis[i][j] ;
323 fhBkgLeadingDeltaEtas[i][j] = jetlc.fhBkgLeadingDeltaEtas[i][j] ;
324 fhBkgFFzs[i][j] = jetlc.fhBkgFFzs[i][j] ;
325 fhBkgFFxis[i][j] = jetlc.fhBkgFFxis[i][j] ;
326 fhBkgFFpts[i][j] = jetlc.fhBkgFFpts[i][j] ;
327 fhBkgNTracksInCones[i][j] = jetlc.fhBkgNTracksInCones[i][j] ;
328 }
329 }
330
331 return *this;
332
333}
334
335//____________________________________________________________________________
336AliAnaParticleJetLeadingConeCorrelation::~AliAnaParticleJetLeadingConeCorrelation()
337{
338 // Remove all pointers except analysis output pointers.
339 delete [] fJetE1;
340 delete [] fJetE2;
341 delete [] fJetSigma1;
342 delete [] fJetSigma2;
343 delete [] fBkgMean;
344 delete [] fBkgRMS;
345 delete [] fJetXMin1;
346 delete [] fJetXMin2;
347 delete [] fJetXMax1;
348 delete [] fJetXMax2;
349 delete [] fJetCones;
350 delete [] fJetNameCones;
351 delete [] fJetPtThres;
352 delete [] fJetNamePtThres;
353}
354
355//____________________________________________________________________________
356Double_t AliAnaParticleJetLeadingConeCorrelation::CalculateJetRatioLimit(const Double_t ptg, const Double_t *par, const Double_t *x) const {
357 //Calculate the ratio of the jet and trigger particle limit for the selection
358 //WARNING: need to check what it does
a3aebfff 359 //printf("CalculateLimit: x1 %f, x2%f\n",x[0],x[1]);
1c5acb87 360 Double_t ePP = par[0] + par[1] * ptg ;
361 Double_t sPP = par[2] + par[3] * ptg ;
362 Double_t f = x[0] + x[1] * ptg ;
363 Double_t ePbPb = ePP + par[4] ;
364 Double_t sPbPb = TMath::Sqrt(sPP*sPP+ par[5]*par[5]) ;
365 Double_t rat = (ePbPb - sPbPb * f) / ptg ;
a3aebfff 366 //printf("CalculateLimit: ePP %f, sPP %f, f %f\n", ePP, sPP, f);
367 //printf("CalculateLimit: ePbPb %f, sPbPb %f, rat %f\n", ePbPb, sPbPb, rat);
1c5acb87 368 return rat ;
369}
370
371//____________________________________________________________________________
372void AliAnaParticleJetLeadingConeCorrelation::FillJetHistos(AliAODPWG4ParticleCorrelation * particle, const TLorentzVector leading, const TLorentzVector jet, const TString type, const TString lastname)
373{
374 //Fill jet and background histograms
375 Double_t ptTrig = particle->Pt();
376 Double_t ptJet = jet.Pt();
377 Double_t ptLead = leading.Pt();
378 Double_t phiTrig = particle->Phi();
379 Double_t phiJet = jet.Phi();
380 Double_t phiLead = leading.Phi();
381 Double_t etaTrig = particle->Eta();
382 Double_t etaJet = jet.Eta();
383 Double_t etaLead = leading.Eta();
384
1e1befe8 385 dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"Pt"+lastname))->
1c5acb87 386 Fill(ptTrig,ptJet);
1e1befe8 387
388 dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"RatioPt"+lastname))->
1c5acb87 389 Fill(ptTrig,ptJet/ptTrig);
1e1befe8 390 dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"LeadingRatioPt"+lastname))->
1c5acb87 391 Fill(ptTrig,ptLead/ptJet);
1e1befe8 392// dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"Phi"+lastname))->
1c5acb87 393// Fill(ptTrig,phiJet);
1e1befe8 394 dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"DeltaPhi"+lastname))->
1c5acb87 395 Fill(ptTrig,phiJet-phiTrig);
1e1befe8 396 dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"LeadingDeltaPhi"+lastname))->
1c5acb87 397 Fill(ptTrig,phiJet-phiLead);
398
1e1befe8 399 // dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"Eta"+lastname))->
1c5acb87 400 // Fill(ptTrig,etaJet);
1e1befe8 401 dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"DeltaEta"+lastname))->
1c5acb87 402 Fill(ptTrig,etaJet-etaTrig);
1e1befe8 403 dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"LeadingDeltaEta"+lastname))->
1c5acb87 404 Fill(ptTrig,etaJet-etaLead);
405
406 //Construct fragmentation function
407 TRefArray * pl = new TRefArray;
a3aebfff 408 if(type == "Jet") pl = particle->GetRefArray(GetAODRefArrayName()+"Tracks");
409 else if(type == "Bkg") particle->GetRefArray(GetAODRefArrayName()+"TracksBkg");
1c5acb87 410
411 //Different pt cut for jet particles in different collisions systems
412 //Only needed when jet is recalculated from AODs
413 Float_t ptcut = fJetPtThreshold;
414 if(fPbPb && !fSeveralConeAndPtCuts && ptTrig > fPtTriggerSelectionCut) ptcut = fJetPtThresPbPb ;
415
416 TVector3 p3;
417 Int_t nTracksInCone = 0;
418 for(Int_t ipr = 0;ipr < pl->GetEntriesFast() ; ipr ++ ){
419 AliAODTrack* track = dynamic_cast<AliAODTrack *>(pl->At(ipr)) ;
420 p3.SetXYZ(track->Px(),track->Py(),track->Pz());
421
422 //Recheck if particle is in jet cone
423 if(fReMakeJet || fSeveralConeAndPtCuts)
424 if(!IsParticleInJetCone(p3.Eta(), p3.Phi(), leading.Eta(), leading.Phi()) ) continue ;
425
426 nTracksInCone++;
427
1e1befe8 428 dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"FFz"+lastname))
1c5acb87 429 ->Fill(ptTrig,p3.Pt()/ptTrig);
1e1befe8 430 dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"FFxi"+lastname))
1c5acb87 431 ->Fill(ptTrig,TMath::Log(ptTrig/p3.Pt()));
1e1befe8 432 dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"FFpt"+lastname))
1c5acb87 433 ->Fill(ptTrig,p3.Pt());
434
435 }//track loop
436
1e1befe8 437 if(nTracksInCone > 0) dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"NTracksInCone"+lastname))
1c5acb87 438 ->Fill(ptTrig, nTracksInCone);
439
440}
441
442//________________________________________________________________________
443TList * AliAnaParticleJetLeadingConeCorrelation::GetCreateOutputObjects()
444{
445 // Create histograms to be saved in output file and
446 // store them in fOutCont
447
477d6cee 448 if(GetDebug()>1) printf("AliAnaParticleJetLeadingConeCorrelation::GetCreateOutputObjects() - Init histograms \n");
449
450 fOutCont = new TList() ;
451 fOutCont->SetName("ParticleJetLeadingInConeCorrelationHistograms") ;
452
453 Int_t nptbins = GetHistoNPtBins();
454 Int_t nphibins = GetHistoNPhiBins();
455 Int_t netabins = GetHistoNEtaBins();
456 Float_t ptmax = GetHistoPtMax();
457 Float_t phimax = GetHistoPhiMax();
458 Float_t etamax = GetHistoEtaMax();
459 Float_t ptmin = GetHistoPtMin();
460 Float_t phimin = GetHistoPhiMin();
461 Float_t etamin = GetHistoEtaMin();
462
463 fhChargedLeadingPt = new TH2F("ChargedLeadingPt","p_{T leading charge} vs p_{T trigger}",nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
464 fhChargedLeadingPt->SetYTitle("p_{T leading charge} /p_{T trigger}");
465 fhChargedLeadingPt->SetXTitle("p_{T trigger} (GeV/c)");
466
467 fhChargedLeadingPhi = new TH2F("ChargedLeadingPhi","#phi_{h^{#pm}} vs p_{T trigger}", nptbins,ptmin,ptmax,nphibins,phimin,phimax);
468 fhChargedLeadingPhi->SetYTitle("#phi_{h^{#pm}} (rad)");
469 fhChargedLeadingPhi->SetXTitle("p_{T trigger} (GeV/c)");
470
471 fhChargedLeadingEta = new TH2F("ChargedLeadingEta","#eta_{h^{#pm}} vs p_{T trigger}",nptbins,ptmin,ptmax,netabins,etamin,etamax);
472 fhChargedLeadingEta->SetYTitle("#eta_{h^{#pm}} ");
473 fhChargedLeadingEta->SetXTitle("p_{T trigger} (GeV/c)");
474
475 fhChargedLeadingDeltaPt = new TH2F("ChargedLeadingDeltaPt","#p_{T trigger} - #p_{T h^{#pm}} vs p_{T trigger}",nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
476 fhChargedLeadingDeltaPt->SetYTitle("#Delta p_{T} (GeV/c)");
477 fhChargedLeadingDeltaPt->SetXTitle("p_{T trigger} (GeV/c)");
478
479 fhChargedLeadingDeltaPhi = new TH2F("ChargedLeadingDeltaPhi","#phi_{trigger} - #phi_{h^{#pm}} vs p_{T trigger}",nptbins,ptmin,ptmax,120,0,TMath::TwoPi());
480 fhChargedLeadingDeltaPhi->SetYTitle("#Delta #phi (rad)");
481 fhChargedLeadingDeltaPhi->SetXTitle("p_{T trigger} (GeV/c)");
482
483 fhChargedLeadingDeltaEta = new TH2F("ChargedLeadingDeltaEta","#eta_{trigger} - #eta_{h^{#pm}} vs p_{T trigger}",nptbins,ptmin,ptmax,120,-2,2);
484 fhChargedLeadingDeltaEta->SetYTitle("#Delta #eta");
485 fhChargedLeadingDeltaEta->SetXTitle("p_{T trigger} (GeV/c)");
486
487 fhChargedLeadingRatioPt = new TH2F("ChargedLeadingRatioPt","p_{T leading charge} /p_{T trigger} vs p_{T trigger}",nptbins,ptmin,ptmax,120,0,2);
488 fhChargedLeadingRatioPt->SetYTitle("p_{T lead charge} /p_{T trigger}");
489 fhChargedLeadingRatioPt->SetXTitle("p_{T trigger} (GeV/c)");
490
491 fOutCont->Add(fhChargedLeadingPt) ;
492 fOutCont->Add(fhChargedLeadingPhi) ;
493 fOutCont->Add(fhChargedLeadingEta) ;
494 fOutCont->Add(fhChargedLeadingDeltaPt) ;
495 fOutCont->Add(fhChargedLeadingDeltaPhi) ;
496 fOutCont->Add(fhChargedLeadingDeltaEta) ;
497 fOutCont->Add(fhChargedLeadingRatioPt) ;
498
499 if(!fJetsOnlyInCTS){
500
501 fhNeutralLeadingPt = new TH2F("NeutralLeadingPt","p_{T leading #pi^{0}} vs p_{T trigger}",nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
502 fhNeutralLeadingPt->SetYTitle("p_{T leading #pi^{0}} /p_{T trigger}");
503 fhNeutralLeadingPt->SetXTitle("p_{T trigger} (GeV/c)");
504
505 fhNeutralLeadingPhi = new TH2F("NeutralLeadingPhi","#phi_{#pi^{0}} vs p_{T trigger}",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
506 fhNeutralLeadingPhi->SetYTitle("#phi_{#pi^{0}} (rad)");
507 fhNeutralLeadingPhi->SetXTitle("p_{T trigger} (GeV/c)");
508
509 fhNeutralLeadingEta = new TH2F("NeutralLeadingEta","#eta_{#pi^{0}} vs p_{T trigger}",nptbins,ptmin,ptmax,netabins,etamin,etamax);
510 fhNeutralLeadingEta->SetYTitle("#eta_{#pi^{0}} ");
511 fhNeutralLeadingEta->SetXTitle("p_{T trigger} (GeV/c)");
512
513 fhNeutralLeadingDeltaPt = new TH2F("NeutralLeadingDeltaPt","#p_{T trigger} - #p_{T #pi^{0}} vs p_{T trigger}",nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
514 fhNeutralLeadingDeltaPt->SetYTitle("#Delta p_{T} (GeV/c)");
515 fhNeutralLeadingDeltaPt->SetXTitle("p_{T trigger} (GeV/c)");
516
517 fhNeutralLeadingDeltaPhi = new TH2F("NeutralLeadingDeltaPhi","#phi_{trigger} - #phi_{#pi^{0}} vs p_{T trigger}",nptbins,ptmin,ptmax,120,0,TMath::TwoPi());
518 fhNeutralLeadingDeltaPhi->SetYTitle("#Delta #phi (rad)");
519 fhNeutralLeadingDeltaPhi->SetXTitle("p_{T trigger} (GeV/c)");
520
521 fhNeutralLeadingDeltaEta = new TH2F("NeutralLeadingDeltaEta","#eta_{trigger} - #eta_{#pi^{0}} vs p_{T trigger}",nptbins,ptmin,ptmax,120,-2,2);
522 fhNeutralLeadingDeltaEta->SetYTitle("#Delta #eta");
523 fhNeutralLeadingDeltaEta->SetXTitle("p_{T trigger} (GeV/c)");
524
525 fhNeutralLeadingRatioPt = new TH2F("NeutralLeadingRatioPt","p_{T leading #pi^{0}} /p_{T trigger} vs p_{T trigger}",nptbins,ptmin,ptmax,120,0,2);
526 fhNeutralLeadingRatioPt->SetYTitle("p_{T lead #pi^{0}} /p_{T trigger}");
527 fhNeutralLeadingRatioPt->SetXTitle("p_{T trigger} (GeV/c)");
528
529 fOutCont->Add(fhNeutralLeadingPt) ;
530 fOutCont->Add(fhNeutralLeadingPhi) ;
531 fOutCont->Add(fhNeutralLeadingEta) ;
532 fOutCont->Add(fhNeutralLeadingDeltaPt) ;
533 fOutCont->Add(fhNeutralLeadingDeltaPhi) ;
534 fOutCont->Add(fhNeutralLeadingDeltaEta) ;
535 fOutCont->Add(fhNeutralLeadingRatioPt) ;
536
537 }
538
539 if(!fSeveralConeAndPtCuts){// not several cones
540
541 //Jet Distributions
542 fhJetPt = new TH2F("JetPt","p_{T jet} vs p_{T trigger}",nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
543 fhJetPt->SetYTitle("p_{T jet}");
544 fhJetPt->SetXTitle("p_{T trigger} (GeV/c)");
545
546 fhJetRatioPt = new TH2F("JetRatioPt","p_{T jet}/p_{T trigger} vs p_{T trigger}",nptbins,ptmin,ptmax,120,0,2);
547 fhJetRatioPt->SetYTitle("p_{T jet}/p_{T trigger}");
548 fhJetRatioPt->SetXTitle("p_{T trigger} (GeV/c)");
549
550 fhJetDeltaPhi = new TH2F("JetDeltaPhi","#phi_{jet} - #phi_{trigger} vs p_{T trigger}",nptbins,ptmin,ptmax,120,0,TMath::TwoPi());
551 fhJetDeltaPhi->SetYTitle("#Delta #phi (rad)");
552 fhJetDeltaPhi->SetXTitle("p_{T trigger} (GeV/c)");
553
554 fhJetDeltaEta = new TH2F("JetDeltaEta","#eta_{jet} - #eta_{trigger} vs p_{T trigger}",nptbins,ptmin,ptmax,120,-2,2);
555 fhJetDeltaEta->SetYTitle("#Delta #eta");
556 fhJetDeltaEta->SetXTitle("p_{T trigger} (GeV/c)");
557
558 fhJetLeadingRatioPt = new TH2F("JetLeadingRatioPt","p_{T jet} vs p_{T trigger}",nptbins,ptmin,ptmax,120,0,2);
559 fhJetLeadingRatioPt->SetYTitle("p_{T leading}/p_{T jet}");
560 fhJetLeadingRatioPt->SetXTitle("p_{T trigger} (GeV/c)");
561
562 fhJetLeadingDeltaPhi = new TH2F("JetLeadingDeltaPhi","#phi_{jet} - #phi_{leading} vs p_{T trigger}",nptbins,ptmin,ptmax,120,0,TMath::TwoPi());
563 fhJetLeadingDeltaPhi->SetYTitle("#Delta #phi (rad)");
564 fhJetLeadingDeltaPhi->SetXTitle("p_{T trigger} (GeV/c)");
565
566 fhJetLeadingDeltaEta = new TH2F("JetLeadingDeltaEta","#eta_{jet} - #eta_{leading} vs p_{T trigger}",nptbins,ptmin,ptmax,120,-2,2);
567 fhJetLeadingDeltaEta->SetYTitle("#Delta #eta");
568 fhJetLeadingDeltaEta->SetXTitle("p_{T trigger} (GeV/c)");
569
570 fhJetFFz = new TH2F("JetFFz","z = p_{T i charged}/p_{T trigger} vs p_{T trigger}",nptbins,ptmin,ptmax,200,0.,2);
571 fhJetFFz->SetYTitle("z");
572 fhJetFFz->SetXTitle("p_{T trigger}");
573
574 fhJetFFxi = new TH2F("JetFFxi","#xi = ln(p_{T trigger}/p_{T i charged}) vs p_{T trigger}",nptbins,ptmin,ptmax,100,0.,10.);
575 fhJetFFxi->SetYTitle("#xi");
576 fhJetFFxi->SetXTitle("p_{T trigger}");
577
578 fhJetFFpt = new TH2F("JetFFpt","#xi = p_{T i charged}) vs p_{T trigger}",nptbins,ptmin,ptmax,200,0.,50.);
579 fhJetFFpt->SetYTitle("p_{T charged hadron}");
580 fhJetFFpt->SetXTitle("p_{T trigger}");
581
582 fhJetNTracksInCone = new TH2F("JetNTracksInCone","N particles in cone vs p_{T trigger}",nptbins,ptmin,ptmax,5000,0, 5000);
583 fhJetNTracksInCone->SetYTitle("N tracks in jet cone");
584 fhJetNTracksInCone->SetXTitle("p_{T trigger} (GeV/c)");
585
586 fOutCont->Add(fhJetPt) ;
587 fOutCont->Add(fhJetRatioPt) ;
588 fOutCont->Add(fhJetDeltaPhi) ;
589 fOutCont->Add(fhJetDeltaEta) ;
590 fOutCont->Add(fhJetLeadingRatioPt) ;
591 fOutCont->Add(fhJetLeadingDeltaPhi) ;
592 fOutCont->Add(fhJetLeadingDeltaEta) ;
593 fOutCont->Add(fhJetFFz) ;
594 fOutCont->Add(fhJetFFxi) ;
595 fOutCont->Add(fhJetFFpt) ;
596 fOutCont->Add(fhJetNTracksInCone) ;
597
598 //Bkg Distributions
599 fhBkgPt = new TH2F("BkgPt","p_{T bkg} vs p_{T trigger}",nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
600 fhBkgPt->SetYTitle("p_{T bkg}");
601 fhBkgPt->SetXTitle("p_{T trigger} (GeV/c)");
602
603 fhBkgRatioPt = new TH2F("BkgRatioPt","p_{T bkg}/p_{T trigger} vs p_{T trigger}",nptbins,ptmin,ptmax,120,0,2);
604 fhBkgRatioPt->SetYTitle("p_{T bkg}/p_{T trigger}");
605 fhBkgRatioPt->SetXTitle("p_{T trigger} (GeV/c)");
606
607 fhBkgDeltaPhi = new TH2F("BkgDeltaPhi","#phi_{bkg} - #phi_{trigger} vs p_{T trigger}",nptbins,ptmin,ptmax,120,0,TMath::TwoPi());
608 fhBkgDeltaPhi->SetYTitle("#Delta #phi (rad)");
609 fhBkgDeltaPhi->SetXTitle("p_{T trigger} (GeV/c)");
610
611 fhBkgDeltaEta = new TH2F("BkgDeltaEta","#eta_{bkg} - #eta_{trigger} vs p_{T trigger}",nptbins,ptmin,ptmax,120,-2,2);
612 fhBkgDeltaEta->SetYTitle("#Delta #eta");
613 fhBkgDeltaEta->SetXTitle("p_{T trigger} (GeV/c)");
614
615 fhBkgLeadingRatioPt = new TH2F("BkgLeadingRatioPt","p_{T bkg} vs p_{T trigger}",nptbins,ptmin,ptmax,120,0,2);
616 fhBkgLeadingRatioPt->SetYTitle("p_{T leading}/p_{T bkg}");
617 fhBkgLeadingRatioPt->SetXTitle("p_{T trigger} (GeV/c)");
618
619 fhBkgLeadingDeltaPhi = new TH2F("BkgLeadingDeltaPhi","#phi_{bkg} - #phi_{leading} vs p_{T trigger}",nptbins,ptmin,ptmax,120,0,TMath::TwoPi());
620 fhBkgLeadingDeltaPhi->SetYTitle("#Delta #phi (rad)");
621 fhBkgLeadingDeltaPhi->SetXTitle("p_{T trigger} (GeV/c)");
622
623 fhBkgLeadingDeltaEta = new TH2F("BkgLeadingDeltaEta","#eta_{bkg} - #eta_{leading} vs p_{T trigger}",nptbins,ptmin,ptmax,120,-2,2);
624 fhBkgLeadingDeltaEta->SetYTitle("#Delta #eta");
625 fhBkgLeadingDeltaEta->SetXTitle("p_{T trigger} (GeV/c)");
626
627 fhBkgFFz = new TH2F("BkgFFz","z = p_{T i charged}/p_{T trigger} vs p_{T trigger}", nptbins,ptmin,ptmax,200,0.,2);
628 fhBkgFFz->SetYTitle("z");
629 fhBkgFFz->SetXTitle("p_{T trigger}");
630
631 fhBkgFFxi = new TH2F("BkgFFxi","#xi = ln(p_{T trigger}/p_{T i charged}) vs p_{T trigger}", nptbins,ptmin,ptmax,100,0.,10.);
632 fhBkgFFxi->SetYTitle("#xi");
633 fhBkgFFxi->SetXTitle("p_{T trigger}");
634
635 fhBkgFFpt = new TH2F("BkgFFpt","p_{T charged hadron } vs p_{T trigger}", nptbins,ptmin,ptmax,200,0.,50.);
636 fhBkgFFpt->SetYTitle("p_{T charged} hadron");
637 fhBkgFFpt->SetXTitle("p_{T trigger}");
638
639 fhBkgNTracksInCone = new TH2F("BkgNTracksInCone","N particles in cone vs p_{T trigger}",nptbins,ptmin,ptmax,5000,0, 5000);
640 fhBkgNTracksInCone->SetYTitle("N tracks in bkg cone");
641 fhBkgNTracksInCone->SetXTitle("p_{T trigger} (GeV/c)");
642
643 fOutCont->Add(fhBkgPt) ;
644 fOutCont->Add(fhBkgRatioPt) ;
645 fOutCont->Add(fhBkgDeltaPhi) ;
646 fOutCont->Add(fhBkgDeltaEta) ;
647 fOutCont->Add(fhBkgLeadingRatioPt) ;
648 fOutCont->Add(fhBkgLeadingDeltaPhi) ;
649 fOutCont->Add(fhBkgLeadingDeltaEta) ;
650 fOutCont->Add(fhBkgFFz) ;
651 fOutCont->Add(fhBkgFFxi) ;
652 fOutCont->Add(fhBkgFFpt) ;
653 fOutCont->Add(fhBkgNTracksInCone) ;
654
655 }//not several cones
656 else{ //If we want to study the jet for different cones and pt
657 for(Int_t icone = 0; icone<fJetNCone; icone++){//icone
658 for(Int_t ipt = 0; ipt<fJetNPt;ipt++){ //ipt
1c5acb87 659
477d6cee 660 TString lastnamehist ="Cone"+ fJetNameCones[icone]+"Pt"+ fJetNamePtThres[ipt];
661 TString lastnametitle =", cone ="+fJetNameCones[icone]+", pt > " +fJetNamePtThres[ipt]+" GeV/c";
1c5acb87 662
477d6cee 663 //Jet Distributions
664 fhJetPts[icone][ipt] = new TH2F("JetPt"+lastnamehist,"p_{T jet} vs p_{T trigger}"+lastnametitle,nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
665 fhJetPts[icone][ipt]->SetYTitle("p_{T jet}");
666 fhJetPts[icone][ipt]->SetXTitle("p_{T trigger} (GeV/c)");
1c5acb87 667
477d6cee 668 fhJetRatioPts[icone][ipt] = new TH2F("JetRatioPt"+lastnamehist,"p_{T jet}/p_{T trigger} vs p_{T trigger}"+lastnametitle,nptbins,ptmin,ptmax,120,0,2);
669 fhJetRatioPts[icone][ipt]->SetYTitle("p_{T jet}/p_{T trigger}");
670 fhJetRatioPts[icone][ipt]->SetXTitle("p_{T trigger} (GeV/c)");
1c5acb87 671
477d6cee 672 fhJetDeltaPhis[icone][ipt] = new TH2F("JetDeltaPhi"+lastnamehist,"#phi_{jet} - #phi_{trigger} vs p_{T trigger}"+lastnametitle,nptbins,ptmin,ptmax,120,0,TMath::TwoPi());
673 fhJetDeltaPhis[icone][ipt]->SetYTitle("#Delta #phi (rad)");
674 fhJetDeltaPhis[icone][ipt]->SetXTitle("p_{T trigger} (GeV/c)");
1c5acb87 675
477d6cee 676 fhJetDeltaEtas[icone][ipt] = new TH2F("JetDeltaEta"+lastnamehist,"#eta_{jet} - #eta_{trigger} vs p_{T trigger}"+lastnametitle,nptbins,ptmin,ptmax,120,-2,2);
677 fhJetDeltaEtas[icone][ipt]->SetYTitle("#Delta #eta");
678 fhJetDeltaEtas[icone][ipt]->SetXTitle("p_{T trigger} (GeV/c)");
1c5acb87 679
477d6cee 680 fhJetLeadingRatioPts[icone][ipt] = new TH2F("JetLeadingRatioPt"+lastnamehist,"p_{T jet} vs p_{T trigger}"+lastnametitle,nptbins,ptmin,ptmax,120,0,2);
681 fhJetLeadingRatioPts[icone][ipt]->SetYTitle("p_{T leading}/p_{T jet}");
682 fhJetLeadingRatioPts[icone][ipt]->SetXTitle("p_{T trigger} (GeV/c)");
1c5acb87 683
477d6cee 684 fhJetLeadingDeltaPhis[icone][ipt] = new TH2F("JetLeadingDeltaPhi"+lastnamehist,"#phi_{jet} - #phi_{leading} vs p_{T trigger}"+lastnametitle,nptbins,ptmin,ptmax,120,0,TMath::TwoPi());
685 fhJetLeadingDeltaPhis[icone][ipt]->SetYTitle("#Delta #phi (rad)");
686 fhJetLeadingDeltaPhis[icone][ipt]->SetXTitle("p_{T trigger} (GeV/c)");
1c5acb87 687
477d6cee 688 fhJetLeadingDeltaEtas[icone][ipt] = new TH2F("JetLeadingDeltaEta"+lastnamehist,"#eta_{jet} - #eta_{leading} vs p_{T trigger}"+lastnametitle,nptbins,ptmin,ptmax,120,-2,2);
689 fhJetLeadingDeltaEtas[icone][ipt]->SetYTitle("#Delta #eta");
690 fhJetLeadingDeltaEtas[icone][ipt]->SetXTitle("p_{T trigger} (GeV/c)");
1c5acb87 691
477d6cee 692 fhJetFFzs[icone][ipt] = new TH2F("JetFFz"+lastnamehist,"z = p_{T i charged}/p_{T trigger} vs p_{T trigger}", 120,0.,120.,200,0.,2);
693 fhJetFFzs[icone][ipt]->SetYTitle("z");
694 fhJetFFzs[icone][ipt]->SetXTitle("p_{T trigger}");
1c5acb87 695
477d6cee 696 fhJetFFxis[icone][ipt] = new TH2F("JetFFxi"+lastnamehist,"#xi = ln(p_{T trigger}/p_{T i charged}) vs p_{T trigger}", 120,0.,120.,100,0.,10.);
697 fhJetFFxis[icone][ipt]->SetYTitle("#xi");
698 fhJetFFxis[icone][ipt]->SetXTitle("p_{T trigger}");
1c5acb87 699
477d6cee 700 fhJetFFpts[icone][ipt] = new TH2F("JetFFpt"+lastnamehist,"p_{T charged hadron } in jet vs p_{T trigger}", 120,0.,120.,200,0.,50.);
701 fhJetFFpts[icone][ipt]->SetYTitle("p_{T charged hadron}");
702 fhJetFFpts[icone][ipt]->SetXTitle("p_{T trigger}");
1c5acb87 703
477d6cee 704 fhJetNTracksInCones[icone][ipt] = new TH2F("JetNTracksInCone"+lastnamehist,"N particles in cone vs p_{T trigger}"+lastnametitle,nptbins,ptmin,ptmax,5000,0, 5000);
705 fhJetNTracksInCones[icone][ipt]->SetYTitle("N tracks in jet cone");
706 fhJetNTracksInCones[icone][ipt]->SetXTitle("p_{T trigger} (GeV/c)");
1c5acb87 707
477d6cee 708 fOutCont->Add(fhJetPts[icone][ipt]) ;
709 fOutCont->Add(fhJetRatioPts[icone][ipt]) ;
710 fOutCont->Add(fhJetDeltaPhis[icone][ipt]) ;
711 fOutCont->Add(fhJetDeltaEtas[icone][ipt]) ;
712 fOutCont->Add(fhJetLeadingRatioPts[icone][ipt]) ;
713 fOutCont->Add(fhJetLeadingDeltaPhis[icone][ipt]) ;
714 fOutCont->Add(fhJetLeadingDeltaEtas[icone][ipt]) ;
715 fOutCont->Add(fhJetFFzs[icone][ipt]) ;
716 fOutCont->Add(fhJetFFxis[icone][ipt]) ;
717 fOutCont->Add(fhJetFFpts[icone][ipt]) ;
718 fOutCont->Add(fhJetNTracksInCones[icone][ipt]) ;
719
720 //Bkg Distributions
721 fhBkgPts[icone][ipt] = new TH2F("BkgPt"+lastnamehist,"p_{T bkg} vs p_{T trigger}"+lastnametitle,nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
722 fhBkgPts[icone][ipt]->SetYTitle("p_{T bkg}");
723 fhBkgPts[icone][ipt]->SetXTitle("p_{T trigger} (GeV/c)");
724
725 fhBkgRatioPts[icone][ipt] = new TH2F("BkgRatioPt"+lastnamehist,"p_{T bkg}/p_{T trigger} vs p_{T trigger}"+lastnametitle,nptbins,ptmin,ptmax,120,0,2);
726 fhBkgRatioPts[icone][ipt]->SetYTitle("p_{T bkg}/p_{T trigger}");
727 fhBkgRatioPts[icone][ipt]->SetXTitle("p_{T trigger} (GeV/c)");
728
729 fhBkgDeltaPhis[icone][ipt] = new TH2F("BkgDeltaPhi"+lastnamehist,"#phi_{bkg} - #phi_{trigger} vs p_{T trigger}"+lastnametitle,nptbins,ptmin,ptmax,120,0,TMath::TwoPi());
730 fhBkgDeltaPhis[icone][ipt]->SetYTitle("#Delta #phi (rad)");
731 fhBkgDeltaPhis[icone][ipt]->SetXTitle("p_{T trigger} (GeV/c)");
732
733 fhBkgDeltaEtas[icone][ipt] = new TH2F("BkgDeltaEta"+lastnamehist,"#eta_{bkg} - #eta_{trigger} vs p_{T trigger}"+lastnametitle,nptbins,ptmin,ptmax,120,-2,2);
734 fhBkgDeltaEtas[icone][ipt]->SetYTitle("#Delta #eta");
735 fhBkgDeltaEtas[icone][ipt]->SetXTitle("p_{T trigger} (GeV/c)");
736
737 fhBkgLeadingRatioPts[icone][ipt] = new TH2F("BkgLeadingRatioPt"+lastnamehist,"p_{T bkg} vs p_{T trigger}"+lastnametitle,nptbins,ptmin,ptmax,120,0,2);
738 fhBkgLeadingRatioPts[icone][ipt]->SetYTitle("p_{T leading}/p_{T bkg}");
739 fhBkgLeadingRatioPts[icone][ipt]->SetXTitle("p_{T trigger} (GeV/c)");
740
741 fhBkgLeadingDeltaPhis[icone][ipt] = new TH2F("BkgLeadingDeltaPhi"+lastnamehist,"#phi_{bkg} - #phi_{leading} vs p_{T trigger}"+lastnametitle,nptbins,ptmin,ptmax,120,0,TMath::TwoPi());
742 fhBkgLeadingDeltaPhis[icone][ipt]->SetYTitle("#Delta #phi (rad)");
743 fhBkgLeadingDeltaPhis[icone][ipt]->SetXTitle("p_{T trigger} (GeV/c)");
744
745 fhBkgLeadingDeltaEtas[icone][ipt] = new TH2F("BkgLeadingDeltaEta"+lastnamehist,"#eta_{bkg} - #eta_{leading} vs p_{T trigger}"+lastnametitle,nptbins,ptmin,ptmax,120,-2,2);
746 fhBkgLeadingDeltaEtas[icone][ipt]->SetYTitle("#Delta #eta");
747 fhBkgLeadingDeltaEtas[icone][ipt]->SetXTitle("p_{T trigger} (GeV/c)");
748
749 fhBkgFFzs[icone][ipt] = new TH2F("BkgFFz"+lastnamehist,"z = p_{T i charged}/p_{T trigger} vs p_{T trigger}", 120,0.,120.,200,0.,2);
750 fhBkgFFzs[icone][ipt]->SetYTitle("z");
751 fhBkgFFzs[icone][ipt]->SetXTitle("p_{T trigger}");
752
753 fhBkgFFxis[icone][ipt] = new TH2F("BkgFFxi"+lastnamehist,"#xi = ln(p_{T trigger}/p_{T i charged}) vs p_{T trigger}", 120,0.,120.,100,0.,10.);
754 fhBkgFFxis[icone][ipt]->SetYTitle("#xi");
755 fhBkgFFxis[icone][ipt]->SetXTitle("p_{T trigger}");
756
757 fhBkgFFpts[icone][ipt] = new TH2F("BkgFFpt"+lastnamehist,"p_{T charged hadron} in jet vs p_{T trigger}", 120,0.,120.,200,0.,50.);
758 fhBkgFFpts[icone][ipt]->SetYTitle("p_{T charged hadron}");
759 fhBkgFFpts[icone][ipt]->SetXTitle("p_{T trigger}");
760
761 fhBkgNTracksInCones[icone][ipt] = new TH2F("BkgNTracksInCone"+lastnamehist,"N particles in cone vs p_{T trigger}"+lastnametitle,nptbins,ptmin,ptmax,5000,0, 5000);
762 fhBkgNTracksInCones[icone][ipt]->SetYTitle("N tracks in bkg cone");
763 fhBkgNTracksInCones[icone][ipt]->SetXTitle("p_{T trigger} (GeV/c)");
764
765 fOutCont->Add(fhBkgPts[icone][ipt]) ;
766 fOutCont->Add(fhBkgRatioPts[icone][ipt]) ;
767 fOutCont->Add(fhBkgDeltaPhis[icone][ipt]) ;
768 fOutCont->Add(fhBkgDeltaEtas[icone][ipt]) ;
769 fOutCont->Add(fhBkgLeadingRatioPts[icone][ipt]) ;
770 fOutCont->Add(fhBkgLeadingDeltaPhis[icone][ipt]) ;
771 fOutCont->Add(fhBkgLeadingDeltaEtas[icone][ipt]) ;
772 fOutCont->Add(fhBkgFFzs[icone][ipt]) ;
773 fOutCont->Add(fhBkgFFxis[icone][ipt]) ;
774 fOutCont->Add(fhBkgFFpts[icone][ipt]) ;
775 fOutCont->Add(fhBkgNTracksInCones[icone][ipt]) ;
1c5acb87 776
477d6cee 777 }//ipt
778 } //icone
779 }//If we want to study any cone or pt threshold
780
1e1befe8 781 if(GetDebug() > 2){
477d6cee 782 printf("AliAnaParticleJetLeadingConeCorrelation::GetCreateOutputObjects() - All histograms names : \n");
477d6cee 783 for(Int_t i = 0 ; i< fOutCont->GetEntries(); i++)
1e1befe8 784 printf("Histo i %d name %s\n",i,((fOutCont->At(i))->GetName()));
477d6cee 785 //cout<< (fOutCont->At(i))->GetName()<<endl;
786 }
787
788 return fOutCont;
a3aebfff 789
1c5acb87 790}
791
792//____________________________________________________________________________
793Bool_t AliAnaParticleJetLeadingConeCorrelation::GetLeadingParticle(AliAODPWG4ParticleCorrelation *particle, TLorentzVector & pLeading)
477d6cee 794 const {
1c5acb87 795 //Search Charged or Neutral leading particle, select the highest one and fill AOD
796
797 TLorentzVector pLeadingCh(0,0,0,0) ;
798 TLorentzVector pLeadingPi0(0,0,0,0) ;
799
800 GetLeadingCharge(particle, pLeadingCh) ;
801 if(!fJetsOnlyInCTS) GetLeadingPi0(particle, pLeadingPi0) ;
802
803 Double_t ptch = pLeadingCh.Pt();
804 Double_t ptpi = pLeadingPi0.Pt();
805
806 if (ptch > 0 || ptpi > 0){
807 if((ptch >= ptpi)){
477d6cee 808 if(GetDebug() > 1)printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingParticle() - Leading found in CTS \n");
1c5acb87 809 pLeading = pLeadingCh;
477d6cee 810 if(GetDebug() > 1) printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingParticle() - Found Leading: pt %f, phi %f deg, eta %f\n",
811 pLeading.Pt(),pLeading.Phi()*TMath::RadToDeg(),pLeading.Eta()) ;
1c5acb87 812 //Put leading in AOD
477d6cee 813 particle->SetLeading(pLeadingCh);
814 particle->SetLeadingDetector("CTS");
815 return kTRUE;
1c5acb87 816 }
817 else{
477d6cee 818 if(GetDebug() > 1)printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingParticle() - Leading found in EMCAL \n");
1c5acb87 819 pLeading = pLeadingPi0;
477d6cee 820 if(GetDebug() > 1) printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingParticle() - Found Leading: pt %f, phi %f, eta %f\n",
821 pLeading.Pt(),pLeading.Phi()*TMath::RadToDeg(),pLeading.Eta()) ;
1c5acb87 822 //Put leading in AOD
823 particle->SetLeading(pLeadingPi0);
824 particle->SetLeadingDetector("EMCAL");
825 return kTRUE;
826 }
827 }
828
477d6cee 829 if(GetDebug() > 1)printf ("AliAnaParticleJetLeadingConeCorrelation::GetLeadingParticle() - NO LEADING PARTICLE FOUND \n");
1c5acb87 830
831 return kFALSE;
832
833}
834
835//____________________________________________________________________________
836void AliAnaParticleJetLeadingConeCorrelation::GetLeadingCharge(AliAODPWG4ParticleCorrelation * particle, TLorentzVector & pLeading) const
837{
838 //Search for the charged particle with highest pt and with
839 //Phi=Phi_trigger-Pi and pT=0.1E_gamma
840
841 if(GetAODCTS()){
842 Double_t ptTrig = particle->Pt();
843 Double_t phiTrig = particle->Phi();
844 Double_t rat = -100 ;
845 Double_t ptl = -100 ;
846 Double_t phil = -100 ;
847 Double_t pt = -100.;
848 Double_t phi = -100.;
849 TVector3 p3;
850
851 for(Int_t ipr = 0;ipr < GetAODCTS()->GetEntriesFast() ; ipr ++ ){
852 AliAODTrack* track = (AliAODTrack *)(GetAODCTS()->At(ipr)) ;
853 p3.SetXYZ(track->Px(),track->Py(),track->Pz());
854 pt = p3.Pt();
855 phi = p3.Phi() ;
856 if(phi<0) phi+=TMath::TwoPi();
857 rat = pt/ptTrig ;
477d6cee 858
1c5acb87 859 //Selection within angular and energy limits
860 if(((phiTrig-phi) > fDeltaPhiMinCut) && ((phiTrig-phi)<fDeltaPhiMaxCut) &&
861 (rat > fLeadingRatioMinCut) && (rat < fLeadingRatioMaxCut) && (pt > ptl)) {
862 phil = phi ;
863 ptl = pt ;
864 pLeading.SetVect(p3);
865 }
866 }// track loop
867
477d6cee 868 if(GetDebug() > 1&& ptl>0 ) printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingCharge() - Leading in CTS: pt %f eta %f phi %f pt/ptTrig %f \n",
869 ptl, pLeading.Eta(), phil,ptl/ptTrig) ;
1c5acb87 870
871 }//CTS list exist
872}
873
874//____________________________________________________________________________
875void AliAnaParticleJetLeadingConeCorrelation::GetLeadingPi0(AliAODPWG4ParticleCorrelation * particle, TLorentzVector & pLeading) const
876{
877 //Search for the neutral pion with highest pt and with
878 //Phi=Phi_trigger-Pi and pT=0.1E_gamma
477d6cee 879
1c5acb87 880 if(GetAODEMCAL()){
881 Double_t ptTrig = particle->Pt();
882 Double_t phiTrig = particle->Phi();
883 Double_t rat = -100 ;
884 Double_t ptl = -100 ;
885 Double_t phil = -100 ;
886 Double_t pt = -100.;
887 Double_t phi = -100.;
888
889 TLorentzVector gammai;
890 TLorentzVector gammaj;
891
892 Double_t vertex[] = {0,0,0};
893 if(!GetReader()->GetDataType()== AliCaloTrackReader::kMC) GetReader()->GetVertex(vertex);
894
895 //Cluster loop, select pairs with good pt, phi and fill AODs or histograms
896 for(Int_t iclus = 0;iclus < GetAODEMCAL()->GetEntriesFast() ; iclus ++ ){
897 AliAODCaloCluster * calo = (AliAODCaloCluster *)(GetAODEMCAL()->At(iclus)) ;
898
899 //Cluster selection, not charged, with photon or pi0 id and in fidutial cut
900 Int_t pdgi=0;
901 if(!SelectCluster(calo,vertex, gammai, pdgi)) continue ;
902
477d6cee 903 if(GetDebug() > 2) printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingPi0() - Neutral cluster: pt %f, phi %f \n",
904 gammai.Pt(),gammai.Phi());
1c5acb87 905
906 //2 gamma overlapped, found with PID
907 if(pdgi == AliCaloPID::kPi0){
908 pt = gammai.Pt();
909 rat = pt/ptTrig;
910 phi = gammai.Phi();
911 if(phi<0) phi+=TMath::TwoPi();
912
913 //Selection within angular and energy limits
914 if(ptl > pt && rat > fLeadingRatioMinCut && rat < fLeadingRatioMaxCut &&
915 (phiTrig-phil) > fDeltaPhiMinCut && (phiTrig-phil) < fDeltaPhiMaxCut )
916 {
917 phi = phil ;
918 pt = ptl ;
919 pLeading.SetPxPyPzE(gammai.Px(),gammai.Py(),gammai.Pz(),gammai.E());
920 }// cuts
921 }// pdg = AliCaloPID::kPi0
477d6cee 922 //Make invariant mass analysis
1c5acb87 923 else if(pdgi == AliCaloPID::kPhoton){
924 //Search the photon companion in case it comes from a Pi0 decay
925 //Apply several cuts to select the good pair
926 for(Int_t jclus = iclus+1; jclus < GetAODEMCAL()->GetEntriesFast() ; jclus ++ ){
927 AliAODCaloCluster * calo2 = (AliAODCaloCluster *) (GetAODEMCAL()->At(jclus)) ;
928
929 //Cluster selection, not charged with photon or pi0 id and in fidutial cut
930 Int_t pdgj=0;
931 if(!SelectCluster(calo2,vertex, gammaj, pdgj)) continue ;
932
933 if(pdgj == AliCaloPID::kPhoton ){
934
935 pt = (gammai+gammaj).Pt();
936 phi = (gammai+gammaj).Phi();
937 rat = pt/ptTrig;
938
939 //Selection within angular and energy limits
940 if(ptl > pt && rat > fLeadingRatioMinCut && rat < fLeadingRatioMaxCut &&
941 (phiTrig-phil) > fDeltaPhiMinCut && (phiTrig-phil) < fDeltaPhiMaxCut ){
942 //Select good pair (aperture and invariant mass)
943 if(GetNeutralMesonSelection()->SelectPair(gammai, gammaj)){
944 phi = phil ;
945 pt = ptl ;
946 pLeading=(gammai+gammaj);
947 }//pi0 selection
948
477d6cee 949 if(GetDebug() > 3 ) printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingPi0() - Neutral Hadron Correlation: Selected gamma pair: pt %2.2f, phi %2.2f, eta %2.2f, M %2.3f\n",
1c5acb87 950 (gammai+gammaj).Pt(),(gammai+gammaj).Phi(),(gammai+gammaj).Eta(), (gammai+gammaj).M());
951 }//Pair selected as leading
952 }//if pair of gammas
953 }//2nd loop
954 }// if pdg = 22
955 }// 1st Loop
956
477d6cee 957 if(GetDebug()>2 && pLeading.Pt() >0 ) printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingPi0() - Leading EMCAL: pt %f eta %f phi %f pt/Eg %f \n", pLeading.Pt(), pLeading.Eta(), pLeading.Phi(), pLeading.Pt()/ptTrig) ;
1c5acb87 958
959 }//EMCAL list exists
960
961}
962
963//____________________________________________________________________________
964void AliAnaParticleJetLeadingConeCorrelation::InitParameters()
965{
477d6cee 966 //Initialize the parameters of the analysis.
a3aebfff 967 SetInputAODName("PWG4Particle");
968 SetAODRefArrayName("JetLeadingCone");
969 AddToHistogramsName("AnaJetLCCorr_");
970
1c5acb87 971 fJetsOnlyInCTS = kFALSE ;
972 fPbPb = kFALSE ;
973 fReMakeJet = kFALSE ;
477d6cee 974
1c5acb87 975 //Leading selection parameters
976 fDeltaPhiMinCut = 2.9 ;
977 fDeltaPhiMaxCut = 3.4 ;
978 fLeadingRatioMinCut = 0.1;
979 fLeadingRatioMaxCut = 1.5;
980
981 //Jet selection parameters
982 //Fixed cut
983 fJetRatioMaxCut = 1.2 ;
984 fJetRatioMinCut = 0.3 ;
985 fJetCTSRatioMaxCut = 1.2 ;
986 fJetCTSRatioMinCut = 0.3 ;
987 fSelect = 0 ; //0, Accept all jets, 1, selection depends on energy, 2 fixed selection
988
477d6cee 989 fSelectIsolated = kFALSE;
990
1c5acb87 991 //Cut depending on gamma energy
992 fPtTriggerSelectionCut = 10.; //For Low pt jets+BKG, another limits applied
993 //Reconstructed jet energy dependence parameters
994 //e_jet = a1+e_gamma b2.
995 //Index 0-> Pt>2 GeV r = 0.3; Index 1-> Pt>0.5 GeV r = 0.3
996 fJetE1[0] = -5.75; fJetE1[1] = -4.1;
997 fJetE2[0] = 1.005; fJetE2[1] = 1.05;
998
999 //Reconstructed sigma of jet energy dependence parameters
1000 //s_jet = a1+e_gamma b2.
1001 //Index 0-> Pt>2 GeV r = 0.3; Index 1-> Pt>0.5 GeV r = 0.3
1002 fJetSigma1[0] = 2.65; fJetSigma1[1] = 2.75;
1003 fJetSigma2[0] = 0.0018; fJetSigma2[1] = 0.033;
1004
1005 //Background mean energy and RMS
1006 //Index 0-> No BKG; Index 1-> BKG > 2 GeV;
1007 //Index 2-> (low pt jets)BKG > 0.5 GeV;
1008 //Index > 2, same for CTS conf
1009 fBkgMean[0] = 0.; fBkgMean[1] = 8.8 ; fBkgMean[2] = 69.5;
1010 fBkgMean[3] = 0.; fBkgMean[4] = 6.4; fBkgMean[5] = 48.6;
1011 fBkgRMS[0] = 0.; fBkgRMS[1] = 7.5; fBkgRMS[2] = 22.0;
1012 fBkgRMS[3] = 0.; fBkgRMS[4] = 5.4; fBkgRMS[5] = 13.2;
1013
1014 //Factor x of min/max = E -+ x * sigma. Obtained after selecting the
1015 //limits for monoenergetic jets.
1016 //Index 0-> No BKG; Index 1-> BKG > 2 GeV;
1017 //Index 2-> (low pt jets) BKG > 0.5 GeV;
1018 //Index > 2, same for CTS conf
1019
1020 fJetXMin1[0] =-0.69 ; fJetXMin1[1] = 0.39 ; fJetXMin1[2] =-0.88 ;
1021 fJetXMin1[3] =-2.0 ; fJetXMin1[4] =-0.442 ; fJetXMin1[5] =-1.1 ;
1022 fJetXMin2[0] = 0.066; fJetXMin2[1] = 0.038; fJetXMin2[2] = 0.034;
1023 fJetXMin2[3] = 0.25 ; fJetXMin2[4] = 0.113; fJetXMin2[5] = 0.077 ;
1024 fJetXMax1[0] =-3.8 ; fJetXMax1[1] =-0.76 ; fJetXMax1[2] =-3.6 ;
1025 fJetXMax1[3] =-2.7 ; fJetXMax1[4] =-1.21 ; fJetXMax1[5] =-3.7 ;
1026 fJetXMax2[0] =-0.012; fJetXMax2[1] =-0.022; fJetXMax2[2] = 0.016;
1027 fJetXMax2[3] =-0.024; fJetXMax2[4] =-0.008; fJetXMax2[5] = 0.027;
1028
1029
1030 //Different cones and pt thresholds to construct the jet
1031
1032 fJetCone = 0.3 ;
1033 fJetPtThreshold = 0.5 ;
1034 fJetPtThresPbPb = 2. ;
1035 fJetNCone = 4 ;
1036 fJetNPt = 4 ;
1037 fJetCones[0] = 0.2 ; fJetNameCones[0] = "02" ;
1038 fJetCones[1] = 0.3 ; fJetNameCones[1] = "03" ;
1039 fJetCones[2] = 0.4 ; fJetNameCones[2] = "04" ;
1040 fJetCones[2] = 0.5 ; fJetNameCones[2] = "05" ;
1041
1042 fJetPtThres[0] = 0.0 ; fJetNamePtThres[0] = "00" ;
1043 fJetPtThres[1] = 0.5 ; fJetNamePtThres[1] = "05" ;
1044 fJetPtThres[2] = 1.0 ; fJetNamePtThres[2] = "10" ;
1045 fJetPtThres[3] = 2.0 ; fJetNamePtThres[3] = "20" ;
1046}
1047
1048//__________________________________________________________________________-
1049Bool_t AliAnaParticleJetLeadingConeCorrelation::IsJetSelected(const Double_t ptTrig, const Double_t ptjet) const {
1050 //Given the pt of the jet and the trigger particle, select the jet or not
1051 //3 options, fSelect=0 accepts all, fSelect=1 selects jets depending on a
1052 //function energy dependent and fSelect=2 selects on simple fixed cuts
1053
1054 if(ptjet == 0) return kFALSE;
477d6cee 1055
1c5acb87 1056 Double_t rat = ptTrig / ptjet ;
1057
1058 //###############################################################
1059 if(fSelect == 0)
1060 return kTRUE; //Accept all jets, no restriction
1061 //###############################################################
1062 else if(fSelect == 1){
1063 //Check if the energy of the reconstructed jet is within an energy window
1064 //WARNING: to be rechecked, don't remember what all the steps mean
1065 Double_t par[6];
1066 Double_t xmax[2];
1067 Double_t xmin[2];
1068
1069 Int_t iCTS = 0;
1070 if(fJetsOnlyInCTS)
1071 iCTS = 3 ;
1072
1073 if(!fPbPb){
1074 //Phythia alone, jets with pt_th > 0.2, r = 0.3
1075 par[0] = fJetE1[0]; par[1] = fJetE2[0];
1076 //Energy of the jet peak
1077 //e_jet = fJetE1[0]+fJetE2[0]*e_gamma, simulation fit
1078 par[2] = fJetSigma1[0]; par[3] = fJetSigma2[0];
1079 //Sigma of the jet peak
1080 //sigma_jet = fJetSigma1[0]+fJetSigma2[0]*e_gamma, simulation fit
1081 par[4] = fBkgMean[0 + iCTS]; par[5] = fBkgRMS[0 + iCTS];
1082 //Parameters reserved for PbPb bkg.
1083 xmax[0] = fJetXMax1[0 + iCTS]; xmax[1] = fJetXMax2[0 + iCTS];
1084 xmin[0] = fJetXMin1[0 + iCTS]; xmin[1] = fJetXMin2[0 + iCTS];
1085 //Factor that multiplies sigma to obtain the best limits,
1086 //by observation, of mono jet ratios (ptjet/ptTrig)
1087 //X_jet = fJetX1[0]+fJetX2[0]*e_gamma
1088
1089 }
1090 else{
1091 if(ptTrig > fPtTriggerSelectionCut){
1092 //Phythia +PbPb with pt_th > 2 GeV/c, r = 0.3
1093 par[0] = fJetE1[0]; par[1] = fJetE2[0];
1094 //Energy of the jet peak, same as in pp
1095 //e_jet = fJetE1[0]+fJetE2[0]*e_gamma, simulation fit
1096 par[2] = fJetSigma1[0]; par[3] = fJetSigma2[0];
1097 //Sigma of the jet peak, same as in pp
1098 //sigma_jet = fJetSigma1[0]+fJetSigma2[0]*e_gamma, simulation fit
1099 par[4] = fBkgMean[1 + iCTS]; par[5] = fBkgRMS[1 + iCTS];
1100 //Mean value and RMS of PbPb Bkg
1101 xmax[0] = fJetXMax1[1 + iCTS]; xmax[1] = fJetXMax2[1 + iCTS];
1102 xmin[0] = fJetXMin1[1 + iCTS]; xmin[1] = fJetXMin2[1 + iCTS];
1103 //Factor that multiplies sigma to obtain the best limits,
1104 //by observation, of mono jet ratios (ptjet/ptTrig) mixed with PbPb Bkg,
1105 //pt_th > 2 GeV, r = 0.3
1106 //X_jet = fJetX1[0]+fJetX2[0]*e_gamma
1107
1108 }
1109 else{
1110 //Phythia + PbPb with pt_th > 0.5 GeV/c, r = 0.3
1111 par[0] = fJetE1[1]; par[1] = fJetE2[1];
1112 //Energy of the jet peak, pt_th > 2 GeV/c, r = 0.3
1113 //e_jet = fJetE1[0]+fJetE2[0]*e_gamma, simulation fit
1114 par[2] = fJetSigma1[1]; par[3] = fJetSigma2[1];
1115 //Sigma of the jet peak, pt_th > 2 GeV/c, r = 0.3
1116 //sigma_jet = fJetSigma1[0]+fJetSigma2[0]*e_gamma, simulation fit
1117 par[4] = fBkgMean[2 + iCTS]; par[5] = fBkgRMS[2 + iCTS];
1118 //Mean value and RMS of PbPb Bkg in a 0.3 cone, pt > 2 GeV.
1119 xmax[0] = fJetXMax1[2 + iCTS]; xmax[1] = fJetXMax2[2 + iCTS];
1120 xmin[0] = fJetXMin1[2 + iCTS]; xmin[1] = fJetXMin2[2 + iCTS];
1121 //Factor that multiplies sigma to obtain the best limits,
1122 //by observation, of mono jet ratios (ptjet/ptTrig) mixed with PbPb Bkg,
1123 //pt_th > 2 GeV, r = 0.3
1124 //X_jet = fJetX1[0]+fJetX2[0]*e_gamma
1125
1126 }//If low pt jet in bkg
1127 }//if Bkg
1128
1129 //Calculate minimum and maximum limits of the jet ratio.
1130 Double_t min = CalculateJetRatioLimit(ptTrig, par, xmin);
1131 Double_t max = CalculateJetRatioLimit(ptTrig, par, xmax);
1132
a3aebfff 1133 if(GetDebug() > 3)printf("AliAnaParticleJetLeadingConeCorrelation::IsJetSelected() - Jet selection? : Limits min %f, max %f, pt_jet %f, pt_gamma %f, pt_jet / pt_gamma %f\n",min,max,ptjet,ptTrig,rat);
1c5acb87 1134
1135 if(( min < rat ) && ( max > ptjet/rat))
1136 return kTRUE;
1137 else
1138 return kFALSE;
1139 }//fSelect = 1
1140 //###############################################################
1141 else if(fSelect == 2){
1142 //Simple selection
1143 if(!fJetsOnlyInCTS){
1144 if((rat < fJetRatioMaxCut) && (rat > fJetRatioMinCut )) return kTRUE;
1145 }
1146 else{
1147 if((rat < fJetCTSRatioMaxCut) && (rat > fJetCTSRatioMinCut )) return kTRUE;
1148 }
1149 }// fSelect = 2
1150 //###############################################################
1151 else{
a3aebfff 1152 printf("AliAnaParticleJetLeadingConeCorrelation::IsJetSelected() - Jet selection option larger than 2, DON'T SELECT JETS\n");
1c5acb87 1153 return kFALSE ;
1154 }
1155
1156 return kFALSE;
1157
1158}
1159
1160//___________________________________________________________________
1161Bool_t AliAnaParticleJetLeadingConeCorrelation::IsParticleInJetCone(const Double_t eta, Double_t phi, const Double_t etal, Double_t phil)
1162 const {
1163 //Check if the particle is inside the cone defined by the leading particle
1164 //WARNING: To be rechecked
1165
1166 if(phi < 0) phi+=TMath::TwoPi();
1167 if(phil < 0) phil+=TMath::TwoPi();
1168 Double_t rad = 10000 + fJetCone;
1169
1170 if(TMath::Abs(phi-phil) <= (TMath::TwoPi() - fJetCone))
1171 rad = TMath::Sqrt(TMath::Power(eta-etal,2)+TMath::Power(phi-phil,2));
1172 else{
1173 if(phi-phil > TMath::TwoPi() - fJetCone)
1174 rad = TMath::Sqrt(TMath::Power(eta-etal,2)+TMath::Power((phi-TMath::TwoPi())-phil,2));
1175 if(phi-phil < -(TMath::TwoPi() - fJetCone))
1176 rad = TMath::Sqrt(TMath::Power(eta-etal,2)+TMath::Power((phi+TMath::TwoPi())-phil,2));
1177 }
1178
1179 if(rad < fJetCone) return kTRUE ;
1180 else return kFALSE ;
1181
1182}
1183
1184//__________________________________________________________________
1185void AliAnaParticleJetLeadingConeCorrelation::MakeAnalysisFillAOD()
1186{
1187 //Particle-Hadron Correlation Analysis, fill AODs
1188
477d6cee 1189 if(!GetInputAODBranch()){
1190 printf("AliAnaParticleJetLeadingConeCorrelation::MakeAnalysisFillAOD() - No input particles in AOD with name branch < %s > \n",
1191 GetInputAODName().Data());
1192 abort();
1193 }
1c5acb87 1194 if(GetDebug() > 1){
477d6cee 1195 printf("AliAnaParticleJetLeadingConeCorrelation::MakeAnalysisFillAOD() - Begin jet leading cone correlation analysis, fill AODs \n");
a3aebfff 1196 printf("AliAnaParticleJetLeadingConeCorrelation::MakeAnalysisFillAOD() - In particle branch aod entries %d\n", GetInputAODBranch()->GetEntriesFast());
1197 printf("AliAnaParticleJetLeadingConeCorrelation::MakeAnalysisFillAOD() - In CTS aod entries %d\n", GetAODCTS()->GetEntriesFast());
1198 printf("AliAnaParticleJetLeadingConeCorrelation::MakeAnalysisFillAOD() - In EMCAL aod entries %d\n", GetAODEMCAL()->GetEntriesFast());
1c5acb87 1199 }
1200
1201 TLorentzVector pLeading(0,0,0,0); //It will contain the kinematics of the found leading particle
1202
1203 //Loop on stored AOD particles, trigger
1204 Int_t naod = GetInputAODBranch()->GetEntriesFast();
1205 for(Int_t iaod = 0; iaod < naod ; iaod++){
a3aebfff 1206 AliAODPWG4ParticleCorrelation* particle = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
1207
1c5acb87 1208 //Search leading particles in CTS and EMCAL
1209 if(GetLeadingParticle(particle, pLeading)){
1210
1211 //Construct the jet around the leading, Fill AOD jet particle list, select jet
1212 //and fill AOD with jet and background
1213 MakeAODJet(particle, pLeading);
a3aebfff 1214
1c5acb87 1215 }//Leading found
1216 }//AOD trigger particle loop
1217
477d6cee 1218 if(GetDebug() >1)printf("AliAnaParticleJetLeadingConeCorrelation::MakeAnalysisFillAOD() - End of jet leading cone analysis, fill AODs \n");
1c5acb87 1219
1220}
1221
1222//__________________________________________________________________
1223void AliAnaParticleJetLeadingConeCorrelation::MakeAnalysisFillHistograms()
1224{
1225
1226 //Particle-Hadron Correlation Analysis, fill histograms
1227
477d6cee 1228 if(!GetInputAODBranch()){
1229 printf("AliAnaParticleJetLeadingConeCorrelation::MakeAnalysisFillHistograms() - No input particles in AOD with name branch < %s > \n",
1230 GetInputAODName().Data());
1231 abort();
1232 }
1c5acb87 1233 if(GetDebug() > 1){
477d6cee 1234 printf("AliAnaParticleJetLeadingConeCorrelation::MakeAnalysisFillHistograms() - Begin jet leading cone correlation analysis, fill histograms \n");
a3aebfff 1235 printf("AliAnaParticleJetLeadingConeCorrelation::MakeAnalysisFillHistograms() - In particle branch aod entries %d\n", GetInputAODBranch()->GetEntriesFast());
1236 printf("AliAnaParticleJetLeadingConeCorrelation::MakeAnalysisFillHistograms() - In CTS aod entries %d\n", GetAODCTS()->GetEntriesFast());
1237 printf("AliAnaParticleJetLeadingConeCorrelation::MakeAnalysisFillHistograms() - In EMCAL aod entries %d\n", GetAODEMCAL()->GetEntriesFast());
1c5acb87 1238 }
477d6cee 1239
1c5acb87 1240 TLorentzVector pLeading(0,0,0,0) ;
1241
1242 //Loop on stored AOD particles, trigger
1243 Int_t naod = GetInputAODBranch()->GetEntriesFast();
1244 for(Int_t iaod = 0; iaod < naod ; iaod++){
1245 AliAODPWG4ParticleCorrelation* particle = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
1246
477d6cee 1247 if(OnlyIsolated() && !particle->IsIsolated()) continue;
1248
1c5acb87 1249 Double_t pt = particle->Pt();
1250 Double_t phi = particle->Phi();
1251 Double_t eta = particle->Eta();
1252
1253 //Get leading particle, fill histograms
1254 pLeading = particle->GetLeading();
1255 TString det = particle->GetLeadingDetector();
477d6cee 1256
1c5acb87 1257 if(det!="" && pLeading.Pt() > 0){
1258 Double_t ptL = pLeading.Pt();
1259 Double_t phiL = pLeading.Phi();
1260 if(phiL < 0 ) phiL+=TMath::TwoPi();
1261 Double_t etaL = pLeading.Eta();
1262
477d6cee 1263 if(GetDebug() > 1) printf("AliAnaParticleJetLeadingConeCorrelation::MakeAnalysisFillHistograms() - Leading found in %s, with pt %3.2f, phi %2.2f, eta %2.2f\n",
1264 det.Data(), ptL, phiL, etaL);
1c5acb87 1265 if(det == "CTS"){
1266 fhChargedLeadingPt->Fill(pt,ptL);
1267 fhChargedLeadingPhi->Fill(pt,phiL);
1268 fhChargedLeadingEta->Fill(pt,etaL);
1269 fhChargedLeadingDeltaPt->Fill(pt,pt-ptL);
1270 fhChargedLeadingDeltaPhi->Fill(pt,phi-phiL);
1271 fhChargedLeadingDeltaEta->Fill(pt,eta-etaL);
1272 fhChargedLeadingRatioPt->Fill(pt,ptL/pt);
1273 }
1274 else if(det== "EMCAL"){
1275 fhNeutralLeadingPt->Fill(pt,ptL);
1276 fhNeutralLeadingPhi->Fill(pt,phiL);
1277 fhNeutralLeadingEta->Fill(pt,etaL);
1278 fhNeutralLeadingDeltaPt->Fill(pt,pt-ptL);
1279 fhNeutralLeadingDeltaPhi->Fill(pt,phi-phiL);
1280 fhNeutralLeadingDeltaEta->Fill(pt,eta-etaL);
1281 fhNeutralLeadingRatioPt->Fill(pt,ptL/pt);
1282 }
1283
1284 //Fill Jet histograms
1285 TLorentzVector bkg(0,0,0,0);
1286 TLorentzVector jet(0,0,0,0);
1287 if(!fSeveralConeAndPtCuts){//just fill histograms
1288 if(!fReMakeJet){
1289 jet=particle->GetCorrelatedJet();
1290 bkg=particle->GetCorrelatedBackground();
1291 }
1292 else MakeJetFromAOD(particle, pLeading, jet,bkg);
1293
1294 if(jet.Pt() > 0){//Jet was found
1295 FillJetHistos(particle, pLeading, jet,"Jet","");
1296 FillJetHistos(particle, pLeading, bkg,"Bkg","");
1297 }
1298 }
1299 else if(fSeveralConeAndPtCuts){
1300 for(Int_t icone = 0; icone<fJetNCone; icone++) {
1301 fJetCone=fJetCones[icone];
1302 for(Int_t ipt = 0; ipt<fJetNPt;ipt++) {
1303 TString lastname ="Cone"+ fJetNameCones[icone]+"Pt"+ fJetNamePtThres[ipt];
1304 fJetPtThreshold=fJetPtThres[ipt];
1305 MakeJetFromAOD(particle, pLeading, jet,bkg);
1306 if(jet.Pt() > 0) {//Jet was found
1307 FillJetHistos(particle, pLeading, jet,"Jet",lastname);
1308 FillJetHistos(particle, pLeading, bkg,"Bkg",lastname);
1309 }
1310 }//icone
1311 }//ipt
1312 }//fSeveralConeAndPtCuts
1313 }//Leading
1314 }//AOD trigger particle loop
1315
477d6cee 1316 if(GetDebug() >1)printf("AliAnaParticleJetLeadingConeCorrelation::MakeAnalysisFillHistograms() - End of jet leading cone analysis, fill histograms \n");
1c5acb87 1317
1318}
1319
1320//____________________________________________________________________________
1321void AliAnaParticleJetLeadingConeCorrelation::MakeAODJet(AliAODPWG4ParticleCorrelation *particle, const TLorentzVector pLeading)
1322const {
1323 //Fill the jet with the particles around the leading particle with
1324 //R=fJetCone and pt_th = fJetPtThres. Calculate the energy of the jet and
1325 //fill aod with found information
1326
1327 TLorentzVector bkg(0,0,0,0);
1328 TLorentzVector jet(0,0,0,0);
1329 TLorentzVector lv (0,0,0,0); //Temporal container for jet particles kinematics
1330
477d6cee 1331 Double_t ptTrig = particle->Pt();
1332 Double_t phiTrig = particle->Phi();
1333 Double_t phil = pLeading.Phi();
1c5acb87 1334 if(phil<0) phil+=TMath::TwoPi();
477d6cee 1335 Double_t etal = pLeading.Eta();
1336 Bool_t first = kTRUE;
1337 Bool_t firstbkg = kTRUE;
1338
1c5acb87 1339 //Different pt cut for jet particles in different collisions systems
1340 Float_t ptcut = fJetPtThreshold;
1341 if(fPbPb && !fSeveralConeAndPtCuts && ptTrig > fPtTriggerSelectionCut) ptcut = fJetPtThresPbPb ;
1342
1343 //Add charged particles to jet if they are in cone around the leading particle
1344 if(!GetAODCTS()) {
477d6cee 1345 printf("AliAnaParticleJetLeadingConeCorrelation::MakeAODJet() - Cannot construct jets without tracks, STOP analysis");
1c5acb87 1346 return;
1347 }
1348
1349 //Fill jet with tracks
1350 TVector3 p3;
a3aebfff 1351 //Initialize reference arrays that will contain jet and background tracks
1352 TRefArray * reftracks = new TRefArray;
1353 TRefArray * reftracksbkg = new TRefArray;
1354
1c5acb87 1355 for(Int_t ipr = 0;ipr < (GetAODCTS())->GetEntriesFast() ; ipr ++ ){
1356 AliAODTrack* track = (AliAODTrack *)((GetAODCTS())->At(ipr)) ;
1357 p3.SetXYZ(track->Px(),track->Py(),track->Pz());
1358
1359 //Particles in jet
1360 if(IsParticleInJetCone(p3.Eta(), p3.Phi(), etal, phil)){
477d6cee 1361
1362 if(first) {
a3aebfff 1363 new (reftracks) TRefArray(TProcessID::GetProcessWithUID(track));
477d6cee 1364 first = kFALSE;
1365 }
1366
a3aebfff 1367 reftracks->Add(track);
1368
1c5acb87 1369 if(p3.Pt() > ptcut ){
1370 lv.SetVect(p3);
1371 jet+=lv;
1372 }
1373 }
1374 //Background around (phi_gamma-pi, eta_leading)
1375 else if(IsParticleInJetCone(p3.Eta(),p3.Phi(),etal, phiTrig)) {
477d6cee 1376
1377 if(firstbkg) {
a3aebfff 1378 new (reftracksbkg) TRefArray(TProcessID::GetProcessWithUID(track));
477d6cee 1379 firstbkg = kFALSE;
1380 }
1381
a3aebfff 1382 reftracksbkg->Add(track);
1383
1c5acb87 1384 if(p3.Pt() > ptcut ){
1385 lv.SetVect(p3);
1386 bkg+=lv;
1387 }
1388 }
1389 }//Track loop
1390
a3aebfff 1391 //Add referenced tracks to AOD
1392 if(reftracks->GetEntriesFast() > 0 ){
1393 reftracks->SetName(GetAODRefArrayName()+"Tracks");
1394 particle->AddRefArray(reftracks);
1395 }
1396 if(reftracksbkg->GetEntriesFast() > 0 ){
1397 reftracksbkg->SetName(GetAODRefArrayName()+"TracksBkg");
1398 particle->AddRefArray(reftracksbkg);
1399 }
1400
1c5acb87 1401 //Add neutral particles to jet
a3aebfff 1402 //Initialize reference arrays that will contain jet and background tracks
1403 TRefArray * refclusters = new TRefArray;
1404 TRefArray * refclustersbkg = new TRefArray;
1c5acb87 1405 if(!fJetsOnlyInCTS && GetAODEMCAL()){
1406
1407 Double_t vertex[] = {0,0,0};
1408 if(!GetReader()->GetDataType()== AliCaloTrackReader::kMC) GetReader()->GetVertex(vertex);
477d6cee 1409
1410 first = kTRUE;
1411 firstbkg = kTRUE;
a3aebfff 1412
1c5acb87 1413 for(Int_t iclus = 0;iclus < (GetAODEMCAL())->GetEntriesFast() ; iclus ++ ){
1414 AliAODCaloCluster * calo = (AliAODCaloCluster *) (GetAODEMCAL()->At(iclus)) ;
1415
1416 //Cluster selection, not charged
1417 if(calo->GetNTracksMatched() > 0) continue ;
1418
1419 calo->GetMomentum(lv,vertex);
1420 //Particles in jet
1421 if(IsParticleInJetCone(lv.Eta(),lv.Phi(), etal, phil)){
a3aebfff 1422
477d6cee 1423 if(first) {
a3aebfff 1424 new (refclusters) TRefArray(TProcessID::GetProcessWithUID(calo));
477d6cee 1425 first = kFALSE;
1426 }
a3aebfff 1427
1428 refclusters->Add(calo);
1429
1c5acb87 1430 if(lv.Pt() > ptcut ) jet+=lv;
1431 }
1432 //Background around (phi_gamma-pi, eta_leading)
1433 else if(IsParticleInJetCone(lv.Eta(),lv.Phi(),etal, phiTrig)){
a3aebfff 1434
477d6cee 1435 if(firstbkg) {
a3aebfff 1436 new (refclustersbkg) TRefArray(TProcessID::GetProcessWithUID(calo));
477d6cee 1437 firstbkg = kFALSE;
1438 }
1439
a3aebfff 1440 refclustersbkg->Add(calo);
1441
1c5acb87 1442 if(lv.Pt() > ptcut ) bkg+=lv;
1443 }
1444 }//cluster loop
1445 }//jets with neutral particles
1446
a3aebfff 1447 //Add referenced clusters to AOD
1448 if(refclusters->GetEntriesFast() > 0 ){
1449 refclusters->SetName(GetAODRefArrayName()+"Clusters");
1450 particle->AddRefArray(refclusters);
1451 }
1452 if(refclustersbkg->GetEntriesFast() > 0 ){
1453 refclustersbkg->SetName(GetAODRefArrayName()+"ClustersBkg");
1454 particle->AddRefArray(refclustersbkg);
1455 }
1456
1c5acb87 1457 //If there is any jet found, select after some criteria and
1458 //and fill AOD with corresponding TLorentzVector kinematics
1459 if(IsJetSelected(particle->Pt(), jet.Pt())) {
1460 particle->SetCorrelatedJet(jet);
1461 particle->SetCorrelatedBackground(bkg);
477d6cee 1462 if(GetDebug()>1) printf("AliAnaParticleJetLeadingConeCorrelation::MakeAODJet() - Found jet: Trigger pt %f, Jet pt %f, Bkg pt %f\n",ptTrig,jet.Pt(),bkg.Pt());
1c5acb87 1463 }
1464
1465}
1466
1467//____________________________________________________________________________
1468void AliAnaParticleJetLeadingConeCorrelation::MakeJetFromAOD(AliAODPWG4ParticleCorrelation *particle, const TLorentzVector pLeading, TLorentzVector & jet, TLorentzVector & bkg)
1469const {
1470 //Fill the jet with the particles around the leading particle with
1471 //R=fJetCone and pt_th = fJetPtThres. Calculate the energy of the jet and
1472 //fill aod tlorentzvectors with jet and bakcground found
1473
1474 TLorentzVector lv (0,0,0,0); //Temporal container for jet particles kinematics
1475
1476 Double_t ptTrig = particle->Pt();
1477 Double_t phiTrig = particle->Phi();
1478 Double_t phil = pLeading.Phi();
1479 Double_t etal = pLeading.Eta();
1480
a3aebfff 1481 TRefArray * refclusters = particle->GetRefArray(GetAODRefArrayName()+"Clusters");
1482 TRefArray * reftracks = particle->GetRefArray(GetAODRefArrayName()+"Tracks");
1483 TRefArray * refclustersbkg = particle->GetRefArray(GetAODRefArrayName()+"ClustersBkg");
1484 TRefArray * reftracksbkg = particle->GetRefArray(GetAODRefArrayName()+"TracksBkg");
1485
1c5acb87 1486 //Different pt cut for jet particles in different collisions systems
1487 Float_t ptcut = fJetPtThreshold;
1488 if(fPbPb && !fSeveralConeAndPtCuts && ptTrig > fPtTriggerSelectionCut) ptcut = fJetPtThresPbPb ;
a3aebfff 1489
1c5acb87 1490 //Fill jet with tracks
1491 //Particles in jet
1492 TVector3 p3;
a3aebfff 1493 if(reftracks){
1494 for(Int_t ipr = 0;ipr < reftracks->GetEntriesFast() ; ipr ++ ){
1495 AliAODTrack* track = (AliAODTrack *) reftracks->At(ipr) ;
1496 p3.SetXYZ(track->Px(),track->Py(),track->Pz());
1497 if(p3.Pt() > ptcut && IsParticleInJetCone(p3.Eta(), p3.Phi(), etal, phil) ){
1498 lv.SetVect(p3);
1499 jet+=lv;
1500 }
1501 }//jet Track loop
1502 }
1c5acb87 1503 //Particles in background
a3aebfff 1504 if(reftracksbkg){
1505 for(Int_t ipr = 0;ipr < reftracksbkg->GetEntriesFast() ; ipr ++ ){
1506 AliAODTrack* track = (AliAODTrack *) reftracksbkg->At(ipr) ;
1507 p3.SetXYZ(track->Px(),track->Py(),track->Pz());
1508 if(p3.Pt() > ptcut && IsParticleInJetCone(p3.Eta(),p3.Phi(),etal, phiTrig) ) {
1509 lv.SetVect(p3);
1510 bkg+=lv;
1511 }
1512 }//background Track loop
1513 }
1514
1c5acb87 1515 //Add neutral particles to jet
a3aebfff 1516 if(!fJetsOnlyInCTS && refclusters){
1c5acb87 1517
1518 Double_t vertex[] = {0,0,0};
1519 if(!GetReader()->GetDataType()== AliCaloTrackReader::kMC) GetReader()->GetVertex(vertex);
1520
1521 //Loop on jet particles
a3aebfff 1522 if(refclusters){
1523 for(Int_t iclus = 0;iclus < refclusters->GetEntriesFast() ; iclus ++ ){
1524 AliAODCaloCluster * calo = (AliAODCaloCluster *) refclusters->At(iclus) ;
1525 calo->GetMomentum(lv,vertex);
1526 if(lv.Pt() > ptcut && IsParticleInJetCone(lv.Eta(),lv.Phi(), etal, phil)) jet+=lv;
1527 }//jet cluster loop
1528 }
1529
1530 //Loop on background particles
1531 if(refclustersbkg){
1532 for(Int_t iclus = 0;iclus < refclustersbkg->GetEntriesFast() ; iclus ++ ){
1533 AliAODCaloCluster * calo = (AliAODCaloCluster *) refclustersbkg->At(iclus) ;
1534 calo->GetMomentum(lv,vertex);
1535 if( lv.Pt() > ptcut && IsParticleInJetCone(lv.Eta(),lv.Phi(),etal, phiTrig)) bkg+=lv;
1536 }//background cluster loop
1537 }
1c5acb87 1538 }//clusters in jet
a3aebfff 1539
1c5acb87 1540 //If there is any jet found, leave jet and bkg as they are,
1541 //if not set them to 0.
1542 if(!IsJetSelected(particle->Pt(), jet.Pt())) {
1543 jet.SetPxPyPzE(0.,0.,0.,0.);
1544 bkg.SetPxPyPzE(0.,0.,0.,0.);
1545 }
1546 else
a3aebfff 1547 if(GetDebug()>1) printf("AliAnaParticleJetLeadingConeCorrelation::MakeJetFromAOD() - Found jet: Trigger pt %f, Jet pt %f, Bkg pt %f\n",ptTrig,jet.Pt(),bkg.Pt());
1c5acb87 1548
1549}
1550
1551//____________________________________________________________________________
1552Bool_t AliAnaParticleJetLeadingConeCorrelation::SelectCluster(AliAODCaloCluster * calo, Double_t *vertex, TLorentzVector & mom, Int_t & pdg) const {
1553 //Select cluster depending on its pid and acceptance selections
1554
1555 //Skip matched clusters with tracks
1556 if(calo->GetNTracksMatched() > 0) return kFALSE;
a3aebfff 1557
1558 //Check PID
1559 calo->GetMomentum(mom,vertex);//Assume that come from vertex in straight line
1560 pdg = AliCaloPID::kPhoton;
1561 if(IsCaloPIDOn()){
1562 //Get most probable PID, 2 options check PID weights (in MC this option is mandatory)
1563 //or redo PID, recommended option for EMCal.
1564 if(!IsCaloPIDRecalculationOn() || GetReader()->GetDataType() == AliCaloTrackReader::kMC )
1565 pdg = GetCaloPID()->GetPdg("EMCAL",calo->PID(),mom.E());//PID with weights
1566 else
1567 pdg = GetCaloPID()->GetPdg("EMCAL",mom,calo);//PID recalculated
1568
1569 if(GetDebug() > 1) printf("AliAnaParticleJetLeadingConeCorrelation::SelectCluster() - PDG of identified particle %d\n",pdg);
1570 //If it does not pass pid, skip
1571 if(pdg != AliCaloPID::kPhoton || pdg != AliCaloPID::kPi0)
1572 return kFALSE ;
1573 }//CaloPID
1574
1c5acb87 1575 //Check acceptance selection
a3aebfff 1576 if(IsFidutialCutOn()){
1577 Bool_t in = GetFidutialCut()->IsInFidutialCut(mom,"EMCAL") ;
1578 if(! in ) return kFALSE ;
1579 }
1580
1581 if(GetDebug() > 1) printf("AliAnaParticleJetLeadingConeCorrelation::SelectCluster() - Cluster selection cuts passed: pT %3.2f, pdg %d\n",mom.Pt(), pdg);
1582
1583 return kTRUE;
1584
1c5acb87 1585}
1586
1587//__________________________________________________________________
1588void AliAnaParticleJetLeadingConeCorrelation::Print(const Option_t * opt) const
1589{
1590
1591 //Print some relevant parameters set for the analysis
1592 if(! opt)
1593 return;
1594
a3aebfff 1595 printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
1596 AliAnaPartCorrBaseClass::Print(" ");
1597
1c5acb87 1598 if(fJetsOnlyInCTS)printf("Jets reconstructed in CTS \n");
1599 else printf("Jets reconstructed in CTS+EMCAL \n");
1600
1601 if(fPbPb) printf("PbPb events, pT cut in jet cone energy reconstruction %2.1f \n", fJetPtThreshold);
1602 else printf("pp events, pT cut in jet cone energy reconstruction %2.1f \n", fJetPtThresPbPb);
1603
1604 printf("If pT of trigger < %f, select jets as in pp? \n", fPtTriggerSelectionCut);
1605
1606 printf("Phi gamma-Leading < %3.2f\n", fDeltaPhiMaxCut) ;
1607 printf("Phi gamma-Leading > %3.2f\n", fDeltaPhiMinCut) ;
1608 printf("pT Leading / pT Trigger < %3.2f\n", fLeadingRatioMaxCut) ;
1609 printf("pT Leading / pT Trigger > %3.2f\n", fLeadingRatioMinCut) ;
1610
1611 if(fSelect == 2){
1612 printf("pT Jet / pT Gamma < %3.2f\n", fJetRatioMaxCut) ;
1613 printf("pT Jet / pT Gamma > %3.2f\n", fJetRatioMinCut) ;
1614 printf("pT Jet (Only CTS)/ pT Trigger < %3.2f\n", fJetCTSRatioMaxCut) ;
1615 printf("pT Jet (Only CTS)/ pT Trigger > %3.2f\n", fJetCTSRatioMinCut) ;
1616 }
1617 else if(fSelect == 0)
1618 printf("Accept all reconstructed jets \n") ;
1619 else if(fSelect == 1)
1620 printf("Accept jets depending on trigger energy \n") ;
1621 else
1622 printf("Wrong jet selection option: %d \n", fSelect) ;
477d6cee 1623
1624 printf("Isolated Trigger? %d\n", fSelectIsolated) ;
1c5acb87 1625
1626}