]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGJE/AliAnalysisTaskJetProperties.cxx
Update for January pA run
[u/mrichter/AliRoot.git] / PWGJE / AliAnalysisTaskJetProperties.cxx
CommitLineData
59a5935a 1// **************************************************************************************
2// * *
3// * Task for Jet properties and jet shape analysis in PWG4 Jet Task Force Train for pp *
4// * *
5// **************************************************************************************
6
7
8/**************************************************************************
9 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
10 * *
11 * Author: The ALICE Off-line Project. *
12 * Contributors are mentioned in the code where appropriate. *
13 * *
14 * Permission to use, copy, modify and distribute this software and its *
15 * documentation strictly for non-commercial purposes is hereby granted *
16 * without fee, provided that the above copyright notice appears in all *
17 * copies and that both the copyright notice and this permission notice *
18 * appear in the supporting documentation. The authors make no claims *
19 * about the suitability of this software for any purpose. It is *
20 * provided "as is" without express or implied warranty. *
21 **************************************************************************/
22
23/* $Id: */
24
25#include "TList.h"
26#include "TH1F.h"
27#include "TH2F.h"
28#include "TH3F.h"
29#include "TString.h"
30#include "THnSparse.h"
31#include "TProfile.h"
32#include "TFile.h"
33#include "TKey.h"
34#include "TRandom3.h"
35
36#include "AliAODInputHandler.h"
37#include "AliAODHandler.h"
38#include "AliESDEvent.h"
39#include "AliAODMCParticle.h"
40#include "AliAODJet.h"
41#include "AliGenPythiaEventHeader.h"
42#include "AliGenHijingEventHeader.h"
43#include "AliInputEventHandler.h"
44
45#include "AliAnalysisHelperJetTasks.h"
46#include "AliAnalysisManager.h"
47#include "AliAnalysisTaskSE.h"
48#include "AliVParticle.h"
49#include "AliAODTrack.h"
50#include "AliVEvent.h"
51
d391281b 52
59a5935a 53#include "AliAnalysisTaskJetProperties.h"
54
55ClassImp(AliAnalysisTaskJetProperties)
56
57//_________________________________________________________________________________//
58
fed89a49 59AliAnalysisTaskJetProperties::AliAnalysisTaskJetProperties()
60: AliAnalysisTaskSE()
61 ,fESD(0)
62 ,fAOD(0)
63 ,fAODJets(0)
64 ,fAODExtension(0)
65 ,fNonStdFile("")
66 ,fBranchJets("jets")
67 ,fTrackType(kTrackAOD)
68 ,fJetRejectType(0)
69 ,fUseAODInputJets(kTRUE)
70 ,fFilterMask(0)
71 ,fUsePhysicsSelection(kTRUE)
72 ,fMaxVertexZ(10)
73 ,fNContributors(2)
74 ,fTrackPtCut(0)
75 ,fTrackEtaMin(0)
76 ,fTrackEtaMax(0)
77 ,fJetPtCut(0)
78 ,fJetEtaMin(0)
79 ,fJetEtaMax(0)
80 ,fAvgTrials(0)
81 ,fJetRadius(0.4)
82 ,fJetList(0)
83 ,fTrackList(0)
84 ,fTrackListUE(0)
85 ,fTrackListJet(0)
86 ,fCommonHistList(0)
87 ,fh1EvtSelection(0)
88 ,fh1VertexNContributors(0)
89 ,fh1VertexZ(0)
90 ,fh1Xsec(0)
91 ,fh1Trials(0)
92 ,fh1PtHard(0)
93 ,fh1PtHardTrials(0)
94 ,fh2EtaJet(0)
95 ,fh2PhiJet(0)
96 ,fh2PtJet(0)
97 ,fh1PtJet(0)
98 ,fh2NtracksJet(0)
99 ,fProNtracksJet(0)
100 ,fh2EtaTrack(0)
101 ,fh2PhiTrack(0)
102 ,fh2PtTrack(0)
103 ,fh2FF(0)
104 ,fh2DelEta(0)
105 ,fh2DelPhi(0)
106 ,fh2DelR(0)
107
108 ,fh1PtLeadingJet(0)
109 ,fh2NtracksLeadingJet(0)
110 ,fProNtracksLeadingJet(0)
111 ,fh2DelR80pcNch(0)
112 ,fProDelR80pcNch(0)
113 ,fh2DelR80pcPt(0)
114 ,fProDelR80pcPt(0)
115 ,fh2AreaCh(0)
116 ,fProAreaCh(0)
117 ,fh3PtDelRNchSum(0)
118 ,fh3PtDelRPtSum(0)
119 ,fProDiffJetShape(0)
120 ,fProIntJetShape(0)
121
122 ,fh1PtSumInJetConeUE(0)
123 ,fh2NtracksLeadingJetUE(0)
124 ,fProNtracksLeadingJetUE(0)
125 ,fh2DelR80pcNchUE(0)
126 ,fProDelR80pcNchUE(0)
127 ,fh2DelR80pcPtUE(0)
128 ,fProDelR80pcPtUE(0)
129 ,fh2AreaChUE(0)
130 ,fProAreaChUE(0)
131 ,fh3PtDelRNchSumUE(0)
132 ,fh3PtDelRPtSumUE(0)
133 ,fProDiffJetShapeUE(0)
134 ,fProIntJetShapeUE(0)
59a5935a 135{
fed89a49 136 for(Int_t ii=0; ii<13; ii++){
efabdd32 137 fProDelRNchSum[ii] = NULL;
138 fProDelRPtSum[ii] = NULL;
139 fProDiffJetShapeA[ii] = NULL;
140 fProIntJetShapeA[ii] = NULL;
fed89a49 141
142 fProDelRNchSumUE[ii] = NULL;
143 fProDelRPtSumUE[ii] = NULL;
144 fProDiffJetShapeAUE[ii] = NULL;
145 fProIntJetShapeAUE[ii] = NULL;
59a5935a 146 }//ii loop
147 // default constructor
148}
149//_________________________________________________________________________________//
150
151AliAnalysisTaskJetProperties::AliAnalysisTaskJetProperties(const char *name)
152 : AliAnalysisTaskSE(name)
153 ,fESD(0)
154 ,fAOD(0)
155 ,fAODJets(0)
156 ,fAODExtension(0)
157 ,fNonStdFile("")
158 ,fBranchJets("jets")
159 ,fTrackType(kTrackAOD)
160 ,fJetRejectType(0)
161 ,fUseAODInputJets(kTRUE)
162 ,fFilterMask(0)
163 ,fUsePhysicsSelection(kTRUE)
164 ,fMaxVertexZ(10)
165 ,fNContributors(2)
166 ,fTrackPtCut(0)
167 ,fTrackEtaMin(0)
168 ,fTrackEtaMax(0)
169 ,fJetPtCut(0)
170 ,fJetEtaMin(0)
171 ,fJetEtaMax(0)
172 ,fAvgTrials(0)
fed89a49 173 ,fJetRadius(0.4)
59a5935a 174 ,fJetList(0)
175 ,fTrackList(0)
fed89a49 176 ,fTrackListUE(0)
177 ,fTrackListJet(0)
59a5935a 178 ,fCommonHistList(0)
179 ,fh1EvtSelection(0)
180 ,fh1VertexNContributors(0)
181 ,fh1VertexZ(0)
182 ,fh1Xsec(0)
183 ,fh1Trials(0)
184 ,fh1PtHard(0)
185 ,fh1PtHardTrials(0)
186 ,fh2EtaJet(0)
187 ,fh2PhiJet(0)
188 ,fh2PtJet(0)
189 ,fh1PtJet(0)
190 ,fh2NtracksJet(0)
191 ,fProNtracksJet(0)
192 ,fh2EtaTrack(0)
193 ,fh2PhiTrack(0)
194 ,fh2PtTrack(0)
195 ,fh2FF(0)
196 ,fh2DelEta(0)
197 ,fh2DelPhi(0)
198 ,fh2DelR(0)
fed89a49 199
59a5935a 200 ,fh1PtLeadingJet(0)
201 ,fh2NtracksLeadingJet(0)
202 ,fProNtracksLeadingJet(0)
203 ,fh2DelR80pcNch(0)
204 ,fProDelR80pcNch(0)
205 ,fh2DelR80pcPt(0)
206 ,fProDelR80pcPt(0)
207 ,fh2AreaCh(0)
208 ,fProAreaCh(0)
209 ,fh3PtDelRNchSum(0)
210 ,fh3PtDelRPtSum(0)
efabdd32 211 ,fProDiffJetShape(0)
212 ,fProIntJetShape(0)
fed89a49 213
214 ,fh1PtSumInJetConeUE(0)
215 ,fh2NtracksLeadingJetUE(0)
216 ,fProNtracksLeadingJetUE(0)
217 ,fh2DelR80pcNchUE(0)
218 ,fProDelR80pcNchUE(0)
219 ,fh2DelR80pcPtUE(0)
220 ,fProDelR80pcPtUE(0)
221 ,fh2AreaChUE(0)
222 ,fProAreaChUE(0)
223 ,fh3PtDelRNchSumUE(0)
224 ,fh3PtDelRPtSumUE(0)
225 ,fProDiffJetShapeUE(0)
226 ,fProIntJetShapeUE(0)
59a5935a 227{
fed89a49 228 for(Int_t ii=0; ii<13; ii++){
efabdd32 229 fProDelRNchSum[ii] = NULL;
230 fProDelRPtSum[ii] = NULL;
231 fProDiffJetShapeA[ii] = NULL;
232 fProIntJetShapeA[ii] = NULL;
fed89a49 233
234 fProDelRNchSumUE[ii] = NULL;
235 fProDelRPtSumUE[ii] = NULL;
236 fProDiffJetShapeAUE[ii] = NULL;
237 fProIntJetShapeAUE[ii] = NULL;
59a5935a 238 }//ii loop
239 // constructor
240 DefineOutput(1,TList::Class());
241}
242//_________________________________________________________________________________//
243
244AliAnalysisTaskJetProperties::~AliAnalysisTaskJetProperties()
245{
246 // destructor
247 if(fJetList) delete fJetList;
248 if(fTrackList) delete fTrackList;
fed89a49 249 if(fTrackListUE) delete fTrackListUE;
250 if(fTrackListJet) delete fTrackListJet;
59a5935a 251}
252//_________________________________________________________________________________//
253
254Bool_t AliAnalysisTaskJetProperties::Notify()
255{
256 //
257 // Implemented Notify() to read the cross sections
258 // and number of trials from pyxsec.root
259 // (taken from AliAnalysisTaskJetSpectrum2)
260 //
261 if(fDebug > 1) Printf("AliAnalysisTaskJetProperties::Notify()");
262
263 TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
264 Float_t xsection = 0;
265 Float_t ftrials = 1;
266
267 fAvgTrials = 1;
268 if(tree){
269 TFile *curfile = tree->GetCurrentFile();
270 if (!curfile) {
271 Error("Notify","No current file");
272 return kFALSE;
273 }
274 if(!fh1Xsec||!fh1Trials){
275 Printf("%s%d No Histogram fh1Xsec",(char*)__FILE__,__LINE__);
276 return kFALSE;
277 }
278 AliAnalysisHelperJetTasks::PythiaInfoFromFile(curfile->GetName(),xsection,ftrials);
279 fh1Xsec->Fill("<#sigma>",xsection);
280 // construct a poor man average trials
281 Float_t nEntries = (Float_t)tree->GetTree()->GetEntries();
282 if(ftrials>=nEntries && nEntries>0.)fAvgTrials = ftrials/nEntries;
283 }
284 return kTRUE;
285}
286//_________________________________________________________________________________//
287
288void AliAnalysisTaskJetProperties::UserCreateOutputObjects()
289{
290 //(Here, event selection part is taken from AliAnalysisTaskFramentationFunction)
291 // create output objects
292 if(fDebug > 1) Printf("AliAnalysisTaskJetProperties::UserCreateOutputObjects()");
293
294 // create list of tracks and jets
295 fJetList = new TList();
296 fJetList->SetOwner(kFALSE);
297 fTrackList = new TList();
298 fTrackList->SetOwner(kFALSE);
fed89a49 299 fTrackListUE = new TList();
300 fTrackListUE->SetOwner(kFALSE);
301 fTrackListJet = new TList();
302 fTrackListJet->SetOwner(kFALSE);
59a5935a 303
304 // Create histograms / output container
305 OpenFile(1);
306 fCommonHistList = new TList();
307 fCommonHistList->SetOwner();
308
309 Bool_t oldStatus = TH1::AddDirectoryStatus();
310 TH1::AddDirectory(kFALSE);
311 // Histograms/TProfiles
312 fh1EvtSelection = new TH1F("fh1EvtSelection", "Event Selection", 6, -0.5, 5.5);
313 fh1EvtSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
314 fh1EvtSelection->GetXaxis()->SetBinLabel(2,"event selection: rejected");
315 fh1EvtSelection->GetXaxis()->SetBinLabel(3,"event class: rejected");
316 fh1EvtSelection->GetXaxis()->SetBinLabel(4,"vertex Ncontr: rejected");
317 fh1EvtSelection->GetXaxis()->SetBinLabel(5,"vertex z: rejected");
318 fh1EvtSelection->GetXaxis()->SetBinLabel(6,"vertex type: rejected");
319
320
321 fh1VertexNContributors = new TH1F("fh1VertexNContributors", "Vertex N contributors", 2500,-.5, 2499.5);
322 fh1VertexZ = new TH1F("fh1VertexZ", "Vertex z distribution", 30, -15., 15.);
323 fh1Xsec = new TProfile("fh1Xsec","xsec from pyxsec.root",1,0,1);
324 fh1Xsec->GetXaxis()->SetBinLabel(1,"<#sigma>");
325 fh1Trials = new TH1F("fh1Trials","trials from pyxsec.root",1,0,1);
326 fh1Trials->GetXaxis()->SetBinLabel(1,"#sum{ntrials}");
327 fh1PtHard = new TH1F("fh1PtHard","PYTHIA Pt hard;p_{T,hard}",350,-.5,349.5);
328 fh1PtHardTrials = new TH1F("fh1PtHardTrials","PYTHIA Pt hard weight with trials;p_{T,hard}",350,-.5,349.5);
329
330
331 Int_t kNbinsPtSlice=20; Float_t xMinPtSlice=0.0; Float_t xMaxPtSlice=100.0;
332 Int_t kNbinsEta=40; Float_t xMinEta=-2.0; Float_t xMaxEta=2.0;
333 Int_t kNbinsPhi=90; Float_t xMinPhi=0.0; Float_t xMaxPhi=TMath::TwoPi();
334 Int_t kNbinsPt=125; Float_t xMinPt=0.0; Float_t xMaxPt=250.0;
335 Int_t kNbinsNtracks=50; Float_t xMinNtracks=0.0; Float_t xMaxNtracks=50.0;
336 Int_t kNbinsFF=80; Float_t xMinFF=0.0; Float_t xMaxFF=2.0;
337 Int_t kNbinsDelR1D=50; Float_t xMinDelR1D=0.0; Float_t xMaxDelR1D=1.0;
338
339
340 fh2EtaJet = new TH2F("EtaJet","EtaJet",
341 kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
342 kNbinsEta, xMinEta, xMaxEta);
343 fh2PhiJet = new TH2F("PhiJet","PhiJet",
344 kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
345 kNbinsPhi, xMinPhi, xMaxPhi);
346 fh2PtJet = new TH2F("PtJet","PtJet",
347 kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
348 kNbinsPt, xMinPt, xMaxPt);
349 fh1PtJet = new TH1F("PtJet1D","PtJet1D",
350 kNbinsPt, xMinPt, xMaxPt);
351 fh2NtracksJet = new TH2F("NtracksJet","NtracksJet",
352 kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
353 kNbinsNtracks, xMinNtracks, xMaxNtracks);
354 fProNtracksJet = new TProfile("AvgNoOfTracksJet","AvgNoOfTracksJet",
355 kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
356 xMinNtracks, xMaxNtracks);
357 fh2EtaTrack = new TH2F("EtaTrack","EtaTrack",
358 kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
359 kNbinsEta, xMinEta, xMaxEta);
360 fh2PhiTrack = new TH2F("PhiTrack","PhiTrack",
361 kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
362 kNbinsPhi, xMinPhi, xMaxPhi);
363 fh2PtTrack = new TH2F("PtTrack","PtTrack",
364 kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
365 kNbinsPt, xMinPt, xMaxPt);
366 fh2FF = new TH2F("FF","FF",
367 kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
368 kNbinsFF, xMinFF, xMaxFF);
369 fh2DelEta = new TH2F("DelEta","DelEta",
370 kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
371 kNbinsEta, xMinEta, xMaxEta);
372 fh2DelPhi = new TH2F("DelPhi","DelPhi",
373 kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
374 kNbinsPhi, xMinPhi, xMaxPhi);
375 fh2DelR = new TH2F("DelR","DelR",
376 kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
377 kNbinsDelR1D, xMinDelR1D, xMaxDelR1D);
378
379
380
381 Int_t kNbinsDelR=100; Float_t xMinDelR=0.0; Float_t xMaxDelR=1.0;
382 Int_t kNbinsPtSliceJS=100; Float_t xMinPtSliceJS=0.0; Float_t xMaxPtSliceJS=250.0;
383
384 fh1PtLeadingJet = new TH1F("PtLeadingJet","PtLeadingJet",
385 kNbinsPt, xMinPt, xMaxPt);
386 fh2NtracksLeadingJet = new TH2F("NtracksLeadingJet","NtracksLeadingJet",
387 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
388 kNbinsNtracks, xMinNtracks, xMaxNtracks);
389 fProNtracksLeadingJet = new TProfile("AvgNoOfTracksLeadingJet","AvgNoOfTracksLeadingJet",
390 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
391 xMinNtracks, xMaxNtracks);
392 fh2DelR80pcNch = new TH2F("delR80pcNch","delR80pcNch",
393 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
394 kNbinsDelR, xMinDelR, xMaxDelR);
395 fProDelR80pcNch = new TProfile("AvgdelR80pcNch","AvgdelR80pcNch",
396 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
397 xMinDelR, xMaxDelR);
398 fh2DelR80pcPt = new TH2F("delR80pcPt","delR80pcPt",
399 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
400 kNbinsDelR, xMinDelR, xMaxDelR);
401 fProDelR80pcPt = new TProfile("AvgdelR80pcPt","AvgdelR80pcPt",
402 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
403 xMinDelR, xMaxDelR);
404 fh2AreaCh = new TH2F("Area","Area",
405 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
406 72, 0.0, TMath::Pi());
407 fProAreaCh = new TProfile("AvgArea","AvgArea",
408 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
409 0.0, TMath::Pi());
410 fh3PtDelRNchSum = new TH3F("PtDelRNchSum","PtDelRNchSum",
411 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
412 kNbinsDelR1D, xMinDelR1D, xMaxDelR1D,
413 kNbinsNtracks, xMinNtracks, xMaxNtracks);
414 fh3PtDelRPtSum = new TH3F("PtDelRPtSum","PtDelRPtSum",
415 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
416 kNbinsDelR1D, xMinDelR1D, xMaxDelR1D,
417 kNbinsPt, xMinPt, xMaxPt);
efabdd32 418 fProDiffJetShape = new TProfile("DiffJetShape","DiffJetShape",
419 10,0.0,1.0,0.0,250.0);
420 fProIntJetShape = new TProfile("IntJetShape","IntJetShape",
421 10,0.0,1.0,0.0,250.0);
422
fed89a49 423
424 fh1PtSumInJetConeUE = new TH1F("PtSumInJetConeUE","PtSumInJetConeUE",
425 500, 0., 100.);
426 fh2NtracksLeadingJetUE = new TH2F("NtracksLeadingJetUE","NtracksLeadingJetUE",
427 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
428 kNbinsNtracks, xMinNtracks, xMaxNtracks);
429 fProNtracksLeadingJetUE = new TProfile("AvgNoOfTracksLeadingJetUE","AvgNoOfTracksLeadingJetUE",
430 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
431 xMinNtracks, xMaxNtracks);
432 fh2DelR80pcNchUE = new TH2F("delR80pcNchUE","delR80pcNchUE",
433 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
434 kNbinsDelR, xMinDelR, xMaxDelR);
435 fProDelR80pcNchUE = new TProfile("AvgdelR80pcNchUE","AvgdelR80pcNchUE",
436 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
437 xMinDelR, xMaxDelR);
438 fh2DelR80pcPtUE = new TH2F("delR80pcPtUE","delR80pcPtUE",
439 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
440 kNbinsDelR, xMinDelR, xMaxDelR);
441 fProDelR80pcPtUE = new TProfile("AvgdelR80pcPtUE","AvgdelR80pcPtUE",
442 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
443 xMinDelR, xMaxDelR);
444 fh2AreaChUE = new TH2F("AreaUE","AreaUE",
445 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
446 72, 0.0, TMath::Pi());
447 fProAreaChUE = new TProfile("AvgAreaUE","AvgAreaUE",
448 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
449 0.0, TMath::Pi());
450 fh3PtDelRNchSumUE = new TH3F("PtDelRNchSumUE","PtDelRNchSumUE",
451 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
452 kNbinsDelR1D, xMinDelR1D, xMaxDelR1D,
453 kNbinsNtracks, xMinNtracks, xMaxNtracks);
454 fh3PtDelRPtSumUE = new TH3F("PtDelRPtSumUE","PtDelRPtSumUE",
455 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
456 kNbinsDelR1D, xMinDelR1D, xMaxDelR1D,
457 kNbinsPt, xMinPt, xMaxPt);
458 fProDiffJetShapeUE = new TProfile("DiffJetShapeUE","DiffJetShapeUE",
459 10,0.0,1.0,0.0,250.0);
460 fProIntJetShapeUE = new TProfile("IntJetShapeUE","IntJetShapeUE",
461 10,0.0,1.0,0.0,250.0);
462
59a5935a 463 TString title;
fed89a49 464 for(Int_t ii=0; ii<13; ii++){
465 if(ii==0)title = "_JetPt20to25";
466 if(ii==1)title = "_JetPt25to30";
467 if(ii==2)title = "_JetPt20to30";
468
469 if(ii==3)title = "_JetPt30to35";
470 if(ii==4)title = "_JetPt35to40";
471 if(ii==5)title = "_JetPt30to40";
472
473 if(ii==6)title = "_JetPt40to50";
474 if(ii==7)title = "_JetPt50to60";
475 if(ii==8)title = "_JetPt40to60";
476
477 if(ii==9) title = "_JetPt60to70";
478 if(ii==10)title = "_JetPt70to80";
479 if(ii==11)title = "_JetPt60to80";
480
481 if(ii==12)title = "_JetPt80to100";
59a5935a 482 fProDelRNchSum[ii] = new TProfile(Form("AvgNchSumDelR%s",title.Data()),Form("AvgNchSumDelR%s",title.Data()),
483 kNbinsDelR1D, xMinDelR1D, xMaxDelR1D,
484 xMinNtracks, xMaxNtracks);
485
486 fProDelRPtSum[ii] = new TProfile(Form("AvgPtSumDelR%s",title.Data()),Form("AvgPtSumDelR%s",title.Data()),
487 kNbinsDelR1D, xMinDelR1D, xMaxDelR1D,
488 xMinPt, xMaxPt);
489 fProDelRNchSum[ii] ->GetXaxis()->SetTitle("R");
490 fProDelRNchSum[ii] ->GetYaxis()->SetTitle("<NchSum>");
491 fProDelRPtSum[ii] ->GetXaxis()->SetTitle("R");
492 fProDelRPtSum[ii] ->GetYaxis()->SetTitle("<PtSum>");
efabdd32 493
494 fProDiffJetShapeA[ii] = new TProfile(Form("DiffJetShape%s",title.Data()),Form("DiffJetShape%s",title.Data()),
495 10,0.0,1.0,0.0,250.0);
496 fProIntJetShapeA[ii] = new TProfile(Form("IntJetShape%s",title.Data()),Form("IntJetShape%s",title.Data()),
497 10,0.0,1.0,0.0,250.0);
498
499 fProDiffJetShapeA[ii]->GetXaxis()->SetTitle("R");
500 fProDiffJetShapeA[ii]->GetYaxis()->SetTitle("Diff jet shape");
501 fProIntJetShapeA[ii]->GetXaxis()->SetTitle("R");
502 fProIntJetShapeA[ii]->GetYaxis()->SetTitle("Integrated jet shape");
503
59a5935a 504 fCommonHistList->Add(fProDelRNchSum[ii]);
505 fCommonHistList->Add(fProDelRPtSum[ii]);
efabdd32 506 fCommonHistList->Add(fProDiffJetShapeA[ii]);
507 fCommonHistList->Add(fProIntJetShapeA[ii]);
fed89a49 508
509 fProDelRNchSumUE[ii] = new TProfile(Form("AvgNchSumDelR%sUE",title.Data()),Form("AvgNchSumDelR%sUE",title.Data()),
510 kNbinsDelR1D, xMinDelR1D, xMaxDelR1D,
511 xMinNtracks, xMaxNtracks);
512
513 fProDelRPtSumUE[ii] = new TProfile(Form("AvgPtSumDelR%sUE",title.Data()),Form("AvgPtSumDelR%sUE",title.Data()),
514 kNbinsDelR1D, xMinDelR1D, xMaxDelR1D,
515 xMinPt, xMaxPt);
516 fProDelRNchSumUE[ii] ->GetXaxis()->SetTitle("R");
517 fProDelRNchSumUE[ii] ->GetYaxis()->SetTitle("<N_{ch}^{Sum,UE}>");
518 fProDelRPtSumUE[ii] ->GetXaxis()->SetTitle("R");
519 fProDelRPtSumUE[ii] ->GetYaxis()->SetTitle("<p_{T}^{sum, UE}>");
520
521 fProDiffJetShapeAUE[ii] = new TProfile(Form("DiffJetShape%sUE",title.Data()),Form("DiffJetShape%sUE",title.Data()),
522 10,0.0,1.0,0.0,250.0);
523 fProIntJetShapeAUE[ii] = new TProfile(Form("IntJetShape%sUE",title.Data()),Form("IntJetShape%sUE",title.Data()),
524 10,0.0,1.0,0.0,250.0);
525
526 fProDiffJetShapeAUE[ii]->GetXaxis()->SetTitle("R");
527 fProDiffJetShapeAUE[ii]->GetYaxis()->SetTitle("Diff jet shape UE");
528 fProIntJetShapeAUE[ii]->GetXaxis()->SetTitle("R");
529 fProIntJetShapeAUE[ii]->GetYaxis()->SetTitle("Integrated jet shape UE");
530
531 fCommonHistList->Add(fProDelRNchSumUE[ii]);
532 fCommonHistList->Add(fProDelRPtSumUE[ii]);
533 fCommonHistList->Add(fProDiffJetShapeAUE[ii]);
534 fCommonHistList->Add(fProIntJetShapeAUE[ii]);
535 }//ii loop
59a5935a 536
537 fh2EtaJet ->GetXaxis()->SetTitle("JetPt"); fh2EtaJet ->GetYaxis()->SetTitle("JetEta");
538 fh2PhiJet ->GetXaxis()->SetTitle("JetPt"); fh2PhiJet ->GetYaxis()->SetTitle("JetPhi");
539 fh2PtJet ->GetXaxis()->SetTitle("JetPt"); fh2PtJet ->GetYaxis()->SetTitle("JetPt");
540 fh1PtJet ->GetXaxis()->SetTitle("JetPt"); fh1PtJet ->GetYaxis()->SetTitle("#jets");
541 fh2NtracksJet ->GetXaxis()->SetTitle("JetPt"); fh2NtracksJet ->GetYaxis()->SetTitle("#tracks");
542 fProNtracksJet->GetXaxis()->SetTitle("JetPt"); fProNtracksJet->GetYaxis()->SetTitle("AgvNtracks");
543 fh2EtaTrack ->GetXaxis()->SetTitle("JetPt"); fh2EtaTrack ->GetYaxis()->SetTitle("TrackEta");
544 fh2PhiTrack ->GetXaxis()->SetTitle("JetPt"); fh2PhiTrack ->GetYaxis()->SetTitle("TrackPhi");
545 fh2PtTrack ->GetXaxis()->SetTitle("JetPt"); fh2PtTrack ->GetYaxis()->SetTitle("TrackPt");
546 fh2FF ->GetXaxis()->SetTitle("JetPt"); fh2FF ->GetYaxis()->SetTitle("FF");
547 fh2DelEta ->GetXaxis()->SetTitle("JetPt"); fh2DelEta ->GetYaxis()->SetTitle("DelEta");
548 fh2DelPhi ->GetXaxis()->SetTitle("JetPt"); fh2DelPhi ->GetYaxis()->SetTitle("DelPhi");
549 fh2DelR ->GetXaxis()->SetTitle("JetPt"); fh2DelR ->GetYaxis()->SetTitle("DelR");
550
551 fh1PtLeadingJet ->GetXaxis()->SetTitle("JetPt"); fh1PtLeadingJet ->GetYaxis()->SetTitle("#leading jets");
552 fh2NtracksLeadingJet ->GetXaxis()->SetTitle("JetPt"); fh2NtracksLeadingJet ->GetYaxis()->SetTitle("#tracks leading jet");
553 fProNtracksLeadingJet ->GetXaxis()->SetTitle("JetPt"); fProNtracksLeadingJet ->GetYaxis()->SetTitle("AvgNtracks leading jet");
554 fh2DelR80pcNch ->GetXaxis()->SetTitle("JetPt"); fh2DelR80pcNch ->GetYaxis()->SetTitle("R containing 80% of tracks");
555 fProDelR80pcNch ->GetXaxis()->SetTitle("JetPt"); fProDelR80pcNch ->GetYaxis()->SetTitle("<R> containing 80% of tracks");
556 fh2DelR80pcPt ->GetXaxis()->SetTitle("JetPt"); fh2DelR80pcPt ->GetYaxis()->SetTitle("R containing 80% of pT");
557 fProDelR80pcPt ->GetXaxis()->SetTitle("JetPt"); fProDelR80pcPt ->GetYaxis()->SetTitle("<R> containing 80% of pT");
558 fh2AreaCh ->GetXaxis()->SetTitle("JetPt"); fh2AreaCh ->GetYaxis()->SetTitle("Jet area");
559 fProAreaCh ->GetXaxis()->SetTitle("JetPt"); fProAreaCh ->GetYaxis()->SetTitle("<jet area>");
560 fh3PtDelRNchSum ->GetXaxis()->SetTitle("JetPt"); fh3PtDelRNchSum ->GetYaxis()->SetTitle("R"); fh3PtDelRNchSum->GetZaxis()->SetTitle("NchSum");
561 fh3PtDelRPtSum ->GetXaxis()->SetTitle("JetPt"); fh3PtDelRPtSum ->GetYaxis()->SetTitle("R"); fh3PtDelRPtSum ->GetZaxis()->SetTitle("PtSum");
efabdd32 562 fProDiffJetShape->GetXaxis()->SetTitle("R");
563 fProDiffJetShape->GetYaxis()->SetTitle("Diff jet shape");
564 fProIntJetShape->GetXaxis()->SetTitle("R");
565 fProIntJetShape->GetYaxis()->SetTitle("Integrated jet shape");
fed89a49 566
567 fh1PtSumInJetConeUE ->GetXaxis()->SetTitle("p_{T}^{sum, UE}(in cone R)"); fh1PtSumInJetConeUE ->GetYaxis()->SetTitle("#leading jets");
568 fh2NtracksLeadingJetUE ->GetXaxis()->SetTitle("JetPt"); fh2NtracksLeadingJetUE ->GetYaxis()->SetTitle("#tracks UE");
569 fProNtracksLeadingJetUE ->GetXaxis()->SetTitle("JetPt"); fProNtracksLeadingJetUE ->GetYaxis()->SetTitle("AvgNtracks UE");
570 fh2DelR80pcNchUE ->GetXaxis()->SetTitle("JetPt"); fh2DelR80pcNchUE ->GetYaxis()->SetTitle("R containing 80% of tracks");
571 fProDelR80pcNchUE ->GetXaxis()->SetTitle("JetPt"); fProDelR80pcNchUE ->GetYaxis()->SetTitle("<R> containing 80% of tracks");
572 fh2DelR80pcPtUE ->GetXaxis()->SetTitle("JetPt"); fh2DelR80pcPtUE ->GetYaxis()->SetTitle("R containing 80% of pT");
573 fProDelR80pcPtUE ->GetXaxis()->SetTitle("JetPt"); fProDelR80pcPtUE ->GetYaxis()->SetTitle("<R> containing 80% of pT");
574 fh2AreaChUE ->GetXaxis()->SetTitle("JetPt"); fh2AreaChUE ->GetYaxis()->SetTitle("UE area");
575 fProAreaChUE ->GetXaxis()->SetTitle("JetPt"); fProAreaChUE ->GetYaxis()->SetTitle("<UE area>");
576 fh3PtDelRNchSumUE ->GetXaxis()->SetTitle("JetPt"); fh3PtDelRNchSumUE ->GetYaxis()->SetTitle("R"); fh3PtDelRNchSumUE->GetZaxis()->SetTitle("NchSumUE");
577 fh3PtDelRPtSumUE ->GetXaxis()->SetTitle("JetPt"); fh3PtDelRPtSumUE ->GetYaxis()->SetTitle("R"); fh3PtDelRPtSumUE ->GetZaxis()->SetTitle("PtSumUE");
578 fProDiffJetShapeUE->GetXaxis()->SetTitle("R");
579 fProDiffJetShapeUE->GetYaxis()->SetTitle("Diff jet shape UE");
580 fProIntJetShapeUE->GetXaxis()->SetTitle("R");
581 fProIntJetShapeUE->GetYaxis()->SetTitle("Integrated jet shape UE");
59a5935a 582
583 fCommonHistList->Add(fh1EvtSelection);
584 fCommonHistList->Add(fh1VertexNContributors);
585 fCommonHistList->Add(fh1VertexZ);
586 fCommonHistList->Add(fh1Xsec);
587 fCommonHistList->Add(fh1Trials);
588 fCommonHistList->Add(fh1PtHard);
589 fCommonHistList->Add(fh1PtHardTrials);
590 fCommonHistList->Add(fh2EtaJet);
591 fCommonHistList->Add(fh2PhiJet);
592 fCommonHistList->Add(fh2PtJet);
593 fCommonHistList->Add(fh1PtJet);
594 fCommonHistList->Add(fh2NtracksJet);
595 fCommonHistList->Add(fProNtracksJet);
596 fCommonHistList->Add(fh2EtaTrack);
597 fCommonHistList->Add(fh2PhiTrack);
598 fCommonHistList->Add(fh2PtTrack);
599 fCommonHistList->Add(fh2FF);
600 fCommonHistList->Add(fh2DelEta);
601 fCommonHistList->Add(fh2DelPhi);
602 fCommonHistList->Add(fh2DelR);
603
604 fCommonHistList->Add(fh1PtLeadingJet);
605 fCommonHistList->Add(fh2NtracksLeadingJet);
606 fCommonHistList->Add(fProNtracksLeadingJet);
607 fCommonHistList->Add(fh2DelR80pcNch);
608 fCommonHistList->Add(fProDelR80pcNch);
609 fCommonHistList->Add(fh2DelR80pcPt);
610 fCommonHistList->Add(fProDelR80pcPt);
611 fCommonHistList->Add(fh2AreaCh);
612 fCommonHistList->Add(fProAreaCh);
613 fCommonHistList->Add(fh3PtDelRNchSum);
614 fCommonHistList->Add(fh3PtDelRPtSum);
efabdd32 615 fCommonHistList->Add(fProDiffJetShape);
616 fCommonHistList->Add(fProIntJetShape);
fed89a49 617
618
619 fCommonHistList->Add(fh1PtSumInJetConeUE);
620 fCommonHistList->Add(fh2NtracksLeadingJetUE);
621 fCommonHistList->Add(fProNtracksLeadingJetUE);
622 fCommonHistList->Add(fh2DelR80pcNchUE);
623 fCommonHistList->Add(fProDelR80pcNchUE);
624 fCommonHistList->Add(fh2DelR80pcPtUE);
625 fCommonHistList->Add(fProDelR80pcPtUE);
626 fCommonHistList->Add(fh2AreaChUE);
627 fCommonHistList->Add(fProAreaChUE);
628 fCommonHistList->Add(fh3PtDelRNchSumUE);
629 fCommonHistList->Add(fh3PtDelRPtSumUE);
630 fCommonHistList->Add(fProDiffJetShapeUE);
631 fCommonHistList->Add(fProIntJetShapeUE);
632
59a5935a 633 // =========== Switch on Sumw2 for all histos ===========
634 for (Int_t i=0; i<fCommonHistList->GetEntries(); ++i){
635 TH1 *h1 = dynamic_cast<TH1*>(fCommonHistList->At(i));
636 if (h1) h1->Sumw2();
637 else{
638 TProfile *hPro = dynamic_cast<TProfile*>(fCommonHistList->At(i));
639 if(hPro) hPro->Sumw2();
640 }
641 }
642
643 TH1::AddDirectory(oldStatus);
644 PostData(1, fCommonHistList);
645}
646//_________________________________________________________________________________//
647
648void AliAnalysisTaskJetProperties::Init()
649{
650 // Initialization
651 if(fDebug > 1) Printf("AliAnalysisTaskJetProperties::Init()");
652
653}
654//_________________________________________________________________________________//
655
656void AliAnalysisTaskJetProperties::UserExec(Option_t *)
657{
658 // Main loop
659 // Called for each event
660 if(fDebug > 1) Printf("AliAnalysisTaskJetProperties::UserExec()");
661
662 if(fDebug > 1) Printf("Analysis event #%5d", (Int_t) fEntry);
663 // Trigger selection
664
665 AliInputEventHandler* inputHandler = (AliInputEventHandler*)
666 ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
667 if(!(inputHandler->IsEventSelected() & AliVEvent::kMB)){
668 if(inputHandler->InheritsFrom("AliESDInputHandler") && fUsePhysicsSelection){ // PhysicsSelection only with ESD input
669 fh1EvtSelection->Fill(1.);
670 if (fDebug > 1 ) Printf(" Trigger Selection: event REJECTED ... ");
671 PostData(1, fCommonHistList);
672 return;
673 }
674 }
675
676 fESD = dynamic_cast<AliESDEvent*>(InputEvent());
677 if(!fESD){
678 if(fDebug>3) Printf("%s:%d ESDEvent not found in the input", (char*)__FILE__,__LINE__);
679 }
680
681 fMCEvent = MCEvent();
682 if(!fMCEvent){
683 if(fDebug>3) Printf("%s:%d MCEvent not found in the input", (char*)__FILE__,__LINE__);
684 }
685
686 // get AOD event from input/ouput
687 TObject* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
688 if( handler && handler->InheritsFrom("AliAODInputHandler") ) {
689 fAOD = ((AliAODInputHandler*)handler)->GetEvent();
690 if(fUseAODInputJets) fAODJets = fAOD;
691 if (fDebug > 1) Printf("%s:%d AOD event from input", (char*)__FILE__,__LINE__);
692 }
693 else {
694 handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
695 if( handler && handler->InheritsFrom("AliAODHandler") ) {
696 fAOD = ((AliAODHandler*)handler)->GetAOD();
697 fAODJets = fAOD;
698 if (fDebug > 1) Printf("%s:%d AOD event from output", (char*)__FILE__,__LINE__);
699 }
700 }
701
702 if(!fAODJets && !fUseAODInputJets){ // case we have AOD in input & output and want jets from output
703 TObject* outHandler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
704 if( outHandler && outHandler->InheritsFrom("AliAODHandler") ) {
705 fAODJets = ((AliAODHandler*)outHandler)->GetAOD();
706 if (fDebug > 1) Printf("%s:%d jets from output AOD", (char*)__FILE__,__LINE__);
707 }
708 }
709
710 if(fNonStdFile.Length()!=0){
711 // case we have an AOD extension - fetch the jets from the extended output
712
713 AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
714 fAODExtension = (aodH?aodH->GetExtension(fNonStdFile.Data()):0);
715 if(!fAODExtension){
716 if(fDebug>1)Printf("AODExtension not found for %s",fNonStdFile.Data());
717 }
718 }
719
720 if(!fAOD){
721 Printf("%s:%d AODEvent not found", (char*)__FILE__,__LINE__);
722 return;
723 }
724 if(!fAODJets){
725 Printf("%s:%d AODEvent with jet branch not found", (char*)__FILE__,__LINE__);
726 return;
727 }
728
729
730 // *** vertex cut ***
731 AliAODVertex* primVtx = fAOD->GetPrimaryVertex();
732 Int_t nTracksPrim = primVtx->GetNContributors();
733 fh1VertexNContributors->Fill(nTracksPrim);
734
735 if (fDebug > 1) Printf("%s:%d primary vertex selection: %d", (char*)__FILE__,__LINE__,nTracksPrim);
736 //if(!nTracksPrim){
737 if(nTracksPrim<fNContributors){
738 if (fDebug > 1) Printf("%s:%d primary vertex selection: event REJECTED...",(char*)__FILE__,__LINE__);
739 fh1EvtSelection->Fill(3.);
740 PostData(1, fCommonHistList);
741 return;
742 }
743
744 fh1VertexZ->Fill(primVtx->GetZ());
745
746 if(TMath::Abs(primVtx->GetZ())>fMaxVertexZ){
747 if (fDebug > 1) Printf("%s:%d primary vertex z = %f: event REJECTED...",(char*)__FILE__,__LINE__,primVtx->GetZ());
748 fh1EvtSelection->Fill(4.);
749 PostData(1, fCommonHistList);
750 return;
751 }
752
753 TString primVtxName(primVtx->GetName());
754
755 if(primVtxName.CompareTo("TPCVertex",TString::kIgnoreCase) == 1){
756 if (fDebug > 1) Printf("%s:%d primary vertex selection: TPC vertex, event REJECTED...",(char*)__FILE__,__LINE__);
757 fh1EvtSelection->Fill(5.);
758 PostData(1, fCommonHistList);
759 return;
760 }
761
762 if (fDebug > 1) Printf("%s:%d event ACCEPTED ...",(char*)__FILE__,__LINE__);
763 fh1EvtSelection->Fill(0.);
764 //___ get MC information __________________________________________________________________
765 Double_t ptHard = 0.;
766 Double_t nTrials = 1; // trials for MC trigger weight for real data
767
768 if(fMCEvent){
769 AliGenEventHeader* genHeader = fMCEvent->GenEventHeader();
770 if(genHeader){
771 AliGenPythiaEventHeader* pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
772 AliGenHijingEventHeader* hijingGenHeader = 0x0;
773
774 if(pythiaGenHeader){
775 if(fDebug>3) Printf("%s:%d pythiaGenHeader found", (char*)__FILE__,__LINE__);
776 nTrials = pythiaGenHeader->Trials();
777 ptHard = pythiaGenHeader->GetPtHard();
778
779 fh1PtHard->Fill(ptHard);
780 fh1PtHardTrials->Fill(ptHard,nTrials);
781
782
783 } else { // no pythia, hijing?
784
785 if(fDebug>3) Printf("%s:%d no pythiaGenHeader found", (char*)__FILE__,__LINE__);
786
787 hijingGenHeader = dynamic_cast<AliGenHijingEventHeader*>(genHeader);
788 if(!hijingGenHeader){
789 Printf("%s:%d no pythiaGenHeader or hjingGenHeader found", (char*)__FILE__,__LINE__);
790 } else {
791 if(fDebug>3) Printf("%s:%d hijingGenHeader found", (char*)__FILE__,__LINE__);
792 }
793 }
794
795 fh1Trials->Fill("#sum{ntrials}",fAvgTrials);
796 }
797 }
798
799 //___ fetch jets __________________________________________________________________________
800 Int_t nJ = GetListOfJets(fJetList);
801 Int_t nJets = 0;
802 if(nJ>=0) nJets = fJetList->GetEntries();
803 if(fDebug>2){
804 Printf("%s:%d Selected jets: %d %d",(char*)__FILE__,__LINE__,nJ,nJets);
805 if(nJ != nJets) Printf("%s:%d Mismatch Selected Jets: %d %d",(char*)__FILE__,__LINE__,nJ,nJets);
806 }
807 FillJetProperties(fJetList);
808 FillJetShape(fJetList);
fed89a49 809 FillJetShapeUE(fJetList);
59a5935a 810 fJetList->Clear();
811 //Post output data.
812 PostData(1, fCommonHistList);
813}
814//_________________________________________________________________________________//
815
816void AliAnalysisTaskJetProperties::Terminate(Option_t *)
817{
818 // terminated
819
820 if(fDebug > 1) printf("AliAnalysisTaskJetProperties::Terminate() \n");
821}
822//_________________________________________________________________________________//
823
824Int_t AliAnalysisTaskJetProperties::GetListOfJets(TList *list)
825{
826 //this functionality is motivated by AliAnalysisTaskFragmentationFunction
827 if(fDebug > 1) printf("AliAnalysisTaskJetProperties::GetListOfJets() \n");
828 // fill list of jets selected according to type
829 if(!list){
830 if(fDebug>1) Printf("%s:%d no input list", (char*)__FILE__,__LINE__);
831 return -1;
832 }
833
834 if(fBranchJets.Length()==0){
835 Printf("%s:%d no jet branch specified", (char*)__FILE__,__LINE__);
836 if(fDebug>1)fAOD->Print();
837 return 0;
838 }
839 TClonesArray *aodJets = 0;
840 if(fBranchJets.Length()) aodJets = dynamic_cast<TClonesArray*>(fAODJets->FindListObject(fBranchJets.Data()));
841 if(!aodJets) aodJets = dynamic_cast<TClonesArray*>(fAODJets->GetList()->FindObject(fBranchJets.Data()));
842 if(fAODExtension&&!aodJets) aodJets = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fBranchJets.Data()));
843
844 if(!aodJets){
845 if(fBranchJets.Length())Printf("%s:%d no jet array with name %s in AOD", (char*)__FILE__,__LINE__,fBranchJets.Data());
846 if(fDebug>1)fAOD->Print();
847 return 0;
848 }
849
850 Int_t nJets = 0;
851 for(Int_t ijet=0; ijet<aodJets->GetEntries(); ijet++){
852 AliAODJet *tmp = dynamic_cast<AliAODJet*>(aodJets->At(ijet));
853 if(!tmp) continue;
854 if( tmp->Pt() < fJetPtCut ) continue;
855 if( tmp->Eta() < fJetEtaMin || tmp->Eta() > fJetEtaMax)continue;
856 if(fJetRejectType==kReject1Track && tmp->GetRefTracks()->GetEntriesFast()==1)continue;//rejecting 1track jet if...
857 list->Add(tmp);
858 nJets++;
859 }//ij loop
860 list->Sort();
861 return nJets;
862}
863//_________________________________________________________________________________//
864
865Int_t AliAnalysisTaskJetProperties::GetListOfJetTracks(TList* list, const AliAODJet* jet)
866{
867 //this functionality is motivated by AliAnalysisTaskFragmentationFunction
868 if(fDebug > 1) printf("AliAnalysisTaskJetProperties::GetListOfJetTracks() \n");
869
870 // list of jet tracks from trackrefs
871 Int_t nTracks = jet->GetRefTracks()->GetEntriesFast();
872 Int_t NoOfTracks=0;
873 for (Int_t itrack=0; itrack<nTracks; itrack++) {
874 if(fTrackType==kTrackUndef){
875 if(fDebug>2)Printf("%s:%d unknown track type %d in AOD", (char*)__FILE__,__LINE__,kTrackUndef);
876 return 0;
877 }//if
878 else if(fTrackType == kTrackAOD){
879 AliAODTrack* trackAod = dynamic_cast<AliAODTrack*>(jet->GetRefTracks()->At(itrack));
880 if(!trackAod){
881 AliError("expected ref track not found ");
882 continue;
883 }
884 list->Add(trackAod);
885 NoOfTracks++;
886 }//if
887
888 else if(fTrackType==kTrackAODMC){
889 AliAODMCParticle* trackmc = dynamic_cast<AliAODMCParticle*>(jet->GetRefTracks()->At(itrack));
890 if(!trackmc){
891 AliError("expected ref trackmc not found ");
892 continue;
893 }
894 list->Add(trackmc);
895 NoOfTracks++;
896 }//if
897
898 else if (fTrackType==kTrackKine){
899 AliVParticle* trackkine = dynamic_cast<AliVParticle*>(jet->GetRefTracks()->At(itrack));
900 if(!trackkine){
901 AliError("expected ref trackkine not found ");
902 continue;
903 }
904 list->Add(trackkine);
905 NoOfTracks++;
906 }//if
907 }//itrack loop
908 list->Sort();
909 return NoOfTracks;
910}//initial loop
911//_________________________________________________________________________________//
912
913void AliAnalysisTaskJetProperties::FillJetProperties(TList *jetlist){
914 //filling up the histograms jets and tracks inside jet
915 if(fDebug > 1) printf("AliAnalysisTaskJetProperties::FillJetProperties() \n");
916
917 for(Int_t iJet=0; iJet < jetlist->GetEntries(); iJet++){
918 Float_t JetPt;Float_t JetPhi;Float_t JetEta;Float_t JetE;
919 AliAODJet *jet = dynamic_cast<AliAODJet*>(jetlist->At(iJet));
920 if(!jet)continue;
921 JetEta = jet->Eta();
922 JetPhi = jet->Phi();
923 JetPt = jet->Pt();
924 JetE = jet->E();
925 fh2EtaJet ->Fill(JetPt,JetEta);
926 fh2PhiJet ->Fill(JetPt,JetPhi);
927 fh2PtJet ->Fill(JetPt,JetPt);
928 fh1PtJet ->Fill(JetPt);
929
fed89a49 930 Int_t nJT = GetListOfJetTracks(fTrackListJet,jet);
59a5935a 931 Int_t nJetTracks = 0;
fed89a49 932 if(nJT>=0) nJetTracks = fTrackListJet->GetEntries();
59a5935a 933 if(fDebug>2){
934 Printf("%s:%d Jet tracks: %d %d",(char*)__FILE__,__LINE__,nJT,nJetTracks);
935 if(nJT != nJetTracks) Printf("%s:%d Mismatch Jet Tracks: %d %d",(char*)__FILE__,__LINE__,nJT,nJetTracks);
936 }
937
fed89a49 938 fh2NtracksJet ->Fill(JetPt,fTrackListJet->GetEntries());
939 fProNtracksJet ->Fill(JetPt,fTrackListJet->GetEntries());
59a5935a 940
fed89a49 941 for (Int_t j =0; j< fTrackListJet->GetEntries(); j++){
59a5935a 942 if(fTrackType==kTrackUndef)continue;
943 Float_t TrackEta=-99.0; Float_t TrackPt=-99.0; Float_t TrackPhi=-99.0;
944 Float_t FF=-99.0; Float_t DelEta=-99.0; Float_t DelPhi=-99.0; Float_t DelR=-99.0; Float_t AreaJ=-99.0;
945 if(fTrackType==kTrackAOD){
fed89a49 946 AliAODTrack *trackaod = dynamic_cast<AliAODTrack*>(fTrackListJet->At(j));
59a5935a 947 if(!trackaod)continue;
948 TrackEta = trackaod->Eta();
949 TrackPhi = trackaod->Phi();
950 TrackPt = trackaod->Pt();
951 }//if kTrackAOD
952 else if(fTrackType==kTrackAODMC){
fed89a49 953 AliAODMCParticle* trackmc = dynamic_cast<AliAODMCParticle*>(fTrackListJet->At(j));
59a5935a 954 if(!trackmc)continue;
955 TrackEta = trackmc->Eta();
956 TrackPhi = trackmc->Phi();
957 TrackPt = trackmc->Pt();
958 }//if kTrackAODMC
959 else if(fTrackType==kTrackKine){
fed89a49 960 AliVParticle* trackkine = dynamic_cast<AliVParticle*>(fTrackListJet->At(j));
59a5935a 961 if(!trackkine)continue;
962 TrackEta = trackkine->Eta();
963 TrackPhi = trackkine->Phi();
964 TrackPt = trackkine->Pt();
965 }//kTrackKine
966 if(JetPt)FF = TrackPt/JetPt;
967 DelEta = TMath::Abs(JetEta - TrackEta);
968 DelPhi = TMath::Abs(JetPhi - TrackPhi);
448afeb9 969 if(DelPhi>TMath::Pi())DelPhi = TMath::Abs(DelPhi-TMath::TwoPi());
59a5935a 970 DelR = TMath::Sqrt(DelEta*DelEta + DelPhi*DelPhi);
971 AreaJ = TMath::Pi()*DelR*DelR;
972 fh2EtaTrack ->Fill(JetPt,TrackEta);
973 fh2PhiTrack ->Fill(JetPt,TrackPhi);
974 fh2PtTrack ->Fill(JetPt,TrackPt);
975 fh2FF ->Fill(JetPt,FF);
976 fh2DelEta ->Fill(JetPt,DelEta);
977 fh2DelPhi ->Fill(JetPt,DelPhi);
978 fh2DelR ->Fill(JetPt,DelR);
979 }//track loop
fed89a49 980 fTrackListJet->Clear();
59a5935a 981 }//iJet loop
982}//FillJetProperties
983//_________________________________________________________________________________//
984
985void AliAnalysisTaskJetProperties::FillJetShape(TList *jetlist){
986 //filling up the histograms
987 if(fDebug > 1) printf("AliAnalysisTaskJetProperties::FillJetShape() \n");
988 Float_t JetEta; Float_t JetPhi; Float_t JetPt;
fed89a49 989 AliAODJet *jet = dynamic_cast<AliAODJet*>(jetlist->At(0));//Leading jet only
59a5935a 990 if(jet){
991 JetEta = jet->Eta();
992 JetPhi = jet->Phi();
993 JetPt = jet->Pt();
994 fh1PtLeadingJet->Fill(JetPt);
efabdd32 995 Float_t NchSumA[50] = {0.};
996 Float_t PtSumA[50] = {0.};
997 Float_t delRPtSum80pc = 0;
998 Float_t delRNtrkSum80pc = 0;
999 Float_t PtSumDiffShape[10] = {0.0};
1000 Float_t PtSumIntShape[10] = {0.0};
1001 Int_t kNbinsR = 10;
1002
fed89a49 1003 Int_t nJT = GetListOfJetTracks(fTrackListJet,jet);
59a5935a 1004 Int_t nJetTracks = 0;
fed89a49 1005 if(nJT>=0) nJetTracks = fTrackListJet->GetEntries();
59a5935a 1006 if(fDebug>2){
1007 Printf("%s:%d Jet tracks: %d %d",(char*)__FILE__,__LINE__,nJT,nJetTracks);
1008 if(nJT != nJetTracks) Printf("%s:%d Mismatch Jet Tracks: %d %d",(char*)__FILE__,__LINE__,nJT,nJetTracks);
1009 }
1010 fh2NtracksLeadingJet->Fill(JetPt,nJetTracks);
1011 fProNtracksLeadingJet->Fill(JetPt,nJetTracks);
1012 Int_t *index = new Int_t [nJetTracks];//dynamic array containing index
1013 Float_t *delRA = new Float_t [nJetTracks];//dynamic array containing delR
1014 Float_t *delEtaA = new Float_t [nJetTracks];//dynamic array containing delEta
1015 Float_t *delPhiA = new Float_t [nJetTracks];//dynamic array containing delPhi
1016 Float_t *trackPtA = new Float_t [nJetTracks];//dynamic array containing pt-track
1017 Float_t *trackEtaA = new Float_t [nJetTracks];//dynamic array containing eta-track
1018 Float_t *trackPhiA = new Float_t [nJetTracks];//dynamic array containing phi-track
1019 for(Int_t ii=0; ii<nJetTracks; ii++){
1020 index[ii] = 0;
1021 delRA[ii] = 0.;
1022 delEtaA[ii] = 0.;
1023 delPhiA[ii] = 0.;
1024 trackPtA[ii] = 0.;
1025 trackEtaA[ii] = 0.;
1026 trackPhiA[ii] = 0.;
1027 }//ii loop
1028
1029 for (Int_t j =0; j< nJetTracks; j++){
1030 if(fTrackType==kTrackUndef)continue;
1031
1032 Float_t TrackEta=-99.0; Float_t TrackPt=-99.0; Float_t TrackPhi=-99.0;
1033 Float_t DelEta=-99.0; Float_t DelPhi=-99.0; Float_t DelR=-99.0; Float_t AreaJ=-99.0;
1034
1035 if(fTrackType==kTrackAOD){
fed89a49 1036 AliAODTrack *trackaod = dynamic_cast<AliAODTrack*>(fTrackListJet->At(j));
59a5935a 1037 if(!trackaod)continue;
1038 TrackEta = trackaod->Eta();
1039 TrackPhi = trackaod->Phi();
1040 TrackPt = trackaod->Pt();
1041 }//if kTrackAOD
1042 else if(fTrackType==kTrackAODMC){
fed89a49 1043 AliAODMCParticle* trackmc = dynamic_cast<AliAODMCParticle*>(fTrackListJet->At(j));
59a5935a 1044 if(!trackmc)continue;
1045 TrackEta = trackmc->Eta();
1046 TrackPhi = trackmc->Phi();
1047 TrackPt = trackmc->Pt();
1048 }//if kTrackAODMC
1049 else if(fTrackType==kTrackKine){
fed89a49 1050 AliVParticle* trackkine = dynamic_cast<AliVParticle*>(fTrackListJet->At(j));
59a5935a 1051 if(!trackkine)continue;
1052 TrackEta = trackkine->Eta();
1053 TrackPhi = trackkine->Phi();
1054 TrackPt = trackkine->Pt();
1055 }//if kTrackKine
1056
1057 DelEta = TMath::Abs(JetEta - TrackEta);
1058 DelPhi = TMath::Abs(JetPhi - TrackPhi);
efabdd32 1059 if(DelPhi>TMath::Pi())DelPhi = TMath::Abs(DelPhi-TMath::TwoPi());
59a5935a 1060 DelR = TMath::Sqrt(DelEta*DelEta + DelPhi*DelPhi);
1061 AreaJ = TMath::Pi()*DelR*DelR;
1062
1063 fh2AreaCh ->Fill(JetPt,AreaJ);
1064 fProAreaCh->Fill(JetPt,AreaJ);
1065 delRA[j] = DelR;
1066 delEtaA[j] = DelEta;
1067 delPhiA[j] = DelPhi;
1068 trackPtA[j] = TrackPt;
1069 trackEtaA[j] = TrackEta;
1070 trackPhiA[j] = TrackPhi;
1071
efabdd32 1072 //calculating diff and integrated jet shapes
1073 Float_t kDeltaR = 0.1;
1074 Float_t RMin = kDeltaR/2.0;
1075 Float_t RMax = kDeltaR/2.0;
1076 Float_t tmpR = 0.05;
1077 for(Int_t ii1=0; ii1<kNbinsR;ii1++){
1078 if((DelR > (tmpR-RMin)) && (DelR <=(tmpR+RMax)))PtSumDiffShape[ii1]+= TrackPt;
1079 if(DelR>0.0 && DelR <=(tmpR+RMax))PtSumIntShape[ii1]+= TrackPt;
1080 tmpR += 0.1;
1081 }//ii1 loop
1082
59a5935a 1083 for(Int_t ibin=1; ibin<=50; ibin++){
1084 Float_t xlow = 0.02*(ibin-1);
1085 Float_t xup = 0.02*ibin;
1086 if( xlow <= DelR && DelR < xup){
1087 NchSumA[ibin-1]++;
1088 PtSumA[ibin-1]+= TrackPt;
1089 }//if loop
1090 }//for ibin loop
1091 }//track loop
fed89a49 1092 fTrackListJet->Clear();
efabdd32 1093
1094 //---------------------
1095 Float_t tmp1R = 0.05;
1096 for(Int_t jj1=0; jj1<kNbinsR;jj1++){
1097 if(JetPt>20 && JetPt<=100){
1098 fProDiffJetShape->Fill(tmp1R,PtSumDiffShape[jj1]/JetPt);
1099 fProIntJetShape ->Fill(tmp1R,PtSumIntShape[jj1]/JetPt);
1100 }
1101 Float_t jetPtMin0=20.0; Float_t jetPtMax0=30.0;
fed89a49 1102 for(Int_t k=0; k<13; k++){
1103 if(k==0) {jetPtMin0=20.0;jetPtMax0=25.0;}
1104 if(k==1) {jetPtMin0=25.0;jetPtMax0=30.0;}
1105 if(k==2) {jetPtMin0=20.0;jetPtMax0=30.0;}
1106 if(k==3) {jetPtMin0=30.0;jetPtMax0=35.0;}
1107 if(k==4) {jetPtMin0=35.0;jetPtMax0=40.0;}
1108 if(k==5) {jetPtMin0=30.0;jetPtMax0=40.0;}
1109 if(k==6) {jetPtMin0=40.0;jetPtMax0=50.0;}
1110 if(k==7) {jetPtMin0=50.0;jetPtMax0=60.0;}
1111 if(k==8) {jetPtMin0=40.0;jetPtMax0=60.0;}
1112 if(k==9) {jetPtMin0=60.0;jetPtMax0=70.0;}
1113 if(k==10){jetPtMin0=70.0;jetPtMax0=80.0;}
1114 if(k==11){jetPtMin0=60.0;jetPtMax0=80.0;}
1115 if(k==12){jetPtMin0=80.0;jetPtMax0=100.0;}
efabdd32 1116 if(JetPt>jetPtMin0 && JetPt<=jetPtMax0){
1117 fProDiffJetShapeA[k]->Fill(tmp1R,PtSumDiffShape[jj1]/JetPt);
1118 fProIntJetShapeA[k] ->Fill(tmp1R,PtSumIntShape[jj1]/JetPt);
1119 }//if
1120 }//k loop
1121 tmp1R +=0.1;
1122 }//jj1 loop
1123 //----------------------//
59a5935a 1124 Float_t PtSum = 0;
1125 Int_t NtrkSum = 0;
1126 Bool_t iflagPtSum = kFALSE;
1127 Bool_t iflagNtrkSum = kFALSE;
1128 TMath::Sort(nJetTracks,delRA,index,0);
1129 for(Int_t ii=0; ii<nJetTracks; ii++){
1130 NtrkSum ++;
1131 PtSum += trackPtA[index[ii]];
1132 /*
1133 cout << index[ii] << "\t" <<
1134 delR[ii] << "\t" <<
1135 delEta[ii]<< "\t" <<
1136 delPhi[ii]<< "\t" <<
1137 trackPt[ii]<< "\t" <<
1138 trackEta[ii]<< "\t" <<
1139 trackPhi[ii]<< "\t DelR " <<
1140 delR[index[ii]] << endl;
1141 */
1142 if(!iflagNtrkSum){
1143 if((Float_t)NtrkSum/(Float_t)nJetTracks > 0.79){
1144 delRNtrkSum80pc = delRA[index[ii]];
1145 iflagNtrkSum = kTRUE;
1146 }//if NtrkSum
1147 }//if iflag
1148 if(!iflagPtSum){
1149 if(PtSum/JetPt >= 0.8000){
1150 delRPtSum80pc = delRA[index[ii]];
1151 iflagPtSum = kTRUE;
1152 }//if PtSum
1153 }//if iflag
1154 }//track loop 2nd
1155 delete [] index;
1156 delete [] delRA;
1157 delete [] delEtaA;
1158 delete [] delPhiA;
1159 delete [] trackPtA;
1160 delete [] trackEtaA;
1161 delete [] trackPhiA;
1162 fh2DelR80pcNch ->Fill(JetPt,delRNtrkSum80pc);
1163 fProDelR80pcNch->Fill(JetPt,delRNtrkSum80pc);
1164 fh2DelR80pcPt ->Fill(JetPt,delRPtSum80pc);
1165 fProDelR80pcPt ->Fill(JetPt,delRPtSum80pc);
1166 for(Int_t ibin=0; ibin<50; ibin++){
1167 Float_t iR = 0.02*ibin + 0.01;
1168 fh3PtDelRNchSum->Fill(JetPt,iR,NchSumA[ibin]);
1169 fh3PtDelRPtSum ->Fill(JetPt,iR,PtSumA[ibin]);
1170 Float_t jetPtMin=20.0; Float_t jetPtMax=30.0;
fed89a49 1171 for(Int_t k=0; k<13; k++){
1172 if(k==0) {jetPtMin=20.0;jetPtMax=25.0;}
1173 if(k==1) {jetPtMin=25.0;jetPtMax=30.0;}
1174 if(k==2) {jetPtMin=20.0;jetPtMax=30.0;}
1175 if(k==3) {jetPtMin=30.0;jetPtMax=35.0;}
1176 if(k==4) {jetPtMin=35.0;jetPtMax=40.0;}
1177 if(k==5) {jetPtMin=30.0;jetPtMax=40.0;}
1178 if(k==6) {jetPtMin=40.0;jetPtMax=50.0;}
1179 if(k==7) {jetPtMin=50.0;jetPtMax=60.0;}
1180 if(k==8) {jetPtMin=40.0;jetPtMax=60.0;}
1181 if(k==9) {jetPtMin=60.0;jetPtMax=70.0;}
1182 if(k==10){jetPtMin=70.0;jetPtMax=80.0;}
1183 if(k==11){jetPtMin=60.0;jetPtMax=80.0;}
1184 if(k==12){jetPtMin=80.0;jetPtMax=100.0;}
59a5935a 1185 if(JetPt>jetPtMin && JetPt<jetPtMax){
1186 fProDelRNchSum[k]->Fill(iR,NchSumA[ibin]);
1187 fProDelRPtSum[k]->Fill(iR,PtSumA[ibin]);
1188 }//if
1189 }//k loop
1190 }//ibin loop
1191 }//if(jet)
1192}//FillJetShape()
1193//_________________________________________________________________________________//
1194
fed89a49 1195void AliAnalysisTaskJetProperties::FillJetShapeUE(TList *jetlist){
1196 //filling up the histograms
1197 if(fDebug > 1) printf("AliAnalysisTaskJetProperties::FillJetShape() \n");
1198 AliAODJet *jet = dynamic_cast<AliAODJet*>(jetlist->At(0));//Leading jet only
1199 if(jet){
1200 Double_t jetMom[3];
1201 jet->PxPyPz(jetMom);
1202
1203 TVector3 jet3mom(jetMom);
1204 // Rotate phi and keep eta unchanged
1205
1206 const Double_t alpha = TMath::Pi()/2.;
1207 Double_t etaTilted = jet3mom.Eta();
1208 Double_t phiTilted = TVector2::Phi_0_2pi(jet3mom.Phi()) + alpha;
1209 if(phiTilted > 2*TMath::Pi()) phiTilted = phiTilted - 2*TMath::Pi();
1210
1211 Double_t JetPt = jet->Pt();
1212
1213 Float_t NchSumA[50] = {0.};
1214 Float_t PtSumA[50] = {0.};
1215 Float_t delRPtSum80pc = 0;
1216 Float_t delRNtrkSum80pc = 0;
1217 Float_t PtSumDiffShape[10] = {0.0};
1218 Float_t PtSumIntShape[10] = {0.0};
1219 Int_t kNbinsR = 10;
1220
1221 //Int_t nTracks = GetListOfTracks(fTrackList,fTrackType);
1222 GetListOfTracks(fTrackList,fTrackType);
1223 Double_t sumPtPerp = 0;
1224 if(fDebug > 100) printf("Cone radius for bckg. = %f Track type = %d\n",fJetRadius, fTrackType);
1225 GetTracksTiltedwrpJetAxis(alpha, fTrackList, fTrackListUE,jet,fJetRadius,sumPtPerp);
1226 fTrackList->Clear();
1227 fh1PtSumInJetConeUE->Fill(sumPtPerp);
1228
1229 Int_t nTracksUE = fTrackListUE->GetEntries();
1230 fh2NtracksLeadingJetUE->Fill(JetPt,nTracksUE);
1231 fProNtracksLeadingJetUE->Fill(JetPt,nTracksUE);
1232
1233 Int_t *index = new Int_t [nTracksUE];//dynamic array containing index
1234 Float_t *delRA = new Float_t [nTracksUE];//dynamic array containing delR
1235 Float_t *delEtaA = new Float_t [nTracksUE];//dynamic array containing delEta
1236 Float_t *delPhiA = new Float_t [nTracksUE];//dynamic array containing delPhi
1237 Float_t *trackPtA = new Float_t [nTracksUE];//dynamic array containing pt-track
1238 Float_t *trackEtaA = new Float_t [nTracksUE];//dynamic array containing eta-track
1239 Float_t *trackPhiA = new Float_t [nTracksUE];//dynamic array containing phi-track
1240 for(Int_t ii=0; ii<nTracksUE; ii++){
1241 index[ii] = 0;
1242 delRA[ii] = 0.;
1243 delEtaA[ii] = 0.;
1244 delPhiA[ii] = 0.;
1245 trackPtA[ii] = 0.;
1246 trackEtaA[ii] = 0.;
1247 trackPhiA[ii] = 0.;
1248 }//ii loop
1249
1250 for (Int_t j =0; j< nTracksUE; j++){
1251 if(fTrackType==kTrackUndef)continue;
1252
1253 Float_t TrackEta=-99.0; Float_t TrackPt=-99.0; Float_t TrackPhi=-99.0;
1254 Float_t DelEta=-99.0; Float_t DelPhi=-99.0; Float_t DelR=-99.0; Float_t AreaJ=-99.0;
1255
1256 if(fTrackType==kTrackAOD){
1257 AliAODTrack *trackaod = dynamic_cast<AliAODTrack*>(fTrackListUE->At(j));
1258 if(!trackaod)continue;
1259 TrackEta = trackaod->Eta();
1260 TrackPhi = trackaod->Phi();
1261 TrackPt = trackaod->Pt();
1262 if(fDebug > 100) printf("FillJetShapeUE itrack, trackPt %d, %f\n",j,TrackPt);
1263 }//if kTrackAOD
1264 else if(fTrackType==kTrackAODMC){
1265 AliAODMCParticle* trackmc = dynamic_cast<AliAODMCParticle*>(fTrackListUE->At(j));
1266 if(!trackmc)continue;
1267 TrackEta = trackmc->Eta();
1268 TrackPhi = trackmc->Phi();
1269 TrackPt = trackmc->Pt();
1270 if(fDebug > 100) printf("FillJetShapeUE itrack, trackPt %d, %f\n",j,TrackPt);
1271 }//if kTrackAODMC
1272 else if(fTrackType==kTrackKine){
1273 AliVParticle* trackkine = dynamic_cast<AliVParticle*>(fTrackListUE->At(j));
1274 if(!trackkine)continue;
1275 TrackEta = trackkine->Eta();
1276 TrackPhi = trackkine->Phi();
1277 TrackPt = trackkine->Pt();
1278 if(fDebug > 100) printf("FillJetShapeUE itrack, trackPt %d, %f\n",j,TrackPt);
1279 }//if kTrackKine
1280
1281 DelEta = TMath::Abs(etaTilted - TrackEta);
1282 DelPhi = TMath::Abs(phiTilted - TrackPhi);
1283 if(DelPhi>TMath::Pi())DelPhi = TMath::Abs(DelPhi-TMath::TwoPi());
1284 DelR = TMath::Sqrt(DelEta*DelEta + DelPhi*DelPhi);
1285 AreaJ = TMath::Pi()*DelR*DelR;
1286
1287 fh2AreaChUE ->Fill(JetPt,AreaJ);
1288 fProAreaChUE->Fill(JetPt,AreaJ);
1289
1290 delRA[j] = DelR;
1291 delEtaA[j] = DelEta;
1292 delPhiA[j] = DelPhi;
1293 trackPtA[j] = TrackPt;
1294 trackEtaA[j] = TrackEta;
1295 trackPhiA[j] = TrackPhi;
1296
1297 //calculating diff and integrated jet shapes
1298 Float_t kDeltaR = 0.1;
1299 Float_t RMin = kDeltaR/2.0;
1300 Float_t RMax = kDeltaR/2.0;
1301 Float_t tmpR = 0.05;
1302 for(Int_t ii1=0; ii1<kNbinsR;ii1++){
1303 if((DelR > (tmpR-RMin)) && (DelR <=(tmpR+RMax)))PtSumDiffShape[ii1]+= TrackPt;
1304 if(DelR>0.0 && DelR <=(tmpR+RMax))PtSumIntShape[ii1]+= TrackPt;
1305 tmpR += 0.1;
1306 }//ii1 loop
1307
1308 for(Int_t ibin=1; ibin<=50; ibin++){
1309 Float_t xlow = 0.02*(ibin-1);
1310 Float_t xup = 0.02*ibin;
1311 if( xlow <= DelR && DelR < xup){
1312 NchSumA[ibin-1]++;
1313 PtSumA[ibin-1]+= TrackPt;
1314 }//if loop
1315 }//for ibin loop
1316 }//track loop
1317 fTrackListUE->Clear();
1318
1319 //---------------------
1320 Float_t tmp1R = 0.05;
1321 for(Int_t jj1=0; jj1<kNbinsR;jj1++){
1322 if(JetPt>20 && JetPt<=100){
1323 fProDiffJetShapeUE->Fill(tmp1R,PtSumDiffShape[jj1]/JetPt);
1324 fProIntJetShapeUE ->Fill(tmp1R,PtSumIntShape[jj1]/JetPt);
1325 }
1326 Float_t jetPtMin0=20.0; Float_t jetPtMax0=30.0;
1327 for(Int_t k=0; k<13; k++){
1328 if(k==0) {jetPtMin0=20.0;jetPtMax0=25.0;}
1329 if(k==1) {jetPtMin0=25.0;jetPtMax0=30.0;}
1330 if(k==2) {jetPtMin0=20.0;jetPtMax0=30.0;}
1331 if(k==3) {jetPtMin0=30.0;jetPtMax0=35.0;}
1332 if(k==4) {jetPtMin0=35.0;jetPtMax0=40.0;}
1333 if(k==5) {jetPtMin0=30.0;jetPtMax0=40.0;}
1334 if(k==6) {jetPtMin0=40.0;jetPtMax0=50.0;}
1335 if(k==7) {jetPtMin0=50.0;jetPtMax0=60.0;}
1336 if(k==8) {jetPtMin0=40.0;jetPtMax0=60.0;}
1337 if(k==9) {jetPtMin0=60.0;jetPtMax0=70.0;}
1338 if(k==10){jetPtMin0=70.0;jetPtMax0=80.0;}
1339 if(k==11){jetPtMin0=60.0;jetPtMax0=80.0;}
1340 if(k==12){jetPtMin0=80.0;jetPtMax0=100.0;}
1341 if(JetPt>jetPtMin0 && JetPt<=jetPtMax0){
1342 fProDiffJetShapeAUE[k]->Fill(tmp1R,PtSumDiffShape[jj1]/JetPt);
1343 fProIntJetShapeAUE[k] ->Fill(tmp1R,PtSumIntShape[jj1]/JetPt);
1344 }//if
1345 }//k loop
1346 tmp1R +=0.1;
1347 }//jj1 loop
1348 //----------------------//
1349 Float_t PtSum = 0;
1350 Int_t NtrkSum = 0;
1351 Bool_t iflagPtSum = kFALSE;
1352 Bool_t iflagNtrkSum = kFALSE;
1353 TMath::Sort(nTracksUE,delRA,index,0);
1354 for(Int_t ii=0; ii<nTracksUE; ii++){
1355 NtrkSum ++;
1356 PtSum += trackPtA[index[ii]];
1357 /*
1358 cout << index[ii] << "\t" <<
1359 delR[ii] << "\t" <<
1360 delEta[ii]<< "\t" <<
1361 delPhi[ii]<< "\t" <<
1362 trackPt[ii]<< "\t" <<
1363 trackEta[ii]<< "\t" <<
1364 trackPhi[ii]<< "\t DelR " <<
1365 delR[index[ii]] << endl;
1366 */
1367 if(!iflagNtrkSum){
1368 if((Float_t)NtrkSum/(Float_t)nTracksUE > 0.79){
1369 delRNtrkSum80pc = delRA[index[ii]];
1370 iflagNtrkSum = kTRUE;
1371 }//if NtrkSum
1372 }//if iflag
1373 if(!iflagPtSum){
1374 if(PtSum/JetPt >= 0.8000){
1375 delRPtSum80pc = delRA[index[ii]];
1376 iflagPtSum = kTRUE;
1377 }//if PtSum
1378 }//if iflag
1379 }//track loop 2nd
1380 delete [] index;
1381 delete [] delRA;
1382 delete [] delEtaA;
1383 delete [] delPhiA;
1384 delete [] trackPtA;
1385 delete [] trackEtaA;
1386 delete [] trackPhiA;
1387 fh2DelR80pcNchUE ->Fill(JetPt,delRNtrkSum80pc);
1388 fProDelR80pcNchUE->Fill(JetPt,delRNtrkSum80pc);
1389 fh2DelR80pcPtUE ->Fill(JetPt,delRPtSum80pc);
1390 fProDelR80pcPtUE ->Fill(JetPt,delRPtSum80pc);
1391 for(Int_t ibin=0; ibin<50; ibin++){
1392 Float_t iR = 0.02*ibin + 0.01;
1393 fh3PtDelRNchSumUE->Fill(JetPt,iR,NchSumA[ibin]);
1394 fh3PtDelRPtSumUE ->Fill(JetPt,iR,PtSumA[ibin]);
1395 Float_t jetPtMin=20.0; Float_t jetPtMax=30.0;
1396 for(Int_t k=0; k<13; k++){
1397 if(k==0) {jetPtMin=20.0;jetPtMax=25.0;}
1398 if(k==1) {jetPtMin=25.0;jetPtMax=30.0;}
1399 if(k==2) {jetPtMin=20.0;jetPtMax=30.0;}
1400 if(k==3) {jetPtMin=30.0;jetPtMax=35.0;}
1401 if(k==4) {jetPtMin=35.0;jetPtMax=40.0;}
1402 if(k==5) {jetPtMin=30.0;jetPtMax=40.0;}
1403 if(k==6) {jetPtMin=40.0;jetPtMax=50.0;}
1404 if(k==7) {jetPtMin=50.0;jetPtMax=60.0;}
1405 if(k==8) {jetPtMin=40.0;jetPtMax=60.0;}
1406 if(k==9) {jetPtMin=60.0;jetPtMax=70.0;}
1407 if(k==10){jetPtMin=70.0;jetPtMax=80.0;}
1408 if(k==11){jetPtMin=60.0;jetPtMax=80.0;}
1409 if(k==12){jetPtMin=80.0;jetPtMax=100.0;}
1410 if(JetPt>jetPtMin && JetPt<jetPtMax){
1411 fProDelRNchSumUE[k]->Fill(iR,NchSumA[ibin]);
1412 fProDelRPtSumUE[k]->Fill(iR,PtSumA[ibin]);
1413 }//if
1414 }//k loop
1415 }//ibin loop
1416 }//if(jet)
1417}//FillJetShapeUE()
1418//_________________________________________________________________________________//
59a5935a 1419
fed89a49 1420void AliAnalysisTaskJetProperties::GetTracksTiltedwrpJetAxis(Float_t alpha, TList* inputlist, TList* outputlist, const AliAODJet* jet, Double_t radius,Double_t& sumPt){
1421 //This part is inherited from AliAnalysisTaskFragmentationFunction.cxx
1422 // List of tracks in cone perpendicular to the jet azimuthal direction
1423 Double_t jetMom[3];
1424 jet->PxPyPz(jetMom);
1425
1426 TVector3 jet3mom(jetMom);
1427 // Rotate phi and keep eta unchanged
1428 Double_t etaTilted = jet3mom.Eta();//no change in eta
1429 Double_t phiTilted = TVector2::Phi_0_2pi(jet3mom.Phi()) + alpha;//rotate phi by alpha
1430 if(phiTilted > 2*TMath::Pi()) phiTilted = phiTilted - 2*TMath::Pi();
1431
1432 for (Int_t itrack=0; itrack<inputlist->GetSize(); itrack++){
1433 if(fTrackType==kTrackUndef)continue;
1434 Double_t trackMom[3];
1435
1436 if(fTrackType==kTrackAOD){
1437 AliAODTrack *trackaod = dynamic_cast<AliAODTrack*>(inputlist->At(itrack));
1438 if(!trackaod)continue;
1439 trackaod->PxPyPz(trackMom);
1440 TVector3 track3mom(trackMom);
1441
1442 Double_t deta = track3mom.Eta() - etaTilted;
1443 Double_t dphi = TMath::Abs(track3mom.Phi() - phiTilted);
1444 if (dphi > TMath::Pi()) dphi = 2. * TMath::Pi() - dphi;
1445 Double_t dR = TMath::Sqrt(deta * deta + dphi * dphi);
1446
1447 if(dR<=radius){
1448 outputlist->Add(trackaod);
1449 sumPt += trackaod->Pt();
1450 if(fDebug > 100) printf("GetTracksTiltewrpJetAxis itrack, trackPt %d, %f\n",itrack,track3mom.Pt());
1451 }//if dR<=radius
1452 }//if kTrackAOD
1453 else if(fTrackType==kTrackAODMC){
1454 AliAODMCParticle* trackmc = dynamic_cast<AliAODMCParticle*>(inputlist->At(itrack));
1455 if(!trackmc)continue;
1456 trackmc->PxPyPz(trackMom);
1457 TVector3 track3mom(trackMom);
1458
1459 Double_t deta = track3mom.Eta() - etaTilted;
1460 Double_t dphi = TMath::Abs(track3mom.Phi() - phiTilted);
1461 if (dphi > TMath::Pi()) dphi = 2. * TMath::Pi() - dphi;
1462 Double_t dR = TMath::Sqrt(deta * deta + dphi * dphi);
1463
1464 if(dR<=radius){
1465 outputlist->Add(trackmc);
1466 sumPt += trackmc->Pt();
1467 if(fDebug > 100) printf("GetTracksTiltewrpJetAxis itrack, trackPt %d, %f\n",itrack,track3mom.Pt());
1468 }//if dR<=radius
1469 }//if kTrackAODMC
1470 else if(fTrackType==kTrackKine){
1471 AliVParticle* trackkine = dynamic_cast<AliVParticle*>(inputlist->At(itrack));
1472 if(!trackkine)continue;
1473 trackkine->PxPyPz(trackMom);
1474 TVector3 track3mom(trackMom);
1475
1476 Double_t deta = track3mom.Eta() - etaTilted;
1477 Double_t dphi = TMath::Abs(track3mom.Phi() - phiTilted);
1478 if (dphi > TMath::Pi()) dphi = 2. * TMath::Pi() - dphi;
1479 Double_t dR = TMath::Sqrt(deta * deta + dphi * dphi);
1480
1481 if(dR<=radius){
1482 outputlist->Add(trackkine);
1483 sumPt += trackkine->Pt();
1484 if(fDebug > 100) printf("GetTracksTiltewrpJetAxis itrack, trackPt %d, %f\n",itrack,track3mom.Pt());
1485 }//if dR<=radius
1486 }//kTrackKine
1487 }//itrack
1488}//initial loop
1489//-----------------------------------------------//
1490Int_t AliAnalysisTaskJetProperties::GetListOfTracks(TList *list, Int_t type)
1491{
1492 //This part is inherited from AliAnalysisTaskFragmentationFunction.cxx (and modified)
1493 // fill list of tracks selected according to type
1494
1495 if(fDebug > 2) Printf("%s:%d Selecting tracks with %d", (char*)__FILE__,__LINE__,type);
1496
1497 if(!list){
1498 if(fDebug>1) Printf("%s:%d no input list", (char*)__FILE__,__LINE__);
1499 return -1;
1500 }
1501
1502 if(!fAOD) return -1;
1503
1504 if(!fAOD->GetTracks()) return 0;
1505
1506 if(type==kTrackUndef) return 0;
1507
1508 Int_t iCount = 0;
1509 if(type==kTrackAOD){
1510 // all rec. tracks, esd filter mask, within acceptance
1511 for(Int_t it=0; it<fAOD->GetNumberOfTracks(); ++it){
1512 AliAODTrack *tr = fAOD->GetTrack(it);
1513 if(!tr)continue;
1514 if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask)))continue;//selecting filtermask
1515 if(tr->Eta() < fTrackEtaMin || tr->Eta() > fTrackEtaMax) continue;
1516 if(tr->Pt() < fTrackPtCut) continue;
1517 list->Add(tr);
1518 iCount++;
1519 }//track loop
1520 }//if type
1521 else if (type==kTrackKine){
1522 // kine particles, primaries, charged only within acceptance
1523 if(!fMCEvent) return iCount;
1524 for(Int_t it=0; it<fMCEvent->GetNumberOfTracks(); ++it){
1525 if(!fMCEvent->IsPhysicalPrimary(it))continue;//selecting only primaries
1526 AliMCParticle* part = (AliMCParticle*) fMCEvent->GetTrack(it);
1527 if(!part)continue;
1528 if(part->Particle()->GetPDG()->Charge()==0)continue;//removing neutrals
1529 if(part->Eta() < fTrackEtaMin || part->Eta() > fTrackEtaMax) continue;//eta cut
1530 if(part->Pt() < fTrackPtCut)continue;//pt cut
1531 list->Add(part);
1532 iCount++;
1533 }//track loop
1534 }//if type
1535 else if (type==kTrackAODMC) {
1536 // MC particles (from AOD), physical primaries, charged only within acceptance
1537 if(!fAOD) return -1;
1538
1539 TClonesArray *tca = dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName()));
1540 if(!tca)return iCount;
1541
1542 for(int it=0; it<tca->GetEntriesFast(); ++it){
1543 AliAODMCParticle *part = dynamic_cast<AliAODMCParticle*>(tca->At(it));
1544 if(!part)continue;
1545 if(!part->IsPhysicalPrimary())continue;
1546 if(part->Charge()==0) continue;
1547 if(part->Eta() > fTrackEtaMax || part->Eta() < fTrackEtaMin)continue;
1548 if(part->Pt() < fTrackPtCut) continue;
1549 list->Add(part);
1550 iCount++;
1551 }//track loop
1552 }//if type
1553
1554 list->Sort();
1555 return iCount;
1556}
1557//_______________________________________________________________________________