Add sigma2 jet shape and fill constituent info. for subtracted jets
[u/mrichter/AliRoot.git] / PWGJE / UserTasks / AliAnalysisTaskJetHadronCorrelation.cxx
CommitLineData
a65a7e70 1//#include <string.h>
2//#include <TStyle.h>
3#include <list>
4#include <string>
5
6#include "TTree.h"
7#include "TCanvas.h"
8#include "AliAnalysisTask.h"
9#include "AliInputEventHandler.h"
10#include "AliESDtrack.h"
11#include "AliAODVertex.h"
12#include "AliAODCluster.h"
13
14#include <TROOT.h>
15#include <TRandom.h>
16#include <TSystem.h>
17#include <TInterpreter.h>
18#include <TChain.h>
19#include <TFile.h>
20#include <TKey.h>
21#include <TH1F.h>
22#include <TH2F.h>
23#include <TH3F.h>
24#include <TProfile.h>
25#include <TList.h>
26#include <TLorentzVector.h>
27#include <TClonesArray.h>
28#include <TRefArray.h>
29
30#include "TDatabasePDG.h"
31#include "AliAnalysisManager.h"
32#include "AliJetFinder.h"
33#include "AliJetHeader.h"
34#include "AliJetReader.h"
35#include "AliJetReaderHeader.h"
a65a7e70 36#include "AliESDEvent.h"
37#include "AliAODEvent.h"
38#include "AliAODHandler.h"
39#include "AliAODInputHandler.h"
40#include "AliAODTrack.h"
41#include "AliAODMCParticle.h"
42#include "AliAODJet.h"
43#include "AliAODJetEventBackground.h"
44#include "AliMCParticle.h"
45#include "AliAODMCParticle.h"
46#include "AliMCEventHandler.h"
47#include "AliMCEvent.h"
48#include "AliStack.h"
49
50#include "AliAODHeader.h"
51#include "AliAODMCHeader.h"
52//#include "AliGenPythiaEventHeader.h"
53#include "AliJetKineReaderHeader.h"
54#include "AliGenCocktailEventHeader.h"
55#include "AliInputEventHandler.h"
56#include "AliGenEventHeader.h"
57#include "AliGenDPMjetEventHeader.h"
58
59#include "AliAnalysisTaskJetHadronCorrelation.h"
60#include "AliAnalysisTaskPhiCorrelations.h"
61//#include "AliAnalysisHelperJetTasks.h"
62#include "AliPWG4HighPtQAMC.h"
63
64using std::cout;
65using std::endl;
66
67ClassImp(AliAnalysisTaskJetHadronCorrelation)
68
69 //________________________________________________________________________
70 AliAnalysisTaskJetHadronCorrelation::AliAnalysisTaskJetHadronCorrelation():
71 AliAnalysisTaskSE(),
72 fUseAODInput(kFALSE),
73 fJetBranch("jets"),
74 fNonStdFile(""),
75 fAODIn(0x0),
76 fAODOut(0x0),
77 fAODExtension(0x0),
78 JFAlg("ANTIKT"),
79 Radius(0.4),
80 Filtermask(272),
81 BackM(0),
82 TrackPtcut(0.15),
83 SkipCone(0),
84 IsMC(kTRUE),
85 JetEScale(1.),
86 TrackEScale(1.),
87 fxsec(0.),
88 ftrial(1.),
89 fHistList(0x0), // Output list
90 fIfiles(0),
91 fH1Events(0x0),
92 fH1Xsec(0x0),
93 fH1Trials(0x0),
94 fH1Track_pt (0x0),
95 fH1Track_phi (0x0),
96 fH1Track_eta (0x0),
97 fH1MCTrack_pt (0x0),
98 fH1MCTrack_phi (0x0),
99 fH1MCTrack_eta (0x0),
100 fH1MCPrimTrack_pt (0x0),
101 fH1MCPrimTrack_phi (0x0),
102 fH1MCPrimTrack_eta (0x0),
103 fH1Jet_pt (0x0),
104 fH1Jet_phi (0x0),
105 fH1Jet_eta (0x0),
106 fH1leadJet_pt (0x0),
107 fH1leadJet_pt_dijet (0x0),
108 fH1subJet_pt_dijet (0x0),
109 fH1JetMC_pt (0x0),
110 fH1leadJetMC_pt (0x0),
111 fH1leadJetMC_pt_dijet (0x0),
112 fH1subJetMC_pt_dijet (0x0),
113 fH2JetsJet_dphi (0x0),
114 fH2JetsJet_deta (0x0),
115 fH2JetsJet_Aj (0x0),
116 fH2JetsJet_pt (0x0),
117 fH2JetsJetMC_dphi (0x0),
118 fH2JetsJetMC_deta (0x0),
119 fH2JetsJetMC_Aj (0x0),
120 fH2JetsJetMC_pt (0x0),
121 fH2Mult_Mtrack (0x0),
122 fH2Mult_Mlead (0x0),
123 fH2Mult_Mjet (0x0),
124 fH2Mult_Njet (0x0),
125 fH2Mult_Aj (0x0),
126 fH2Mlead_Aj (0x0),
127 fH2Jet_pt_Mlead (0x0),
128 fH2Jet_pt_Munder (0x0),
129 fH2leadJetMCptResolution (0x0),
130 fH2TrackMCptResolution (0x0),
131 fH2TrackMCptEfficiency (0x0),
132 fH2AjCorrelation_MCRec (0x0),
133 fH2MleadCorrelation_MCRec(0x0)
134{
135 for(int j=0;j<5;j++){
136 fH1ndiJ_ediv [j]=0;
137 fH1Aj [j]=0;
138 fH1Mlead [j]=0;
139 fH1leadJetMC_dphiResolution [j]=0;
140 fH1subJetMC_dphiResolution [j]=0;
141 fH1leadJetMC_Efficiency [j]=0;
142 fH1subJetMC_Efficiency [j]=0;
143 for(int k=0;k<5;k++){
144 fH1JetHadron_dphi_ediv [j][k]=0;
145 fH1JetHadron_dphi_tptweight_ediv [j][k]=0;
146 fH1JetHadron_dphi_tJptweight_ediv[j][k]=0;
147 fH1JetHadronMC_dphi_ediv [j][k]=0;
148 fH1JetHadronMC_dphi_tptweight_ediv [j][k]=0;
149 fH1JetHadronMC_dphi_tJptweight_ediv[j][k]=0;
150 fH1JetHadronMCPrim_dphi_ediv [j][k]=0;
151 fH1JetHadronMCPrim_dphi_tptweight_ediv [j][k]=0;
152 fH1JetHadronMCPrim_dphi_tJptweight_ediv[j][k]=0;
153 }
154 }
155 for(int j=0;j<3;j++){
156 fH1ndiJ_2040Mlead [j]=0;
157 fH1ndiJ_2040Aj [j]=0;
158 for(int k=0;k<5;k++){
159 fH1JetHadron_dphi_tptweight2040_Mleaddep[j][k]=0;
160 fH1JetHadron_dphi_tptweight2040_Ajdep [j][k]=0;
161 fH1JetHadronMC_dphi_tptweight2040_Mleaddep[j][k]=0;
162 fH1JetHadronMC_dphi_tptweight2040_Ajdep [j][k]=0;
163 fH1JetHadronMCPrim_dphi_tptweight2040_Mleaddep[j][k]=0;
164 fH1JetHadronMCPrim_dphi_tptweight2040_Ajdep [j][k]=0;
165 }
166 }
167 // Default constructor
168}
169
170//________________________________________________________________________
171AliAnalysisTaskJetHadronCorrelation::AliAnalysisTaskJetHadronCorrelation(const char *name):
172 AliAnalysisTaskSE(name),
173 fUseAODInput(kFALSE),
174 fJetBranch("jets"),
175 fNonStdFile(""),
176 fAODIn(0x0),
177 fAODOut(0x0),
178 fAODExtension(0x0),
179 JFAlg("ANTIKT"),
180 Radius(0.4),
181 Filtermask(272),
182 BackM(0),
183 TrackPtcut(0.15),
184 SkipCone(0),
185 IsMC(kTRUE),
186 JetEScale(1.),
187 TrackEScale(1.),
188 fxsec(0.),
189 ftrial(1.),
190 fHistList(0x0), // Output list
191 fIfiles(0),
192 fH1Events(0x0),
193 fH1Xsec(0x0),
194 fH1Trials(0x0),
195 fH1Track_pt (0x0),
196 fH1Track_phi (0x0),
197 fH1Track_eta (0x0),
198 fH1MCTrack_pt (0x0),
199 fH1MCTrack_phi (0x0),
200 fH1MCTrack_eta (0x0),
201 fH1MCPrimTrack_pt (0x0),
202 fH1MCPrimTrack_phi (0x0),
203 fH1MCPrimTrack_eta (0x0),
204 fH1Jet_pt (0x0),
205 fH1Jet_phi (0x0),
206 fH1Jet_eta (0x0),
207 fH1leadJet_pt (0x0),
208 fH1leadJet_pt_dijet (0x0),
209 fH1subJet_pt_dijet (0x0),
210 fH1JetMC_pt (0x0),
211 fH1leadJetMC_pt (0x0),
212 fH1leadJetMC_pt_dijet (0x0),
213 fH1subJetMC_pt_dijet (0x0),
214 fH2JetsJet_dphi (0x0),
215 fH2JetsJet_deta (0x0),
216 fH2JetsJet_Aj (0x0),
217 fH2JetsJet_pt (0x0),
218 fH2JetsJetMC_dphi (0x0),
219 fH2JetsJetMC_deta (0x0),
220 fH2JetsJetMC_Aj (0x0),
221 fH2JetsJetMC_pt (0x0),
222 fH2Mult_Mtrack (0x0),
223 fH2Mult_Mlead (0x0),
224 fH2Mult_Mjet (0x0),
225 fH2Mult_Njet (0x0),
226 fH2Mult_Aj (0x0),
227 fH2Mlead_Aj (0x0),
228 fH2Jet_pt_Mlead (0x0),
229 fH2Jet_pt_Munder (0x0),
230 fH2leadJetMCptResolution (0x0),
231 fH2TrackMCptResolution (0x0),
232 fH2TrackMCptEfficiency (0x0),
233 fH2AjCorrelation_MCRec (0x0),
234 fH2MleadCorrelation_MCRec(0x0)
235{
236
237 for(int j=0;j<5;j++){
238 fH1ndiJ_ediv [j]=0;
239 fH1Aj [j]=0;
240 fH1Mlead [j]=0;
241 fH1leadJetMC_dphiResolution [j]=0;
242 fH1subJetMC_dphiResolution [j]=0;
243 fH1leadJetMC_Efficiency [j]=0;
244 fH1subJetMC_Efficiency [j]=0;
245 for(int k=0;k<5;k++){
246 fH1JetHadron_dphi_ediv [j][k]=0;
247 fH1JetHadron_dphi_tptweight_ediv [j][k]=0;
248 fH1JetHadron_dphi_tJptweight_ediv[j][k]=0;
249 fH1JetHadronMC_dphi_ediv [j][k]=0;
250 fH1JetHadronMC_dphi_tptweight_ediv [j][k]=0;
251 fH1JetHadronMC_dphi_tJptweight_ediv[j][k]=0;
252 fH1JetHadronMCPrim_dphi_ediv [j][k]=0;
253 fH1JetHadronMCPrim_dphi_tptweight_ediv [j][k]=0;
254 fH1JetHadronMCPrim_dphi_tJptweight_ediv[j][k]=0;
255 }
256 }
257 for(int j=0;j<3;j++){
258 fH1ndiJ_2040Mlead [j]=0;
259 fH1ndiJ_2040Aj [j]=0;
260 for(int k=0;k<5;k++){
261 fH1JetHadron_dphi_tptweight2040_Mleaddep[j][k]=0;
262 fH1JetHadron_dphi_tptweight2040_Ajdep [j][k]=0;
263 fH1JetHadronMC_dphi_tptweight2040_Mleaddep[j][k]=0;
264 fH1JetHadronMC_dphi_tptweight2040_Ajdep [j][k]=0;
265 fH1JetHadronMCPrim_dphi_tptweight2040_Mleaddep[j][k]=0;
266 fH1JetHadronMCPrim_dphi_tptweight2040_Ajdep [j][k]=0;
267 }
268 }
269
270 // Default constructor
271 // Constructor
272 // Define input and output slots here
273 // Input slot #0 works with a TChain
274 DefineInput(0, TChain::Class());
275 // Output slot #0 id reserved by the base class for AOD
276 // Output slot #1 writes into a TH1 container
277 DefineOutput(1, TList::Class());
278}
279
280//________________________________________________________________________
281void AliAnalysisTaskJetHadronCorrelation::UserCreateOutputObjects()
282{
283 // Create histograms
284 // Called once
285
286 fHistList = new TList();fHistList->SetOwner(kTRUE); cout<<"TList is created for output "<<endl;
287 //if (!fHistList){ fHistList = new TList();fHistList->SetOwner(kTRUE); cout<<"TList is created for output "<<endl;}
288
289 Bool_t oldStatus = TH1::AddDirectoryStatus();
290 TH1::AddDirectory(kFALSE);
291 Float_t pi=TMath::Pi();
292
293
294 char *histname;
295
296 fH1Events = new TH1F ("Events" ,"Events" ,1,0,1);
297 fH1Xsec = new TProfile("Xsec" ,"Xsec" ,1,0,1);
298 fH1Trials = new TH1F ("Trials" ,"Trials" ,1,0,1);
299
300 fH1Track_pt = new TH1F("Track_pt" ,"Track_pt" ,200,0,200);
301 fH1Track_phi = new TH1F("Track_phi" ,"Track_phi" ,100,0,2*pi);
302 fH1Track_eta = new TH1F("Track_eta" ,"Track_eta" ,100,-1.,1);
303 fH1MCTrack_pt = new TH1F("MCTrack_pt" ,"MCTrack_pt" ,200,0,200);
304 fH1MCTrack_phi = new TH1F("MCTrack_phi" ,"MCTrack_phi" ,100,0,2*pi);
305 fH1MCTrack_eta = new TH1F("MCTrack_eta" ,"MCTrack_eta" ,100,-1.,1);
306 fH1MCPrimTrack_pt = new TH1F("MCPrimTrack_pt" ,"MCPrimTrack_pt" ,200,0,200);
307 fH1MCPrimTrack_phi = new TH1F("MCPrimTrack_phi" ,"MCPrimTrack_phi" ,100,0,2*pi);
308 fH1MCPrimTrack_eta = new TH1F("MCPrimTrack_eta" ,"MCPrimTrack_eta" ,100,-1.,1);
309 fH1Jet_pt = new TH1F("Jet_pt" ,"Jet_pt" ,200,0,200);
310 fH1Jet_phi = new TH1F("Jet_phi" ,"Jet_pt" ,100,0,2*pi);
311 fH1Jet_eta = new TH1F("Jet_eta" ,"Jet_pt" ,100,-1.,1);
312 fH1leadJet_pt = new TH1F("leadJet_pt" ,"leadJet_pt" ,200,0,200);
313 fH1leadJet_pt_dijet = new TH1F("leadJet_pt_dijet" ,"leadJet_pt_dijet" ,200,0,200);
314 fH1subJet_pt_dijet = new TH1F("subJet_pt_dijet" ,"subJet_pt_dijet" ,200,0,200);
315 fH1JetMC_pt = new TH1F("JetMC_pt" ,"JetMC_pt" ,200,0,200);
316 fH1leadJetMC_pt = new TH1F("leadJetMC_pt" ,"leadJetMC_pt" ,200,0,200);
317 fH1leadJetMC_pt_dijet = new TH1F("leadJetMC_pt_dijet","leadJetMC_pt_dijet",200,0,200);
318 fH1subJetMC_pt_dijet = new TH1F("subJetMC_pt_dijet" ,"subJetMC_pt_dijet" ,200,0,200);
319 fH2JetsJet_dphi = new TH2F("JetsJet_dphi" ,"JetsJet_dphi" ,200,0,200,100,-2*pi,2*pi);
320 fH2JetsJet_deta = new TH2F("JetsJet_deta" ,"JetsJet_deta" ,200,0,200,100,-1.5,1.5);
321 fH2JetsJet_Aj = new TH2F("JetsJet_Aj" ,"JetsJet_Aj" ,200,0,200,100,0,1.2);
322 fH2JetsJet_pt = new TH2F("JetsJet_pt" ,"JetsJet_pt" ,200,0,200,200,0,200);
323 fH2JetsJetMC_dphi = new TH2F("JetsJetMC_dphi" ,"JetsJetMC_dphi" ,200,0,200,100,-2*pi,2*pi);
324 fH2JetsJetMC_deta = new TH2F("JetsJetMC_deta" ,"JetsJetMC_deta" ,200,0,200,100,-1.5,1.5);
325 fH2JetsJetMC_Aj = new TH2F("JetsJetMC_Aj" ,"JetsJetMC_Aj" ,200,0,200,100,0,1.2);
326 fH2JetsJetMC_pt = new TH2F("JetsJetMC_pt" ,"JetsJetMC_pt" ,200,0,200,200,0,200);
327 fH2Mult_Mtrack = new TH2F("Mult_Mtrack" ,"Mult_Mtrack" ,50,0,250,50,0,250);
328 fH2Mult_Mlead = new TH2F("Mult_Mlead" ,"Mult_Mlead" ,50,0,250,25,0,25);
329 fH2Mult_Mjet = new TH2F("Mult_Mjet" ,"Mult_Mjet" ,50,0,250,50,0,100);
330 fH2Mult_Njet = new TH2F("Mult_Njet" ,"Mult_Njet" ,50,0,250,50,0,50);
331 fH2Mult_Aj = new TH2F("Mult_Aj" ,"Mult_Aj" ,50,0,250,25,0,1.2);
332 fH2Mlead_Aj = new TH2F("Mlead_Aj" ,"Mlead_Aj" ,25,0,25,25,0,1.2);
333 fH2Jet_pt_Mlead = new TH2F("Jet_pt_Mlead" ,"Jet_pt_Mlead" ,50,0,200,25,0,25);
334 fH2Jet_pt_Munder = new TH2F("Jet_pt_Munder" ,"Jet_pt_Munder" ,50,0,200,25,0,5);
335 fH2leadJetMCptResolution = new TH2F("leadJetMCptResolution" ,"leadJetMCptResolution" ,200,0,200,200,0,200);
336 fH2TrackMCptResolution = new TH2F("TrackMCptResolution" ,"TrackMCptResolution" ,200,0,200,200,0,200);
337 fH2TrackMCptEfficiency = new TH2F("TrackMCptEfficiency" ,"TrackMCptEfficiency" ,200,0,200,100,0,1.2);
338 fH2AjCorrelation_MCRec = new TH2F("AjCorrelation_MCRec" ,"AjCorrelation_MCRec" ,60,0,1.2,60,0,1.2);
339 fH2MleadCorrelation_MCRec = new TH2F("MleadCorrelation_MCRec","MleadCorrelation_MCRec",60,0,60,60,0,60);
340
341 for(int j=0;j<5;j++){
342 histname = Form("ndiJ_ediv%d",j);
343 fH1ndiJ_ediv[j]= new TH1F(histname,histname,1,1,2);
344 histname = Form("Aj%d",j);
345 fH1Aj[j] = new TH1F(histname,histname,50,0,1.2);
346 histname = Form("Mlead%d",j);
347 fH1Mlead[j] = new TH1F(histname,histname,50,0,50);
348 histname = Form("leadJetMC_dphiResolution%d",j);
349 fH1leadJetMC_dphiResolution[j] = new TH1F(histname,histname,200,-2*pi,2*pi);
350 histname = Form("subJetMC_dphiResolution%d",j);
351 fH1subJetMC_dphiResolution[j] = new TH1F(histname,histname,200,-2*pi,2*pi);
352 histname = Form("leadJetMC_Efficiency%d",j);
353 fH1leadJetMC_Efficiency[j] = new TH1F(histname,histname,100,0,1.2);
354 histname = Form("subJetMC_Efficiency%d",j);
355 fH1subJetMC_Efficiency[j] = new TH1F(histname,histname,100,0,1.2);
356 for(int k=0;k<5;k++){
357 histname = Form("JetHadron_dphi_ediv%d%d",j,k);
358 fH1JetHadron_dphi_ediv [j][k]= new TH1F(histname,histname,200,-1./2.*pi,3./2.*pi);
359 histname = Form("JetHadron_dphi_tptweight_ediv%d%d",j,k);
360 fH1JetHadron_dphi_tptweight_ediv [j][k]= new TH1F(histname,histname,200,-1./2.*pi,3./2.*pi);
361 histname = Form("JetHadron_dphi_tJptweight_ediv%d%d",j,k);
362 fH1JetHadron_dphi_tJptweight_ediv [j][k]= new TH1F(histname,histname,200,-1./2.*pi,3./2.*pi);
363 histname = Form("JetHadronMC_dphi_ediv%d%d",j,k);
364 fH1JetHadronMC_dphi_ediv [j][k]= new TH1F(histname,histname,200,-1./2.*pi,3./2.*pi);
365 histname = Form("JetHadronMC_dphi_tptweight_ediv%d%d",j,k);
366 fH1JetHadronMC_dphi_tptweight_ediv [j][k]= new TH1F(histname,histname,200,-1./2.*pi,3./2.*pi);
367 histname = Form("JetHadronMC_dphi_tJptweight_ediv%d%d",j,k);
368 fH1JetHadronMC_dphi_tJptweight_ediv [j][k]= new TH1F(histname,histname,200,-1./2.*pi,3./2.*pi);
369 histname = Form("JetHadronMCPrim_dphi_ediv%d%d",j,k);
370 fH1JetHadronMCPrim_dphi_ediv [j][k]= new TH1F(histname,histname,200,-1./2.*pi,3./2.*pi);
371 histname = Form("JetHadronMCPrim_dphi_tptweight_ediv%d%d",j,k);
372 fH1JetHadronMCPrim_dphi_tptweight_ediv [j][k]= new TH1F(histname,histname,200,-1./2.*pi,3./2.*pi);
373 histname = Form("JetHadronMCPrim_dphi_tJptweight_ediv%d%d",j,k);
374 fH1JetHadronMCPrim_dphi_tJptweight_ediv [j][k]= new TH1F(histname,histname,200,-1./2.*pi,3./2.*pi);
375 }
376 }
377 for(int j=0;j<3;j++){
378 histname = Form("ndiJ_2040Mlead%d",j);
379 fH1ndiJ_2040Mlead[j]= new TH1F(histname,histname,1,1,2);
380 histname = Form("ndiJ_2040Aj%d",j);
381 fH1ndiJ_2040Aj[j]= new TH1F(histname,histname,1,1,2);
382 for(int k=0;k<5;k++){
383 histname = Form("JetHadron_dphi_tptweight2040_Mleaddep%d%d",j,k);
384 fH1JetHadron_dphi_tptweight2040_Mleaddep [j][k]= new TH1F(histname,histname,200,-1./2.*pi,3./2.*pi);
385 histname = Form("JetHadron_dphi_tptweight2040_Ajdep%d%d",j,k);
386 fH1JetHadron_dphi_tptweight2040_Ajdep [j][k]= new TH1F(histname,histname,200,-1./2.*pi,3./2.*pi);
387 histname = Form("JetHadronMC_dphi_tptweight2040_Mleaddep%d%d",j,k);
388 fH1JetHadronMC_dphi_tptweight2040_Mleaddep [j][k]= new TH1F(histname,histname,200,-1./2.*pi,3./2.*pi);
389 histname = Form("JetHadronMC_dphi_tptweight2040_Ajdep%d%d",j,k);
390 fH1JetHadronMC_dphi_tptweight2040_Ajdep [j][k]= new TH1F(histname,histname,200,-1./2.*pi,3./2.*pi);
391 histname = Form("JetHadronMCPrim_dphi_tptweight2040_Mleaddep%d%d",j,k);
392 fH1JetHadronMCPrim_dphi_tptweight2040_Mleaddep [j][k]= new TH1F(histname,histname,200,-1./2.*pi,3./2.*pi);
393 histname = Form("JetHadronMCPrim_dphi_tptweight2040_Ajdep%d%d",j,k);
394 fH1JetHadronMCPrim_dphi_tptweight2040_Ajdep [j][k]= new TH1F(histname,histname,200,-1./2.*pi,3./2.*pi);
395 }
396 }
397
398
399 if(IsMC){
400 fHistList->Add(fH1Events );
401 fHistList->Add(fH1Xsec );
402 fHistList->Add(fH1Trials );
403 fHistList->Add(fH1Track_pt );
404 fHistList->Add(fH1Track_phi );
405 fHistList->Add(fH1Track_eta );
406 fHistList->Add(fH1MCTrack_pt );
407 fHistList->Add(fH1MCTrack_phi );
408 fHistList->Add(fH1MCTrack_eta );
409 fHistList->Add(fH1MCPrimTrack_pt );
410 fHistList->Add(fH1MCPrimTrack_phi );
411 fHistList->Add(fH1MCPrimTrack_eta );
412 fHistList->Add(fH1Jet_pt );
413 fHistList->Add(fH1Jet_phi );
414 fHistList->Add(fH1Jet_eta );
415 fHistList->Add(fH1leadJet_pt );
416 fHistList->Add(fH1leadJet_pt_dijet );
417 fHistList->Add(fH1subJet_pt_dijet );
418 fHistList->Add(fH1JetMC_pt );
419 fHistList->Add(fH1leadJetMC_pt );
420 fHistList->Add(fH1leadJetMC_pt_dijet);
421 fHistList->Add(fH1subJetMC_pt_dijet );
422 fHistList->Add(fH2JetsJet_dphi );
423 fHistList->Add(fH2JetsJet_deta );
424 fHistList->Add(fH2JetsJet_Aj );
425 fHistList->Add(fH2JetsJet_pt );
426 fHistList->Add(fH2JetsJetMC_dphi );
427 fHistList->Add(fH2JetsJetMC_deta );
428 fHistList->Add(fH2JetsJetMC_Aj );
429 fHistList->Add(fH2JetsJetMC_pt );
430 fHistList->Add(fH2Mult_Mtrack );
431 fHistList->Add(fH2Mult_Mlead );
432 fHistList->Add(fH2Mult_Mjet );
433 fHistList->Add(fH2Mult_Njet );
434 fHistList->Add(fH2Mult_Aj );
435 fHistList->Add(fH2Mlead_Aj );
436 fHistList->Add(fH2Jet_pt_Mlead );
437 fHistList->Add(fH2Jet_pt_Munder );
438 fHistList->Add(fH2leadJetMCptResolution );
439 fHistList->Add(fH2TrackMCptResolution );
440 fHistList->Add(fH2TrackMCptEfficiency );
441 fHistList->Add(fH2AjCorrelation_MCRec );
442 fHistList->Add(fH2MleadCorrelation_MCRec);
443 for(int j=0;j<5;j++){
444 fHistList->Add(fH1ndiJ_ediv [j]);
445 fHistList->Add(fH1Aj [j]);
446 fHistList->Add(fH1Mlead [j]);
447 fHistList->Add(fH1leadJetMC_dphiResolution [j]);
448 fHistList->Add(fH1subJetMC_dphiResolution [j]);
449 fHistList->Add(fH1leadJetMC_Efficiency [j]);
450 fHistList->Add(fH1subJetMC_Efficiency [j]);
451 for(int k=0;k<5;k++){
452 fHistList->Add(fH1JetHadron_dphi_ediv [j][k]);
453 fHistList->Add(fH1JetHadron_dphi_tptweight_ediv [j][k]);
454 fHistList->Add(fH1JetHadron_dphi_tJptweight_ediv [j][k]);
455 fHistList->Add(fH1JetHadronMC_dphi_ediv [j][k]);
456 fHistList->Add(fH1JetHadronMC_dphi_tptweight_ediv [j][k]);
457 fHistList->Add(fH1JetHadronMC_dphi_tJptweight_ediv [j][k]);
458 fHistList->Add(fH1JetHadronMCPrim_dphi_ediv [j][k]);
459 fHistList->Add(fH1JetHadronMCPrim_dphi_tptweight_ediv [j][k]);
460 fHistList->Add(fH1JetHadronMCPrim_dphi_tJptweight_ediv [j][k]);
461 }
462 }
463 for(int j=0;j<3;j++){
464 fHistList->Add(fH1ndiJ_2040Mlead [j]);
465 fHistList->Add(fH1ndiJ_2040Aj [j]);
466 for(int k=0;k<5;k++){
467 fHistList->Add(fH1JetHadron_dphi_tptweight2040_Mleaddep [j][k]);
468 fHistList->Add(fH1JetHadron_dphi_tptweight2040_Ajdep [j][k]);
469 fHistList->Add(fH1JetHadronMC_dphi_tptweight2040_Mleaddep [j][k]);
470 fHistList->Add(fH1JetHadronMC_dphi_tptweight2040_Ajdep [j][k]);
471 fHistList->Add(fH1JetHadronMCPrim_dphi_tptweight2040_Mleaddep [j][k]);
472 fHistList->Add(fH1JetHadronMCPrim_dphi_tptweight2040_Ajdep [j][k]);
473 }
474 }
475 }
476 else{
477 fHistList->Add(fH1Events );
478 fHistList->Add(fH1Track_pt );
479 fHistList->Add(fH1Track_phi );
480 fHistList->Add(fH1Track_eta );
481 fHistList->Add(fH1Jet_pt );
482 fHistList->Add(fH1Jet_phi );
483 fHistList->Add(fH1Jet_eta );
484 fHistList->Add(fH1leadJet_pt );
485 fHistList->Add(fH1leadJet_pt_dijet );
486 fHistList->Add(fH1subJet_pt_dijet );
487 fHistList->Add(fH2JetsJet_dphi );
488 fHistList->Add(fH2JetsJet_deta );
489 fHistList->Add(fH2JetsJet_Aj );
490 fHistList->Add(fH2JetsJet_pt );
491 fHistList->Add(fH2Mult_Mtrack );
492 fHistList->Add(fH2Mult_Mlead );
493 fHistList->Add(fH2Mult_Mjet );
494 fHistList->Add(fH2Mult_Njet );
495 fHistList->Add(fH2Mult_Aj );
496 fHistList->Add(fH2Mlead_Aj );
497 fHistList->Add(fH2Jet_pt_Mlead );
498 fHistList->Add(fH2Jet_pt_Munder );
499 for(int j=0;j<5;j++){
500 fHistList->Add(fH1ndiJ_ediv [j]);
501 fHistList->Add(fH1Aj [j]);
502 fHistList->Add(fH1Mlead [j]);
503 for(int k=0;k<5;k++){
504 fHistList->Add(fH1JetHadron_dphi_ediv [j][k]);
505 fHistList->Add(fH1JetHadron_dphi_tptweight_ediv [j][k]);
506 fHistList->Add(fH1JetHadron_dphi_tJptweight_ediv [j][k]);
507 }
508 }
509 for(int j=0;j<3;j++){
510 fHistList->Add(fH1ndiJ_2040Mlead [j]);
511 fHistList->Add(fH1ndiJ_2040Aj [j]);
512 for(int k=0;k<5;k++){
513 fHistList->Add(fH1JetHadron_dphi_tptweight2040_Mleaddep [j][k]);
514 fHistList->Add(fH1JetHadron_dphi_tptweight2040_Ajdep [j][k]);
515 }
516 }
517 }
518
519 // =========== Switch on Sumw2 for all histos ===========
520 for (Int_t i=0; i<fHistList->GetEntries(); ++i)
521 {
522 TH1 *h1 = dynamic_cast<TH1*>(fHistList->At(i));
523 if (h1)
524 {
525 h1->Sumw2();
526 continue;
527 }
528 THnSparse *hn = dynamic_cast<THnSparse*>(fHistList->At(i));
529 if(hn)hn->Sumw2();
530 }
531 TH1::AddDirectory(oldStatus);
532
533 PostData(1,fHistList);
534}
535
536
537//----------------------------------------------------------------------
538void AliAnalysisTaskJetHadronCorrelation::Init()
539{
540 // Initialization
541 if (fDebug) printf("AnalysisTaskJetHadronCorrelation::Init() \n");
542
543}
544
545Bool_t AliAnalysisTaskJetHadronCorrelation::Notify()
546{
547
548
549 fAODIn = dynamic_cast<AliAODEvent*>(InputEvent());
550 fAODOut = AODEvent();
551 if(fNonStdFile.Length()!=0){
552 AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
553 fAODExtension = (aodH?aodH->GetExtension(fNonStdFile.Data()):0);
554 if(fAODExtension){
555 if(fDebug>1)Printf("AODExtension found for %s ",fNonStdFile.Data());
556 }
557 }
558
559 TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
560 fxsec=0;
561 ftrial=1;
562
563 if(tree){
564 TFile *curfile = tree->GetCurrentFile();
565 if(!curfile){
566 Error("Notify","No current file");
567 return kFALSE;
568 }
569 if(IsMC){
570 AliPWG4HighPtQAMC::PythiaInfoFromFile(curfile->GetName(),fxsec,ftrial);
571 fH1Xsec ->Fill(0.,fxsec);
572 fH1Trials->Fill(0.,ftrial);
573 }
574
575 }
576
577 printf("Reading File %s ",fInputHandler->GetTree()->GetCurrentFile()->GetName());
578 return kTRUE;
579}
580void AliAnalysisTaskJetHadronCorrelation::FinishTaskOutput()
581{
582}
583
584
585
586//________________________________________________________________________
587void AliAnalysisTaskJetHadronCorrelation::UserExec(Option_t *)
588{
589
590
591 // Main loop (called each event)
592 // Execute analysis for current event
593
594 AliAODEvent *fAOD;
595 TObject* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
596 if( handler && handler->InheritsFrom("AliAODInputHandler") ) {
597 fAOD = ((AliAODInputHandler*)handler)->GetEvent();
598 if(fUseAODInput) fAODIn = fAOD;
599 if (fDebug > 1) Printf("%s:%d AOD event from input", (char*)__FILE__,__LINE__);
600 }
601 else {
602 handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
603 if( handler && handler->InheritsFrom("AliAODHandler") ) {
604 fAOD = ((AliAODHandler*)handler)->GetAOD();
605 fAODIn = fAOD;
606 if (fDebug > 1) Printf("%s:%d AOD event from output", (char*)__FILE__,__LINE__);
607 }
608 }
609
610 if(!fAODIn && !fUseAODInput){ // case we have AOD in input & output and want jets from output
611 TObject* outHandler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
612 if( outHandler && outHandler->InheritsFrom("AliAODHandler") ) {
613 fAODIn = ((AliAODHandler*)outHandler)->GetAOD();
614 if (fDebug > 1) Printf("%s:%d jets from output AOD", (char*)__FILE__,__LINE__);
615 }
616 }
617 if (!fAODIn) {
618 Printf("ERROR %s : fAODIn not available",(char*)__FILE__);
619 return;
620 }
621
622 AliInputEventHandler* inputHandler = (AliInputEventHandler*)
623 ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
624 if(!(inputHandler->IsEventSelected() & AliVEvent::kMB)){
625 if (fDebug > 1 ) Printf(" Trigger Selection: event REJECTED ... ");
626 return;
627 }
628 fH1Events->Fill(0);
629
630 AliAODHeader* aliH = dynamic_cast <AliAODHeader*> (fAODIn->GetHeader());
631 if(!aliH){
632 Printf("ERROR: AliAODHeader not available");
633 return;
634 }
635 double Mult = aliH->GetRefMultiplicity();
636
637 // start jet analysis -------------------------Init.
638 Float_t pi=TMath::Pi();
639
640 Double_t Jet_pt [20][5000];
641 Double_t Jet_phi [20][5000];
642 Double_t Jet_eta [20][5000];
643 Double_t Jet_area [20][5000];
644 Double_t subJet_pt [20][5000];
645 Double_t subJet_eta[20][5000];
646 Double_t Track_n ;
647 Double_t Track_pt [5000];
648 Double_t Track_eta[5000];
649 Double_t Track_phi[5000];
650 Double_t MCTrack_n ;
651 Double_t MCTrack_pt [5000];
652 Double_t MCTrack_eta[5000];
653 Double_t MCTrack_phi[5000];
654
655 Track_n=0;MCTrack_n=0;
656 for(int i=0;i<20;i++){
657 for(int j=0;j<1000;j++){
658 Jet_pt[i][j]=0.;
659 Jet_phi[i][j]=999.;
660 Jet_eta[i][j]=999.;
661 Jet_area[i][j]=999.;
662 subJet_pt[i][j]=0.;
663 subJet_eta[i][j]=999.;
664 Track_pt [j]=0.;
665 Track_phi[j]=999.;
666 Track_eta[j]=999.;
667 MCTrack_pt [j]=0.;
668 MCTrack_phi[j]=999.;
669 MCTrack_eta[j]=999.;
670 }
671 }
672
673 int nLJetAOD=999; double ptLJetAOD=0;double phiLJetAOD=999.;double etaLJetAOD=999.;double ptsLJetAOD=0;double phisLJetAOD=900.;double etasLJetAOD=900.;
674 int nLJetMC2=999; double ptLJetMC2=0;double phiLJetMC2=999.;double etaLJetMC2=999.;double ptsLJetMC2=0;double phisLJetMC2=900.;double etasLJetMC2=900.;
675 int nLJetMC =999; double ptLJetMC =0;double phiLJetMC =999.;double etaLJetMC =999.;double ptsLJetMC =0;double phisLJetMC =900.;double etasLJetMC =900.;
676 bool findLJetAOD=false;
677 bool findLJetMC2=false;
678 bool findDiJet=false,findDiJetMC=false;
679 int nLJet = 999;
680 int Mjet_tot =0;
681 int Njet_tot =0;
682
683 double Aj=99.,AjMC=99.;
684 double Mlead=99.,MleadMC=99.;
685 int Munder=99.;
686
687 ////--------------------------------------------------------------------Init.
688
689 TString cAdd = "";
690 TString Branchname_gen,Branchname_gen2,Branchname_rec;
691 if((JFAlg=="ANTIKT")||(JFAlg=="KT")){
692 cAdd += Form("%02d_",(int)((Radius+0.01)*10.));
693 cAdd += Form("B%d",(int)BackM);
694 cAdd += Form("_Filter%05d",Filtermask);
695 cAdd += Form("_Cut%05d",(int)(1000.*TrackPtcut));
696 cAdd += Form("_Skip%02d",SkipCone);
697 Branchname_gen = Form("clustersAODMC_%s%s",JFAlg.Data(),cAdd.Data());
698 Branchname_gen2 = Form("clustersAODMC2_%s%s",JFAlg.Data(),cAdd.Data());
699 Branchname_rec = Form("clustersAOD_%s%s",JFAlg.Data(),cAdd.Data());
700 }
701 else{
702 cAdd += Form("%02d_",(int)((Radius+0.01)*10.));
703 cAdd += Form("B%d",(int)BackM);
704 cAdd += Form("_Filter%05d",Filtermask);
705 cAdd += Form("_Cut%05d",(int)(1000.*TrackPtcut));
706 Branchname_gen = Form("jetsAODMC_%s%s",JFAlg.Data(),cAdd.Data());
707 Branchname_gen2 = Form("jetsAODMC2_%s%s",JFAlg.Data(),cAdd.Data());
708 Branchname_rec = Form("jetsAOD_%s%s",JFAlg.Data(),cAdd.Data());
709 }
710
711
712
713 //count number of tracks@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
714 //Reconstructed Track
715 TClonesArray* tracks = dynamic_cast <TClonesArray*> (fAODIn->GetTracks());
716 if(!tracks){
717 if (fDebug > 1) Printf("%s:%d could not get AODtracks", (char*)__FILE__,__LINE__);
718 return;
719 }
720
721 Bool_t TrackEff[5000];
722 for(int i=0;i<5000;i++){
723 TrackEff[i]=false;
724 }
725 Int_t nt = fAODIn->GetNumberOfTracks();
726 AliAODTrack* trackAOD=NULL;
727 for(int ntrack =0;ntrack<nt;ntrack++){
728 trackAOD = (AliAODTrack*) (tracks->At(ntrack));
729 Bool_t bgoodT=false;
730 if(Filtermask!=272){if(trackAOD->TestFilterMask(Filtermask))bgoodT=true;}
731 else {if(trackAOD->IsHybridGlobalConstrainedGlobal())bgoodT=true;} //for hybrid Track cuts
732 if(!bgoodT)continue;
733 if(TMath::Abs(trackAOD->Eta())<0.9){
734 Track_n++;
735 fH1Track_pt ->Fill(trackAOD->Pt()*TrackEScale);
736 fH1Track_phi->Fill(trackAOD->Phi());
737 fH1Track_eta->Fill(trackAOD->Eta());
738 //cout<<"Scale "<<TrackEScale<<" org pt "<<trackAOD->Pt()<< " scaled pt "<< trackAOD->Pt()*TrackEScale <<endl;
739 if(IsMC){
740 // track pt resplution-------------------
741 Int_t MCID = TMath::Abs(trackAOD->GetLabel());
742 TClonesArray* mctracks = dynamic_cast <TClonesArray*> (fAODIn->GetList()->FindObject(AliAODMCParticle::StdBranchName()));
743 if(!mctracks){
744 if (fDebug > 1) Printf("%s:%d could not get AODMCtracks", (char*)__FILE__,__LINE__);
745 continue;
746 }
747 AliAODMCParticle *trackMCAOD = (AliAODMCParticle*) mctracks->At(MCID);
748 fH2TrackMCptResolution->Fill(trackMCAOD->Pt(),trackAOD->Pt());
749 TrackEff[MCID]=true;
750 // --------------------------------------
751 }
752 }
753 }
754 if(IsMC){
755 //MC Track
756 TClonesArray* mctracks = dynamic_cast <TClonesArray*> (fAODIn->GetList()->FindObject(AliAODMCParticle::StdBranchName()));
757 if(!mctracks){
758 if (fDebug > 1) Printf("%s:%d could not get AODMCtracks", (char*)__FILE__,__LINE__);
759 return;
760 }
761 Int_t ntmc = mctracks->GetEntriesFast();
762 AliAODMCParticle* trackMCAOD;
763 int lastprim=0;
764 for(int ntrack =0;ntrack<ntmc;ntrack++){
765 trackMCAOD = (AliAODMCParticle*) (mctracks->At(ntrack));
766 if((trackMCAOD->IsPhysicalPrimary())==1)lastprim=ntrack;
767 }
768 for(int ntrack =0;ntrack<ntmc;ntrack++){
769 trackMCAOD = (AliAODMCParticle*) (mctracks->At(ntrack));
770 if((trackMCAOD->GetPdgCode()>10)&&((trackMCAOD->GetMother())>1)&&(ntrack>lastprim)&&(trackMCAOD->Charge())){// for Decay particles
771 if(TMath::Abs(trackMCAOD->Eta())<0.9){
772 fH1MCTrack_pt ->Fill(trackMCAOD->Pt());
773 fH1MCTrack_phi->Fill(trackMCAOD->Phi());
774 fH1MCTrack_eta->Fill(trackMCAOD->Eta());
775 if(TrackEff[ntrack])fH2TrackMCptEfficiency->Fill(trackMCAOD->Pt(),1);
776 else fH2TrackMCptEfficiency->Fill(trackMCAOD->Pt(),0);
777 }
778 }
779 if((trackMCAOD->IsPhysicalPrimary())&&(trackMCAOD->Charge())){// for Physical particles
780 if(TMath::Abs(trackMCAOD->Eta())<0.9){
781 MCTrack_n++;
782 fH1MCTrack_pt ->Fill(trackMCAOD->Pt());
783 fH1MCTrack_phi->Fill(trackMCAOD->Phi());
784 fH1MCTrack_eta->Fill(trackMCAOD->Eta());
785 fH1MCPrimTrack_pt ->Fill(trackMCAOD->Pt());
786 fH1MCPrimTrack_phi->Fill(trackMCAOD->Phi());
787 fH1MCPrimTrack_eta->Fill(trackMCAOD->Eta());
788 if(TrackEff[ntrack])fH2TrackMCptEfficiency->Fill(trackMCAOD->Pt(),1);
789 else fH2TrackMCptEfficiency->Fill(trackMCAOD->Pt(),0);
790 }
791 }
792 }
793 }
794 //@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ count number of tracks
795
796
797
798
799 for(int algorithm=0;algorithm<3;algorithm++){
800
801 if(algorithm==0)fJetBranch = Branchname_rec.Data();
802 if(algorithm==1)fJetBranch = Branchname_gen2.Data();
803 if(algorithm==2)fJetBranch = Branchname_gen.Data();
804
805 if((!IsMC&&(algorithm==1||algorithm==2)))continue;
806
807 TClonesArray* jets = dynamic_cast <TClonesArray*> (fAODIn->FindListObject(fJetBranch.Data()));
808 if(!jets){
809 printf(" Tere are no Branch named %s \n",fJetBranch.Data());
810 continue;
811 }
812 Int_t nj = jets->GetEntriesFast();
813 if (fDebug) printf("There are %5d jets in the event \n", nj);
814 AliAODJet* jetsAOD;
815 //Find Leading Jet -------------------------------------------------------
816 for(int njet =0;njet<nj;njet++){
817 jetsAOD = (AliAODJet*) (jets->At(njet));
818 Jet_pt [algorithm][njet] = jetsAOD->Pt()*JetEScale;
819 Jet_phi [algorithm][njet] = jetsAOD->Phi();
820 Jet_eta [algorithm][njet] = jetsAOD->Eta();
821 Jet_area [algorithm][njet] = jetsAOD->EffectiveAreaCharged();
822
823
824 TRefArray *reftracks = jetsAOD->GetRefTracks();
825 if(algorithm==0){if(Jet_pt[algorithm][njet]>1.)Mjet_tot += reftracks->GetEntriesFast();Njet_tot++;}
826 double eta_cut_Jet=0.5;
827 if((TMath::Abs(Jet_eta[algorithm][njet])<eta_cut_Jet)&&(Jet_pt[algorithm][njet]>10.)){
828 if(algorithm==0){
829 fH1Jet_pt ->Fill(Jet_pt [algorithm][njet]);
830 fH1Jet_phi->Fill(Jet_phi[algorithm][njet]);
831 fH1Jet_eta->Fill(Jet_eta[algorithm][njet]);
832 if(Jet_pt[algorithm][njet]>ptLJetAOD){
833 findLJetAOD=true;
834 nLJetAOD=njet;ptLJetAOD=Jet_pt[algorithm][njet];phiLJetAOD=Jet_phi[algorithm][njet];etaLJetAOD=Jet_eta[algorithm][njet];
835 }
836 }
837 if(algorithm==1){
838 fH1JetMC_pt->Fill(Jet_pt[algorithm][njet]);
839 if(Jet_pt[algorithm][njet]>ptLJetMC2){
840 findLJetMC2=true;
841 nLJetMC2=njet;ptLJetMC2=Jet_pt[algorithm][njet];phiLJetMC2=Jet_phi[algorithm][njet];etaLJetMC2=Jet_eta[algorithm][njet];
842 }
843 }
844 if(algorithm==2){
845 if(Jet_pt[algorithm][njet]>ptLJetMC){
846 nLJetMC=njet;ptLJetMC=Jet_pt[algorithm][njet];phiLJetMC=Jet_phi[algorithm][njet];etaLJetMC=Jet_eta[algorithm][njet];
847 }
848 }
849 }
850 }//njet loop
851 if(algorithm==0){nLJet=nLJetAOD;fH1leadJet_pt ->Fill(Jet_pt[algorithm][nLJet]);}
852 if(algorithm==1){nLJet=nLJetMC2;fH1leadJetMC_pt->Fill(Jet_pt[algorithm][nLJet]);}
853 if(algorithm==2){nLJet=nLJetMC;}
854 if(findLJetAOD&&(algorithm==0)){
855 jetsAOD = (AliAODJet*) (jets->At(nLJet));
856 TRefArray *reftracks = jetsAOD->GetRefTracks();
857 Mlead = reftracks->GetEntriesFast();
858 }
859 if(findLJetMC2&&(algorithm==1)){
860 jetsAOD = (AliAODJet*) (jets->At(nLJetMC2));
861 TRefArray *reftracks = jetsAOD->GetRefTracks();
862 MleadMC = reftracks->GetEntriesFast();
863 }
864 //----------------------------------------------------------- Leading Jet
865 if(nj<2)continue;
866 //Find Sub leading Jet ==================================================
867 for(int njet=0;njet<nj;njet++){
868 if(njet==nLJet)continue;
869 jetsAOD = (AliAODJet *)jets->At(njet);
870 subJet_pt [algorithm][njet] = jetsAOD->Pt()*JetEScale;
871 subJet_eta[algorithm][njet] = jetsAOD->Eta();
872 double eta_cut_Jet=0.5;
873 if((TMath::Abs(subJet_eta[algorithm][njet])<eta_cut_Jet) && (subJet_pt[algorithm][njet]>10.)){
874 if(subJet_pt[algorithm][njet]>ptsLJetAOD&&algorithm==0){
875 ptsLJetAOD=Jet_pt[algorithm][njet];phisLJetAOD=Jet_phi[algorithm][njet];etasLJetAOD=Jet_eta[algorithm][njet];
876 }
877 if(subJet_pt[algorithm][njet]>ptsLJetMC2 &&algorithm==1){
878 ptsLJetMC2=Jet_pt[algorithm][njet];phisLJetMC2=Jet_phi[algorithm][njet];etasLJetMC2=Jet_eta[algorithm][njet];
879 }
880 if(subJet_pt[algorithm][njet]>ptsLJetMC &&algorithm==2){
881 ptsLJetMC =Jet_pt[algorithm][njet];phisLJetMC =Jet_phi[algorithm][njet];etasLJetMC =Jet_eta[algorithm][njet];
882 }
883 }
884 }
885 //====================================================== Sub leading Jet
886
887 double Leading_pt=0.;double Leading_phi=999.;double Leading_eta=999.;double sLeading_pt=0.;double sLeading_phi=999.;double sLeading_eta=999.;
888 if(algorithm==0){Leading_pt=ptLJetAOD;Leading_phi=phiLJetAOD;Leading_eta=etaLJetAOD;sLeading_pt=ptsLJetAOD;sLeading_phi=phisLJetAOD;sLeading_eta=etasLJetAOD;}
889 if(algorithm==1){Leading_pt=ptLJetMC2;Leading_phi=phiLJetMC2;Leading_eta=etaLJetMC2;sLeading_pt=ptsLJetMC2;sLeading_phi=phisLJetMC2;sLeading_eta=etasLJetMC2;}
890 if(algorithm==2){Leading_pt=ptLJetMC ;Leading_phi=phiLJetMC ;Leading_eta=etaLJetMC ;sLeading_pt=ptsLJetMC ;sLeading_phi=phisLJetMC ;sLeading_eta=etasLJetMC ;}
891
892 ////Di-Jet event trigger +++++++++++++++++++++++++++++++++++++++++++++++
893 double DPhi = DeltaPhi(Leading_phi,sLeading_phi);
894 double DEta = Leading_eta-sLeading_eta;
895 if(algorithm==0){
896 fH2JetsJet_dphi->Fill(Leading_pt,DPhi);
897 fH2JetsJet_deta->Fill(Leading_pt,DEta);
898 }
899 if(algorithm==1){
900 fH2JetsJetMC_dphi->Fill(Leading_pt,DPhi);
901 fH2JetsJetMC_deta->Fill(Leading_pt,DEta);
902 }
903 if((TMath::Cos(DPhi)<-0.5)&&(Leading_pt>10.)&&(sLeading_pt>10.)){
904 if(algorithm==0)Aj = (Leading_pt-sLeading_pt)/(Leading_pt+sLeading_pt);
905 if(algorithm==1)AjMC = (Leading_pt-sLeading_pt)/(Leading_pt+sLeading_pt);
906 if(algorithm==0){
907 fH1subJet_pt_dijet ->Fill(sLeading_pt);
908 fH1leadJet_pt_dijet->Fill(Leading_pt);
909 fH2JetsJet_Aj ->Fill(Leading_pt,Aj);
910 fH2JetsJet_pt ->Fill(Leading_pt,sLeading_pt);
911 fH2Mult_Aj ->Fill(Mult,Aj);
912 fH2Mlead_Aj ->Fill(Mlead,Aj);
913 for(int eb=0;eb<5;eb++){
914 if(TMath::Abs(Leading_pt -10.-20.*(eb))<10.){
915 fH1Aj[eb] ->Fill(Aj);
916 }
917 }
918 }
919 if(algorithm==1){
920 fH1leadJetMC_pt_dijet->Fill(Leading_pt);
921 fH1subJetMC_pt_dijet ->Fill(sLeading_pt);
922 fH2JetsJetMC_Aj ->Fill(Leading_pt,AjMC);
923 fH2JetsJetMC_pt ->Fill(Leading_pt,sLeading_pt);
924 findDiJetMC=true;
925 }
926 findDiJet=true;
927
928 }
929 ////+++++++++++++++++++++++++++++++++++++++++++++++ Di-Jet event trigger
930
931 if(algorithm!=0)continue;// for only data & reconstructed Jets
932
933
934 //Jet-Hadron Correlation###############################################
935 if((findDiJet)&&(Leading_pt>10.)&&(sLeading_pt>10.)){
936 double eta_cut_Jet=0.5;
937 if(TMath::Abs(Leading_eta)<eta_cut_Jet){
938 for(int eb=0;eb<5;eb++){
939 if(TMath::Abs(Leading_pt -10.-20.*(eb))<10.){
940 fH1ndiJ_ediv[eb]->Fill(1);
941 if(eb==1){
942 if((0<Mlead)&&Mlead<7) {fH1ndiJ_2040Mlead[0]->Fill(1);}
943 else if((7<=Mlead)&&(Mlead<10)) {fH1ndiJ_2040Mlead[1]->Fill(1);}
944 else {fH1ndiJ_2040Mlead[2]->Fill(1);}
945 if((0<Aj)&&(Aj<0.19)) {fH1ndiJ_2040Aj [0]->Fill(1);}
946 else if((0.19<=Aj)&&(Aj<0.38)) {fH1ndiJ_2040Aj [1]->Fill(1);}
947 else {fH1ndiJ_2040Aj [2]->Fill(1);}
948 }
949 fH1Mlead[eb]->Fill(Mlead);
950 for(int ntrack =0;ntrack<nt;ntrack++){
951 trackAOD = (AliAODTrack*) (fAODIn->GetTrack(ntrack));
952 Bool_t bgoodT=false;
953 if(Filtermask!=272){if(trackAOD->TestFilterMask(Filtermask))bgoodT=true;}
954 else{ if(trackAOD->IsHybridGlobalConstrainedGlobal())bgoodT=true;} //for hybrid Track cuts
955 if(!bgoodT)continue;
956 Track_pt [ntrack] = (trackAOD->Pt()*TrackEScale);
957 Track_phi [ntrack] = trackAOD->Phi();
958 Track_eta [ntrack] = trackAOD->Eta();
959
960 //cout<<"Scale "<<TrackEScale<<" org pt "<<trackAOD->Pt()<< " scaled pt "<< trackAOD->Pt()*TrackEScale <<endl;
961
962 double DelPhi = DeltaPhi(Leading_phi,Track_phi[ntrack]);
963 if(TMath::Abs(Track_eta[ntrack])<0.9){
964 if((TMath::Abs(DelPhi-pi/2.)<pi/8.)||((DelPhi+pi/2.)<pi/8.)||(TMath::Abs(DelPhi-3./2.*pi)<pi/8.))Munder++;
965 for(int teb=0;teb<5;teb++){
966 if(teb==0){if(!( Track_pt[ntrack]>0.15))continue;}
967 if(teb==1){if(!((Track_pt[ntrack]<1.5)&&(Track_pt[ntrack]>0.15)))continue;}
968 if(teb==2){if(!((Track_pt[ntrack]<3.0)&&(Track_pt[ntrack]>1.5)))continue;}
969 if(teb==3){if(!((Track_pt[ntrack]<4.5)&&(Track_pt[ntrack]>3.0)))continue;}
970 if(teb==4){if(!( Track_pt[ntrack]>4.5))continue;}
971 fH1JetHadron_dphi_ediv [eb][teb]->Fill(DelPhi);
972 fH1JetHadron_dphi_tptweight_ediv [eb][teb]->Fill(DelPhi,Track_pt[ntrack]);
973 fH1JetHadron_dphi_tJptweight_ediv [eb][teb]->Fill(DelPhi,Track_pt[ntrack]/Leading_pt);
974 if(eb==1){
975 if((0<Mlead)&&Mlead<7) {fH1JetHadron_dphi_tptweight2040_Mleaddep[0][teb]->Fill(DelPhi,Track_pt[ntrack]);}
976 else if((7<=Mlead)&&(Mlead<10)){fH1JetHadron_dphi_tptweight2040_Mleaddep[1][teb]->Fill(DelPhi,Track_pt[ntrack]);}
977 else {fH1JetHadron_dphi_tptweight2040_Mleaddep[2][teb]->Fill(DelPhi,Track_pt[ntrack]);}
978 if((0<Aj)&&(Aj<0.19)) {fH1JetHadron_dphi_tptweight2040_Ajdep [0][teb]->Fill(DelPhi,Track_pt[ntrack]);}
979 else if((0.19<=Aj)&&(Aj<0.38)) {fH1JetHadron_dphi_tptweight2040_Ajdep [1][teb]->Fill(DelPhi,Track_pt[ntrack]);}
980 else {fH1JetHadron_dphi_tptweight2040_Ajdep [2][teb]->Fill(DelPhi,Track_pt[ntrack]);}
981 }
982 }
983 }
984 }//Track Loop
985 if(IsMC){
986 //MC Track
987 TClonesArray* mctracks = dynamic_cast <TClonesArray*> (fAODIn->GetList()->FindObject(AliAODMCParticle::StdBranchName()));
988 if(!mctracks){
989 if (fDebug > 1) Printf("%s:%d could not get AODMCtracks", (char*)__FILE__,__LINE__);
990 continue;
991 }
992 Int_t ntmc = mctracks->GetEntriesFast();
993 AliAODMCParticle* trackMCAOD;
994 int lastprim=0;
995 for(int ntrack =0;ntrack<ntmc;ntrack++){
996 trackMCAOD = (AliAODMCParticle*) (mctracks->At(ntrack));
997 if((trackMCAOD->IsPhysicalPrimary())==1)lastprim=ntrack;
998 }
999 for(int ntrack =0;ntrack<ntmc;ntrack++){
1000 trackMCAOD = (AliAODMCParticle*) (mctracks->At(ntrack));
1001 if((trackMCAOD->GetPdgCode()>10)&&((trackMCAOD->GetMother())>1)&&(ntrack>lastprim)&&(trackMCAOD->Charge())){// for Decay particles
1002 MCTrack_pt [ntrack] = trackMCAOD->Pt();
1003 MCTrack_phi [ntrack] = trackMCAOD->Phi();
1004 MCTrack_eta [ntrack] = trackMCAOD->Eta();
1005 double DelPhi = DeltaPhi(Leading_phi,MCTrack_phi[ntrack]);
1006 if(TMath::Abs(MCTrack_eta[ntrack])<0.9){
1007 for(int teb=0;teb<5;teb++){
1008 if(teb==0){if(!( MCTrack_pt[ntrack]>0.15))continue;}
1009 if(teb==1){if(!((MCTrack_pt[ntrack]<1.5)&&(MCTrack_pt[ntrack]>0.15)))continue;}
1010 if(teb==2){if(!((MCTrack_pt[ntrack]<3.0)&&(MCTrack_pt[ntrack]>1.5)))continue;}
1011 if(teb==3){if(!((MCTrack_pt[ntrack]<4.5)&&(MCTrack_pt[ntrack]>3.0)))continue;}
1012 if(teb==4){if(!( MCTrack_pt[ntrack]>4.5))continue;}
1013 fH1JetHadronMC_dphi_ediv [eb][teb]->Fill(DelPhi);
1014 fH1JetHadronMC_dphi_tptweight_ediv [eb][teb]->Fill(DelPhi,MCTrack_pt[ntrack]);
1015 fH1JetHadronMC_dphi_tJptweight_ediv [eb][teb]->Fill(DelPhi,MCTrack_pt[ntrack]/Leading_pt);
1016 if(eb==1){
1017 if((0<Mlead)&&Mlead<7) {fH1JetHadronMC_dphi_tptweight2040_Mleaddep[0][teb]->Fill(DelPhi,MCTrack_pt[ntrack]);}
1018 else if((7<=Mlead)&&(Mlead<10)){fH1JetHadronMC_dphi_tptweight2040_Mleaddep[1][teb]->Fill(DelPhi,MCTrack_pt[ntrack]);}
1019 else {fH1JetHadronMC_dphi_tptweight2040_Mleaddep[2][teb]->Fill(DelPhi,MCTrack_pt[ntrack]);}
1020 if((0<Aj)&&(Aj<0.19)) {fH1JetHadronMC_dphi_tptweight2040_Ajdep [0][teb]->Fill(DelPhi,MCTrack_pt[ntrack]);}
1021 else if((0.19<=Aj)&&(Aj<0.38)) {fH1JetHadronMC_dphi_tptweight2040_Ajdep [1][teb]->Fill(DelPhi,MCTrack_pt[ntrack]);}
1022 else {fH1JetHadronMC_dphi_tptweight2040_Ajdep [2][teb]->Fill(DelPhi,MCTrack_pt[ntrack]);}
1023 }
1024 }
1025 }
1026 }
1027 if((trackMCAOD->IsPhysicalPrimary())&&(trackMCAOD->Charge())){// for Physical particles
1028 MCTrack_pt [ntrack] = trackMCAOD->Pt();
1029 MCTrack_phi [ntrack] = trackMCAOD->Phi();
1030 MCTrack_eta [ntrack] = trackMCAOD->Eta();
1031 double DelPhi = DeltaPhi(Leading_phi,MCTrack_phi[ntrack]);
1032 if(TMath::Abs(MCTrack_eta[ntrack])<0.9){
1033 for(int teb=0;teb<5;teb++){
1034 if(teb==0){if(!( MCTrack_pt[ntrack]>0.15))continue;}
1035 if(teb==1){if(!((MCTrack_pt[ntrack]<1.5)&&(MCTrack_pt[ntrack]>0.15)))continue;}
1036 if(teb==2){if(!((MCTrack_pt[ntrack]<3.0)&&(MCTrack_pt[ntrack]>1.5)))continue;}
1037 if(teb==3){if(!((MCTrack_pt[ntrack]<4.5)&&(MCTrack_pt[ntrack]>3.0)))continue;}
1038 if(teb==4){if(!( MCTrack_pt[ntrack]>4.5))continue;}
1039 fH1JetHadronMC_dphi_ediv [eb][teb]->Fill(DelPhi);
1040 fH1JetHadronMC_dphi_tptweight_ediv [eb][teb]->Fill(DelPhi,MCTrack_pt[ntrack]);
1041 fH1JetHadronMC_dphi_tJptweight_ediv [eb][teb]->Fill(DelPhi,MCTrack_pt[ntrack]/Leading_pt);
1042 fH1JetHadronMCPrim_dphi_ediv [eb][teb]->Fill(DelPhi);
1043 fH1JetHadronMCPrim_dphi_tptweight_ediv [eb][teb]->Fill(DelPhi,MCTrack_pt[ntrack]);
1044 fH1JetHadronMCPrim_dphi_tJptweight_ediv [eb][teb]->Fill(DelPhi,MCTrack_pt[ntrack]/Leading_pt);
1045 if(eb==1){
1046 if((0<Mlead)&&Mlead<7) {fH1JetHadronMC_dphi_tptweight2040_Mleaddep[0][teb]->Fill(DelPhi,MCTrack_pt[ntrack]);}
1047 else if((7<=Mlead)&&(Mlead<10)){fH1JetHadronMC_dphi_tptweight2040_Mleaddep[1][teb]->Fill(DelPhi,MCTrack_pt[ntrack]);}
1048 else {fH1JetHadronMC_dphi_tptweight2040_Mleaddep[2][teb]->Fill(DelPhi,MCTrack_pt[ntrack]);}
1049 if((0<Aj)&&(Aj<0.19)) {fH1JetHadronMC_dphi_tptweight2040_Ajdep [0][teb]->Fill(DelPhi,MCTrack_pt[ntrack]);}
1050 else if((0.19<=Aj)&&(Aj<0.38)) {fH1JetHadronMC_dphi_tptweight2040_Ajdep [1][teb]->Fill(DelPhi,MCTrack_pt[ntrack]);}
1051 else {fH1JetHadronMC_dphi_tptweight2040_Ajdep [2][teb]->Fill(DelPhi,MCTrack_pt[ntrack]);}
1052
1053 if((0<Mlead)&&Mlead<7) {fH1JetHadronMCPrim_dphi_tptweight2040_Mleaddep[0][teb]->Fill(DelPhi,MCTrack_pt[ntrack]);}
1054 else if((7<=Mlead)&&(Mlead<10)){fH1JetHadronMCPrim_dphi_tptweight2040_Mleaddep[1][teb]->Fill(DelPhi,MCTrack_pt[ntrack]);}
1055 else {fH1JetHadronMCPrim_dphi_tptweight2040_Mleaddep[2][teb]->Fill(DelPhi,MCTrack_pt[ntrack]);}
1056 if((0<Aj)&&(Aj<0.19)) {fH1JetHadronMCPrim_dphi_tptweight2040_Ajdep [0][teb]->Fill(DelPhi,MCTrack_pt[ntrack]);}
1057 else if((0.19<=Aj)&&(Aj<0.38)) {fH1JetHadronMCPrim_dphi_tptweight2040_Ajdep [1][teb]->Fill(DelPhi,MCTrack_pt[ntrack]);}
1058 else {fH1JetHadronMCPrim_dphi_tptweight2040_Ajdep [2][teb]->Fill(DelPhi,MCTrack_pt[ntrack]);}
1059 }
1060 }
1061 }
1062 }
1063 }
1064 }
1065 }
1066 }// Momentum Loop Jet
1067 fH2Jet_pt_Munder ->Fill(Leading_pt,(double)Munder/(1.8*pi/2.)*Jet_area[0][nLJet]);
1068 fH2Jet_pt_Mlead ->Fill(Leading_pt,Mlead);
1069 }//eta cut
1070 }// Di-Jet
1071 //############################################### Jet-Hadron Correlation
1072 }// algorithm LOOP
1073 if(IsMC){
1074 for(int eb=0;eb<5;eb++){
1075 double DPhi,DEta;
1076 if(TMath::Abs(ptLJetAOD -10.-20.*(eb))<10.){
1077 DPhi = DeltaPhi(phiLJetMC,phiLJetAOD);
1078 DEta = etaLJetMC-etaLJetAOD;
1079 fH1leadJetMC_dphiResolution[eb]->Fill(DPhi);
1080 if(sqrt(pow(DPhi,2)+pow(DEta,2))<0.4)fH1leadJetMC_Efficiency[eb]->Fill(1);
1081 else fH1leadJetMC_Efficiency[eb]->Fill(0);
1082 DPhi = DeltaPhi(phisLJetMC,phisLJetAOD);
1083 DEta = etasLJetMC-etasLJetAOD;
1084 fH1subJetMC_dphiResolution[eb]->Fill(DPhi);
1085 if(sqrt(pow(DPhi,2)+pow(DEta,2))<0.4)fH1subJetMC_Efficiency[eb]->Fill(1);
1086 else fH1subJetMC_Efficiency[eb]->Fill(0);
1087 DPhi = DeltaPhi(phiLJetMC2,phiLJetAOD);
1088 DEta = etaLJetMC2-etaLJetAOD;
1089
1090 if(sqrt(pow(DPhi,2)+pow(DEta,2))<0.4)fH2leadJetMCptResolution->Fill(ptLJetMC2,ptLJetAOD);
1091 if(findDiJetMC)fH2AjCorrelation_MCRec ->Fill(AjMC,Aj);
1092 if(findDiJetMC)fH2MleadCorrelation_MCRec->Fill(MleadMC,Mlead);
1093 }
1094 }
1095 fH2Mult_Mtrack->Fill(Mult,Track_n);
1096 fH2Mult_Mjet ->Fill(Mult,Mjet_tot);
1097 fH2Mult_Njet ->Fill(Mult,Njet_tot);
1098 if(findLJetAOD)fH2Mult_Mlead ->Fill(Mult,Mlead);
1099 }
1100 else{
1101 fH2Mult_Mtrack->Fill(Mult,Track_n);
1102 fH2Mult_Mjet ->Fill(Mult,Mjet_tot);
1103 fH2Mult_Njet ->Fill(Mult,Njet_tot);
1104 if(findLJetAOD)fH2Mult_Mlead ->Fill(Mult,Mlead);
1105 }
1106
1107 PostData(1, fHistList);
1108 return;
1109}
1110
1111//________________________________________________________________________
1112void AliAnalysisTaskJetHadronCorrelation::Terminate(Option_t *){
1113 // Terminate analysis
1114 if (fDebug) printf("AnalysisTaskPt: Terminate() \n");
1115}
1116
1117Double_t AliAnalysisTaskJetHadronCorrelation::DeltaPhi(Double_t phi1,Double_t phi2){
1118 Float_t pi=TMath::Pi();
1119 Double_t dphi = phi1-phi2;
1120 if (dphi<(-1./2*pi))dphi = dphi +2*pi;
1121 else if(dphi>( 3./2*pi))dphi = dphi -2*pi;
1122 return dphi;
1123}