]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGJE/AliAnalysisTaskJetProperties.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
take only good quality tracks as trigger
[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
59 AliAnalysisTaskJetProperties::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 ,fJetList(0)
82 ,fTrackList(0)
83 ,fCommonHistList(0)
84 ,fh1EvtSelection(0)
85 ,fh1VertexNContributors(0)
86 ,fh1VertexZ(0)
87 ,fh1Xsec(0)
88 ,fh1Trials(0)
89 ,fh1PtHard(0)
90 ,fh1PtHardTrials(0)
91 ,fh2EtaJet(0)
92 ,fh2PhiJet(0)
93 ,fh2PtJet(0)
94 ,fh1PtJet(0)
95 ,fh2NtracksJet(0)
96 ,fProNtracksJet(0)
97 ,fh2EtaTrack(0)
98 ,fh2PhiTrack(0)
99 ,fh2PtTrack(0)
100 ,fh2FF(0)
101 ,fh2DelEta(0)
102 ,fh2DelPhi(0)
103 ,fh2DelR(0)
104 ,fh1PtLeadingJet(0)
105 ,fh2NtracksLeadingJet(0)
106 ,fProNtracksLeadingJet(0)
107 ,fh2DelR80pcNch(0)
108 ,fProDelR80pcNch(0)
109 ,fh2DelR80pcPt(0)
110 ,fProDelR80pcPt(0)
111 ,fh2AreaCh(0)
112 ,fProAreaCh(0)
113 ,fh3PtDelRNchSum(0)
114 ,fh3PtDelRPtSum(0)
efabdd32 115 ,fProDiffJetShape(0)
116 ,fProIntJetShape(0)
59a5935a 117{
118 for(Int_t ii=0; ii<5; ii++){
efabdd32 119 fProDelRNchSum[ii] = NULL;
120 fProDelRPtSum[ii] = NULL;
121 fProDiffJetShapeA[ii] = NULL;
122 fProIntJetShapeA[ii] = NULL;
59a5935a 123 }//ii loop
124 // default constructor
125}
126//_________________________________________________________________________________//
127
128AliAnalysisTaskJetProperties::AliAnalysisTaskJetProperties(const char *name)
129 : AliAnalysisTaskSE(name)
130 ,fESD(0)
131 ,fAOD(0)
132 ,fAODJets(0)
133 ,fAODExtension(0)
134 ,fNonStdFile("")
135 ,fBranchJets("jets")
136 ,fTrackType(kTrackAOD)
137 ,fJetRejectType(0)
138 ,fUseAODInputJets(kTRUE)
139 ,fFilterMask(0)
140 ,fUsePhysicsSelection(kTRUE)
141 ,fMaxVertexZ(10)
142 ,fNContributors(2)
143 ,fTrackPtCut(0)
144 ,fTrackEtaMin(0)
145 ,fTrackEtaMax(0)
146 ,fJetPtCut(0)
147 ,fJetEtaMin(0)
148 ,fJetEtaMax(0)
149 ,fAvgTrials(0)
150 ,fJetList(0)
151 ,fTrackList(0)
152 ,fCommonHistList(0)
153 ,fh1EvtSelection(0)
154 ,fh1VertexNContributors(0)
155 ,fh1VertexZ(0)
156 ,fh1Xsec(0)
157 ,fh1Trials(0)
158 ,fh1PtHard(0)
159 ,fh1PtHardTrials(0)
160 ,fh2EtaJet(0)
161 ,fh2PhiJet(0)
162 ,fh2PtJet(0)
163 ,fh1PtJet(0)
164 ,fh2NtracksJet(0)
165 ,fProNtracksJet(0)
166 ,fh2EtaTrack(0)
167 ,fh2PhiTrack(0)
168 ,fh2PtTrack(0)
169 ,fh2FF(0)
170 ,fh2DelEta(0)
171 ,fh2DelPhi(0)
172 ,fh2DelR(0)
173 ,fh1PtLeadingJet(0)
174 ,fh2NtracksLeadingJet(0)
175 ,fProNtracksLeadingJet(0)
176 ,fh2DelR80pcNch(0)
177 ,fProDelR80pcNch(0)
178 ,fh2DelR80pcPt(0)
179 ,fProDelR80pcPt(0)
180 ,fh2AreaCh(0)
181 ,fProAreaCh(0)
182 ,fh3PtDelRNchSum(0)
183 ,fh3PtDelRPtSum(0)
efabdd32 184 ,fProDiffJetShape(0)
185 ,fProIntJetShape(0)
59a5935a 186{
187 for(Int_t ii=0; ii<5; ii++){
efabdd32 188 fProDelRNchSum[ii] = NULL;
189 fProDelRPtSum[ii] = NULL;
190 fProDiffJetShapeA[ii] = NULL;
191 fProIntJetShapeA[ii] = NULL;
59a5935a 192 }//ii loop
193 // constructor
194 DefineOutput(1,TList::Class());
195}
196//_________________________________________________________________________________//
197
198AliAnalysisTaskJetProperties::~AliAnalysisTaskJetProperties()
199{
200 // destructor
201 if(fJetList) delete fJetList;
202 if(fTrackList) delete fTrackList;
203}
204//_________________________________________________________________________________//
205
206Bool_t AliAnalysisTaskJetProperties::Notify()
207{
208 //
209 // Implemented Notify() to read the cross sections
210 // and number of trials from pyxsec.root
211 // (taken from AliAnalysisTaskJetSpectrum2)
212 //
213 if(fDebug > 1) Printf("AliAnalysisTaskJetProperties::Notify()");
214
215 TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
216 Float_t xsection = 0;
217 Float_t ftrials = 1;
218
219 fAvgTrials = 1;
220 if(tree){
221 TFile *curfile = tree->GetCurrentFile();
222 if (!curfile) {
223 Error("Notify","No current file");
224 return kFALSE;
225 }
226 if(!fh1Xsec||!fh1Trials){
227 Printf("%s%d No Histogram fh1Xsec",(char*)__FILE__,__LINE__);
228 return kFALSE;
229 }
230 AliAnalysisHelperJetTasks::PythiaInfoFromFile(curfile->GetName(),xsection,ftrials);
231 fh1Xsec->Fill("<#sigma>",xsection);
232 // construct a poor man average trials
233 Float_t nEntries = (Float_t)tree->GetTree()->GetEntries();
234 if(ftrials>=nEntries && nEntries>0.)fAvgTrials = ftrials/nEntries;
235 }
236 return kTRUE;
237}
238//_________________________________________________________________________________//
239
240void AliAnalysisTaskJetProperties::UserCreateOutputObjects()
241{
242 //(Here, event selection part is taken from AliAnalysisTaskFramentationFunction)
243 // create output objects
244 if(fDebug > 1) Printf("AliAnalysisTaskJetProperties::UserCreateOutputObjects()");
245
246 // create list of tracks and jets
247 fJetList = new TList();
248 fJetList->SetOwner(kFALSE);
249 fTrackList = new TList();
250 fTrackList->SetOwner(kFALSE);
251
252 // Create histograms / output container
253 OpenFile(1);
254 fCommonHistList = new TList();
255 fCommonHistList->SetOwner();
256
257 Bool_t oldStatus = TH1::AddDirectoryStatus();
258 TH1::AddDirectory(kFALSE);
259 // Histograms/TProfiles
260 fh1EvtSelection = new TH1F("fh1EvtSelection", "Event Selection", 6, -0.5, 5.5);
261 fh1EvtSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
262 fh1EvtSelection->GetXaxis()->SetBinLabel(2,"event selection: rejected");
263 fh1EvtSelection->GetXaxis()->SetBinLabel(3,"event class: rejected");
264 fh1EvtSelection->GetXaxis()->SetBinLabel(4,"vertex Ncontr: rejected");
265 fh1EvtSelection->GetXaxis()->SetBinLabel(5,"vertex z: rejected");
266 fh1EvtSelection->GetXaxis()->SetBinLabel(6,"vertex type: rejected");
267
268
269 fh1VertexNContributors = new TH1F("fh1VertexNContributors", "Vertex N contributors", 2500,-.5, 2499.5);
270 fh1VertexZ = new TH1F("fh1VertexZ", "Vertex z distribution", 30, -15., 15.);
271 fh1Xsec = new TProfile("fh1Xsec","xsec from pyxsec.root",1,0,1);
272 fh1Xsec->GetXaxis()->SetBinLabel(1,"<#sigma>");
273 fh1Trials = new TH1F("fh1Trials","trials from pyxsec.root",1,0,1);
274 fh1Trials->GetXaxis()->SetBinLabel(1,"#sum{ntrials}");
275 fh1PtHard = new TH1F("fh1PtHard","PYTHIA Pt hard;p_{T,hard}",350,-.5,349.5);
276 fh1PtHardTrials = new TH1F("fh1PtHardTrials","PYTHIA Pt hard weight with trials;p_{T,hard}",350,-.5,349.5);
277
278
279 Int_t kNbinsPtSlice=20; Float_t xMinPtSlice=0.0; Float_t xMaxPtSlice=100.0;
280 Int_t kNbinsEta=40; Float_t xMinEta=-2.0; Float_t xMaxEta=2.0;
281 Int_t kNbinsPhi=90; Float_t xMinPhi=0.0; Float_t xMaxPhi=TMath::TwoPi();
282 Int_t kNbinsPt=125; Float_t xMinPt=0.0; Float_t xMaxPt=250.0;
283 Int_t kNbinsNtracks=50; Float_t xMinNtracks=0.0; Float_t xMaxNtracks=50.0;
284 Int_t kNbinsFF=80; Float_t xMinFF=0.0; Float_t xMaxFF=2.0;
285 Int_t kNbinsDelR1D=50; Float_t xMinDelR1D=0.0; Float_t xMaxDelR1D=1.0;
286
287
288 fh2EtaJet = new TH2F("EtaJet","EtaJet",
289 kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
290 kNbinsEta, xMinEta, xMaxEta);
291 fh2PhiJet = new TH2F("PhiJet","PhiJet",
292 kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
293 kNbinsPhi, xMinPhi, xMaxPhi);
294 fh2PtJet = new TH2F("PtJet","PtJet",
295 kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
296 kNbinsPt, xMinPt, xMaxPt);
297 fh1PtJet = new TH1F("PtJet1D","PtJet1D",
298 kNbinsPt, xMinPt, xMaxPt);
299 fh2NtracksJet = new TH2F("NtracksJet","NtracksJet",
300 kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
301 kNbinsNtracks, xMinNtracks, xMaxNtracks);
302 fProNtracksJet = new TProfile("AvgNoOfTracksJet","AvgNoOfTracksJet",
303 kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
304 xMinNtracks, xMaxNtracks);
305 fh2EtaTrack = new TH2F("EtaTrack","EtaTrack",
306 kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
307 kNbinsEta, xMinEta, xMaxEta);
308 fh2PhiTrack = new TH2F("PhiTrack","PhiTrack",
309 kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
310 kNbinsPhi, xMinPhi, xMaxPhi);
311 fh2PtTrack = new TH2F("PtTrack","PtTrack",
312 kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
313 kNbinsPt, xMinPt, xMaxPt);
314 fh2FF = new TH2F("FF","FF",
315 kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
316 kNbinsFF, xMinFF, xMaxFF);
317 fh2DelEta = new TH2F("DelEta","DelEta",
318 kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
319 kNbinsEta, xMinEta, xMaxEta);
320 fh2DelPhi = new TH2F("DelPhi","DelPhi",
321 kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
322 kNbinsPhi, xMinPhi, xMaxPhi);
323 fh2DelR = new TH2F("DelR","DelR",
324 kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
325 kNbinsDelR1D, xMinDelR1D, xMaxDelR1D);
326
327
328
329 Int_t kNbinsDelR=100; Float_t xMinDelR=0.0; Float_t xMaxDelR=1.0;
330 Int_t kNbinsPtSliceJS=100; Float_t xMinPtSliceJS=0.0; Float_t xMaxPtSliceJS=250.0;
331
332 fh1PtLeadingJet = new TH1F("PtLeadingJet","PtLeadingJet",
333 kNbinsPt, xMinPt, xMaxPt);
334 fh2NtracksLeadingJet = new TH2F("NtracksLeadingJet","NtracksLeadingJet",
335 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
336 kNbinsNtracks, xMinNtracks, xMaxNtracks);
337 fProNtracksLeadingJet = new TProfile("AvgNoOfTracksLeadingJet","AvgNoOfTracksLeadingJet",
338 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
339 xMinNtracks, xMaxNtracks);
340 fh2DelR80pcNch = new TH2F("delR80pcNch","delR80pcNch",
341 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
342 kNbinsDelR, xMinDelR, xMaxDelR);
343 fProDelR80pcNch = new TProfile("AvgdelR80pcNch","AvgdelR80pcNch",
344 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
345 xMinDelR, xMaxDelR);
346 fh2DelR80pcPt = new TH2F("delR80pcPt","delR80pcPt",
347 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
348 kNbinsDelR, xMinDelR, xMaxDelR);
349 fProDelR80pcPt = new TProfile("AvgdelR80pcPt","AvgdelR80pcPt",
350 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
351 xMinDelR, xMaxDelR);
352 fh2AreaCh = new TH2F("Area","Area",
353 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
354 72, 0.0, TMath::Pi());
355 fProAreaCh = new TProfile("AvgArea","AvgArea",
356 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
357 0.0, TMath::Pi());
358 fh3PtDelRNchSum = new TH3F("PtDelRNchSum","PtDelRNchSum",
359 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
360 kNbinsDelR1D, xMinDelR1D, xMaxDelR1D,
361 kNbinsNtracks, xMinNtracks, xMaxNtracks);
362 fh3PtDelRPtSum = new TH3F("PtDelRPtSum","PtDelRPtSum",
363 kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
364 kNbinsDelR1D, xMinDelR1D, xMaxDelR1D,
365 kNbinsPt, xMinPt, xMaxPt);
efabdd32 366 fProDiffJetShape = new TProfile("DiffJetShape","DiffJetShape",
367 10,0.0,1.0,0.0,250.0);
368 fProIntJetShape = new TProfile("IntJetShape","IntJetShape",
369 10,0.0,1.0,0.0,250.0);
370
59a5935a 371 TString title;
372 for(Int_t ii=0; ii<5; ii++){
373 if(ii==0)title = "_JetPt20to30";
374 if(ii==1)title = "_JetPt30to40";
375 if(ii==2)title = "_JetPt40to60";
376 if(ii==3)title = "_JetPt60to80";
377 if(ii==4)title = "_JetPt80to100";
378 fProDelRNchSum[ii] = new TProfile(Form("AvgNchSumDelR%s",title.Data()),Form("AvgNchSumDelR%s",title.Data()),
379 kNbinsDelR1D, xMinDelR1D, xMaxDelR1D,
380 xMinNtracks, xMaxNtracks);
381
382 fProDelRPtSum[ii] = new TProfile(Form("AvgPtSumDelR%s",title.Data()),Form("AvgPtSumDelR%s",title.Data()),
383 kNbinsDelR1D, xMinDelR1D, xMaxDelR1D,
384 xMinPt, xMaxPt);
385 fProDelRNchSum[ii] ->GetXaxis()->SetTitle("R");
386 fProDelRNchSum[ii] ->GetYaxis()->SetTitle("<NchSum>");
387 fProDelRPtSum[ii] ->GetXaxis()->SetTitle("R");
388 fProDelRPtSum[ii] ->GetYaxis()->SetTitle("<PtSum>");
efabdd32 389
390 fProDiffJetShapeA[ii] = new TProfile(Form("DiffJetShape%s",title.Data()),Form("DiffJetShape%s",title.Data()),
391 10,0.0,1.0,0.0,250.0);
392 fProIntJetShapeA[ii] = new TProfile(Form("IntJetShape%s",title.Data()),Form("IntJetShape%s",title.Data()),
393 10,0.0,1.0,0.0,250.0);
394
395 fProDiffJetShapeA[ii]->GetXaxis()->SetTitle("R");
396 fProDiffJetShapeA[ii]->GetYaxis()->SetTitle("Diff jet shape");
397 fProIntJetShapeA[ii]->GetXaxis()->SetTitle("R");
398 fProIntJetShapeA[ii]->GetYaxis()->SetTitle("Integrated jet shape");
399
59a5935a 400 fCommonHistList->Add(fProDelRNchSum[ii]);
401 fCommonHistList->Add(fProDelRPtSum[ii]);
efabdd32 402 fCommonHistList->Add(fProDiffJetShapeA[ii]);
403 fCommonHistList->Add(fProIntJetShapeA[ii]);
404
405}//ii loop
59a5935a 406
407 fh2EtaJet ->GetXaxis()->SetTitle("JetPt"); fh2EtaJet ->GetYaxis()->SetTitle("JetEta");
408 fh2PhiJet ->GetXaxis()->SetTitle("JetPt"); fh2PhiJet ->GetYaxis()->SetTitle("JetPhi");
409 fh2PtJet ->GetXaxis()->SetTitle("JetPt"); fh2PtJet ->GetYaxis()->SetTitle("JetPt");
410 fh1PtJet ->GetXaxis()->SetTitle("JetPt"); fh1PtJet ->GetYaxis()->SetTitle("#jets");
411 fh2NtracksJet ->GetXaxis()->SetTitle("JetPt"); fh2NtracksJet ->GetYaxis()->SetTitle("#tracks");
412 fProNtracksJet->GetXaxis()->SetTitle("JetPt"); fProNtracksJet->GetYaxis()->SetTitle("AgvNtracks");
413 fh2EtaTrack ->GetXaxis()->SetTitle("JetPt"); fh2EtaTrack ->GetYaxis()->SetTitle("TrackEta");
414 fh2PhiTrack ->GetXaxis()->SetTitle("JetPt"); fh2PhiTrack ->GetYaxis()->SetTitle("TrackPhi");
415 fh2PtTrack ->GetXaxis()->SetTitle("JetPt"); fh2PtTrack ->GetYaxis()->SetTitle("TrackPt");
416 fh2FF ->GetXaxis()->SetTitle("JetPt"); fh2FF ->GetYaxis()->SetTitle("FF");
417 fh2DelEta ->GetXaxis()->SetTitle("JetPt"); fh2DelEta ->GetYaxis()->SetTitle("DelEta");
418 fh2DelPhi ->GetXaxis()->SetTitle("JetPt"); fh2DelPhi ->GetYaxis()->SetTitle("DelPhi");
419 fh2DelR ->GetXaxis()->SetTitle("JetPt"); fh2DelR ->GetYaxis()->SetTitle("DelR");
420
421 fh1PtLeadingJet ->GetXaxis()->SetTitle("JetPt"); fh1PtLeadingJet ->GetYaxis()->SetTitle("#leading jets");
422 fh2NtracksLeadingJet ->GetXaxis()->SetTitle("JetPt"); fh2NtracksLeadingJet ->GetYaxis()->SetTitle("#tracks leading jet");
423 fProNtracksLeadingJet ->GetXaxis()->SetTitle("JetPt"); fProNtracksLeadingJet ->GetYaxis()->SetTitle("AvgNtracks leading jet");
424 fh2DelR80pcNch ->GetXaxis()->SetTitle("JetPt"); fh2DelR80pcNch ->GetYaxis()->SetTitle("R containing 80% of tracks");
425 fProDelR80pcNch ->GetXaxis()->SetTitle("JetPt"); fProDelR80pcNch ->GetYaxis()->SetTitle("<R> containing 80% of tracks");
426 fh2DelR80pcPt ->GetXaxis()->SetTitle("JetPt"); fh2DelR80pcPt ->GetYaxis()->SetTitle("R containing 80% of pT");
427 fProDelR80pcPt ->GetXaxis()->SetTitle("JetPt"); fProDelR80pcPt ->GetYaxis()->SetTitle("<R> containing 80% of pT");
428 fh2AreaCh ->GetXaxis()->SetTitle("JetPt"); fh2AreaCh ->GetYaxis()->SetTitle("Jet area");
429 fProAreaCh ->GetXaxis()->SetTitle("JetPt"); fProAreaCh ->GetYaxis()->SetTitle("<jet area>");
430 fh3PtDelRNchSum ->GetXaxis()->SetTitle("JetPt"); fh3PtDelRNchSum ->GetYaxis()->SetTitle("R"); fh3PtDelRNchSum->GetZaxis()->SetTitle("NchSum");
431 fh3PtDelRPtSum ->GetXaxis()->SetTitle("JetPt"); fh3PtDelRPtSum ->GetYaxis()->SetTitle("R"); fh3PtDelRPtSum ->GetZaxis()->SetTitle("PtSum");
efabdd32 432 fProDiffJetShape->GetXaxis()->SetTitle("R");
433 fProDiffJetShape->GetYaxis()->SetTitle("Diff jet shape");
434 fProIntJetShape->GetXaxis()->SetTitle("R");
435 fProIntJetShape->GetYaxis()->SetTitle("Integrated jet shape");
59a5935a 436
437 fCommonHistList->Add(fh1EvtSelection);
438 fCommonHistList->Add(fh1VertexNContributors);
439 fCommonHistList->Add(fh1VertexZ);
440 fCommonHistList->Add(fh1Xsec);
441 fCommonHistList->Add(fh1Trials);
442 fCommonHistList->Add(fh1PtHard);
443 fCommonHistList->Add(fh1PtHardTrials);
444 fCommonHistList->Add(fh2EtaJet);
445 fCommonHistList->Add(fh2PhiJet);
446 fCommonHistList->Add(fh2PtJet);
447 fCommonHistList->Add(fh1PtJet);
448 fCommonHistList->Add(fh2NtracksJet);
449 fCommonHistList->Add(fProNtracksJet);
450 fCommonHistList->Add(fh2EtaTrack);
451 fCommonHistList->Add(fh2PhiTrack);
452 fCommonHistList->Add(fh2PtTrack);
453 fCommonHistList->Add(fh2FF);
454 fCommonHistList->Add(fh2DelEta);
455 fCommonHistList->Add(fh2DelPhi);
456 fCommonHistList->Add(fh2DelR);
457
458 fCommonHistList->Add(fh1PtLeadingJet);
459 fCommonHistList->Add(fh2NtracksLeadingJet);
460 fCommonHistList->Add(fProNtracksLeadingJet);
461 fCommonHistList->Add(fh2DelR80pcNch);
462 fCommonHistList->Add(fProDelR80pcNch);
463 fCommonHistList->Add(fh2DelR80pcPt);
464 fCommonHistList->Add(fProDelR80pcPt);
465 fCommonHistList->Add(fh2AreaCh);
466 fCommonHistList->Add(fProAreaCh);
467 fCommonHistList->Add(fh3PtDelRNchSum);
468 fCommonHistList->Add(fh3PtDelRPtSum);
efabdd32 469 fCommonHistList->Add(fProDiffJetShape);
470 fCommonHistList->Add(fProIntJetShape);
59a5935a 471 // =========== Switch on Sumw2 for all histos ===========
472 for (Int_t i=0; i<fCommonHistList->GetEntries(); ++i){
473 TH1 *h1 = dynamic_cast<TH1*>(fCommonHistList->At(i));
474 if (h1) h1->Sumw2();
475 else{
476 TProfile *hPro = dynamic_cast<TProfile*>(fCommonHistList->At(i));
477 if(hPro) hPro->Sumw2();
478 }
479 }
480
481 TH1::AddDirectory(oldStatus);
482 PostData(1, fCommonHistList);
483}
484//_________________________________________________________________________________//
485
486void AliAnalysisTaskJetProperties::Init()
487{
488 // Initialization
489 if(fDebug > 1) Printf("AliAnalysisTaskJetProperties::Init()");
490
491}
492//_________________________________________________________________________________//
493
494void AliAnalysisTaskJetProperties::UserExec(Option_t *)
495{
496 // Main loop
497 // Called for each event
498 if(fDebug > 1) Printf("AliAnalysisTaskJetProperties::UserExec()");
499
500 if(fDebug > 1) Printf("Analysis event #%5d", (Int_t) fEntry);
501 // Trigger selection
502
503 AliInputEventHandler* inputHandler = (AliInputEventHandler*)
504 ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
505 if(!(inputHandler->IsEventSelected() & AliVEvent::kMB)){
506 if(inputHandler->InheritsFrom("AliESDInputHandler") && fUsePhysicsSelection){ // PhysicsSelection only with ESD input
507 fh1EvtSelection->Fill(1.);
508 if (fDebug > 1 ) Printf(" Trigger Selection: event REJECTED ... ");
509 PostData(1, fCommonHistList);
510 return;
511 }
512 }
513
514 fESD = dynamic_cast<AliESDEvent*>(InputEvent());
515 if(!fESD){
516 if(fDebug>3) Printf("%s:%d ESDEvent not found in the input", (char*)__FILE__,__LINE__);
517 }
518
519 fMCEvent = MCEvent();
520 if(!fMCEvent){
521 if(fDebug>3) Printf("%s:%d MCEvent not found in the input", (char*)__FILE__,__LINE__);
522 }
523
524 // get AOD event from input/ouput
525 TObject* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
526 if( handler && handler->InheritsFrom("AliAODInputHandler") ) {
527 fAOD = ((AliAODInputHandler*)handler)->GetEvent();
528 if(fUseAODInputJets) fAODJets = fAOD;
529 if (fDebug > 1) Printf("%s:%d AOD event from input", (char*)__FILE__,__LINE__);
530 }
531 else {
532 handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
533 if( handler && handler->InheritsFrom("AliAODHandler") ) {
534 fAOD = ((AliAODHandler*)handler)->GetAOD();
535 fAODJets = fAOD;
536 if (fDebug > 1) Printf("%s:%d AOD event from output", (char*)__FILE__,__LINE__);
537 }
538 }
539
540 if(!fAODJets && !fUseAODInputJets){ // case we have AOD in input & output and want jets from output
541 TObject* outHandler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
542 if( outHandler && outHandler->InheritsFrom("AliAODHandler") ) {
543 fAODJets = ((AliAODHandler*)outHandler)->GetAOD();
544 if (fDebug > 1) Printf("%s:%d jets from output AOD", (char*)__FILE__,__LINE__);
545 }
546 }
547
548 if(fNonStdFile.Length()!=0){
549 // case we have an AOD extension - fetch the jets from the extended output
550
551 AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
552 fAODExtension = (aodH?aodH->GetExtension(fNonStdFile.Data()):0);
553 if(!fAODExtension){
554 if(fDebug>1)Printf("AODExtension not found for %s",fNonStdFile.Data());
555 }
556 }
557
558 if(!fAOD){
559 Printf("%s:%d AODEvent not found", (char*)__FILE__,__LINE__);
560 return;
561 }
562 if(!fAODJets){
563 Printf("%s:%d AODEvent with jet branch not found", (char*)__FILE__,__LINE__);
564 return;
565 }
566
567
568 // *** vertex cut ***
569 AliAODVertex* primVtx = fAOD->GetPrimaryVertex();
570 Int_t nTracksPrim = primVtx->GetNContributors();
571 fh1VertexNContributors->Fill(nTracksPrim);
572
573 if (fDebug > 1) Printf("%s:%d primary vertex selection: %d", (char*)__FILE__,__LINE__,nTracksPrim);
574 //if(!nTracksPrim){
575 if(nTracksPrim<fNContributors){
576 if (fDebug > 1) Printf("%s:%d primary vertex selection: event REJECTED...",(char*)__FILE__,__LINE__);
577 fh1EvtSelection->Fill(3.);
578 PostData(1, fCommonHistList);
579 return;
580 }
581
582 fh1VertexZ->Fill(primVtx->GetZ());
583
584 if(TMath::Abs(primVtx->GetZ())>fMaxVertexZ){
585 if (fDebug > 1) Printf("%s:%d primary vertex z = %f: event REJECTED...",(char*)__FILE__,__LINE__,primVtx->GetZ());
586 fh1EvtSelection->Fill(4.);
587 PostData(1, fCommonHistList);
588 return;
589 }
590
591 TString primVtxName(primVtx->GetName());
592
593 if(primVtxName.CompareTo("TPCVertex",TString::kIgnoreCase) == 1){
594 if (fDebug > 1) Printf("%s:%d primary vertex selection: TPC vertex, event REJECTED...",(char*)__FILE__,__LINE__);
595 fh1EvtSelection->Fill(5.);
596 PostData(1, fCommonHistList);
597 return;
598 }
599
600 if (fDebug > 1) Printf("%s:%d event ACCEPTED ...",(char*)__FILE__,__LINE__);
601 fh1EvtSelection->Fill(0.);
602 //___ get MC information __________________________________________________________________
603 Double_t ptHard = 0.;
604 Double_t nTrials = 1; // trials for MC trigger weight for real data
605
606 if(fMCEvent){
607 AliGenEventHeader* genHeader = fMCEvent->GenEventHeader();
608 if(genHeader){
609 AliGenPythiaEventHeader* pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
610 AliGenHijingEventHeader* hijingGenHeader = 0x0;
611
612 if(pythiaGenHeader){
613 if(fDebug>3) Printf("%s:%d pythiaGenHeader found", (char*)__FILE__,__LINE__);
614 nTrials = pythiaGenHeader->Trials();
615 ptHard = pythiaGenHeader->GetPtHard();
616
617 fh1PtHard->Fill(ptHard);
618 fh1PtHardTrials->Fill(ptHard,nTrials);
619
620
621 } else { // no pythia, hijing?
622
623 if(fDebug>3) Printf("%s:%d no pythiaGenHeader found", (char*)__FILE__,__LINE__);
624
625 hijingGenHeader = dynamic_cast<AliGenHijingEventHeader*>(genHeader);
626 if(!hijingGenHeader){
627 Printf("%s:%d no pythiaGenHeader or hjingGenHeader found", (char*)__FILE__,__LINE__);
628 } else {
629 if(fDebug>3) Printf("%s:%d hijingGenHeader found", (char*)__FILE__,__LINE__);
630 }
631 }
632
633 fh1Trials->Fill("#sum{ntrials}",fAvgTrials);
634 }
635 }
636
637 //___ fetch jets __________________________________________________________________________
638 Int_t nJ = GetListOfJets(fJetList);
639 Int_t nJets = 0;
640 if(nJ>=0) nJets = fJetList->GetEntries();
641 if(fDebug>2){
642 Printf("%s:%d Selected jets: %d %d",(char*)__FILE__,__LINE__,nJ,nJets);
643 if(nJ != nJets) Printf("%s:%d Mismatch Selected Jets: %d %d",(char*)__FILE__,__LINE__,nJ,nJets);
644 }
645 FillJetProperties(fJetList);
646 FillJetShape(fJetList);
647 fJetList->Clear();
648 //Post output data.
649 PostData(1, fCommonHistList);
650}
651//_________________________________________________________________________________//
652
653void AliAnalysisTaskJetProperties::Terminate(Option_t *)
654{
655 // terminated
656
657 if(fDebug > 1) printf("AliAnalysisTaskJetProperties::Terminate() \n");
658}
659//_________________________________________________________________________________//
660
661Int_t AliAnalysisTaskJetProperties::GetListOfJets(TList *list)
662{
663 //this functionality is motivated by AliAnalysisTaskFragmentationFunction
664 if(fDebug > 1) printf("AliAnalysisTaskJetProperties::GetListOfJets() \n");
665 // fill list of jets selected according to type
666 if(!list){
667 if(fDebug>1) Printf("%s:%d no input list", (char*)__FILE__,__LINE__);
668 return -1;
669 }
670
671 if(fBranchJets.Length()==0){
672 Printf("%s:%d no jet branch specified", (char*)__FILE__,__LINE__);
673 if(fDebug>1)fAOD->Print();
674 return 0;
675 }
676 TClonesArray *aodJets = 0;
677 if(fBranchJets.Length()) aodJets = dynamic_cast<TClonesArray*>(fAODJets->FindListObject(fBranchJets.Data()));
678 if(!aodJets) aodJets = dynamic_cast<TClonesArray*>(fAODJets->GetList()->FindObject(fBranchJets.Data()));
679 if(fAODExtension&&!aodJets) aodJets = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fBranchJets.Data()));
680
681 if(!aodJets){
682 if(fBranchJets.Length())Printf("%s:%d no jet array with name %s in AOD", (char*)__FILE__,__LINE__,fBranchJets.Data());
683 if(fDebug>1)fAOD->Print();
684 return 0;
685 }
686
687 Int_t nJets = 0;
688 for(Int_t ijet=0; ijet<aodJets->GetEntries(); ijet++){
689 AliAODJet *tmp = dynamic_cast<AliAODJet*>(aodJets->At(ijet));
690 if(!tmp) continue;
691 if( tmp->Pt() < fJetPtCut ) continue;
692 if( tmp->Eta() < fJetEtaMin || tmp->Eta() > fJetEtaMax)continue;
693 if(fJetRejectType==kReject1Track && tmp->GetRefTracks()->GetEntriesFast()==1)continue;//rejecting 1track jet if...
694 list->Add(tmp);
695 nJets++;
696 }//ij loop
697 list->Sort();
698 return nJets;
699}
700//_________________________________________________________________________________//
701
702Int_t AliAnalysisTaskJetProperties::GetListOfJetTracks(TList* list, const AliAODJet* jet)
703{
704 //this functionality is motivated by AliAnalysisTaskFragmentationFunction
705 if(fDebug > 1) printf("AliAnalysisTaskJetProperties::GetListOfJetTracks() \n");
706
707 // list of jet tracks from trackrefs
708 Int_t nTracks = jet->GetRefTracks()->GetEntriesFast();
709 Int_t NoOfTracks=0;
710 for (Int_t itrack=0; itrack<nTracks; itrack++) {
711 if(fTrackType==kTrackUndef){
712 if(fDebug>2)Printf("%s:%d unknown track type %d in AOD", (char*)__FILE__,__LINE__,kTrackUndef);
713 return 0;
714 }//if
715 else if(fTrackType == kTrackAOD){
716 AliAODTrack* trackAod = dynamic_cast<AliAODTrack*>(jet->GetRefTracks()->At(itrack));
717 if(!trackAod){
718 AliError("expected ref track not found ");
719 continue;
720 }
721 list->Add(trackAod);
722 NoOfTracks++;
723 }//if
724
725 else if(fTrackType==kTrackAODMC){
726 AliAODMCParticle* trackmc = dynamic_cast<AliAODMCParticle*>(jet->GetRefTracks()->At(itrack));
727 if(!trackmc){
728 AliError("expected ref trackmc not found ");
729 continue;
730 }
731 list->Add(trackmc);
732 NoOfTracks++;
733 }//if
734
735 else if (fTrackType==kTrackKine){
736 AliVParticle* trackkine = dynamic_cast<AliVParticle*>(jet->GetRefTracks()->At(itrack));
737 if(!trackkine){
738 AliError("expected ref trackkine not found ");
739 continue;
740 }
741 list->Add(trackkine);
742 NoOfTracks++;
743 }//if
744 }//itrack loop
745 list->Sort();
746 return NoOfTracks;
747}//initial loop
748//_________________________________________________________________________________//
749
750void AliAnalysisTaskJetProperties::FillJetProperties(TList *jetlist){
751 //filling up the histograms jets and tracks inside jet
752 if(fDebug > 1) printf("AliAnalysisTaskJetProperties::FillJetProperties() \n");
753
754 for(Int_t iJet=0; iJet < jetlist->GetEntries(); iJet++){
755 Float_t JetPt;Float_t JetPhi;Float_t JetEta;Float_t JetE;
756 AliAODJet *jet = dynamic_cast<AliAODJet*>(jetlist->At(iJet));
757 if(!jet)continue;
758 JetEta = jet->Eta();
759 JetPhi = jet->Phi();
760 JetPt = jet->Pt();
761 JetE = jet->E();
762 fh2EtaJet ->Fill(JetPt,JetEta);
763 fh2PhiJet ->Fill(JetPt,JetPhi);
764 fh2PtJet ->Fill(JetPt,JetPt);
765 fh1PtJet ->Fill(JetPt);
766
767 Int_t nJT = GetListOfJetTracks(fTrackList,jet);
768 Int_t nJetTracks = 0;
769 if(nJT>=0) nJetTracks = fTrackList->GetEntries();
770 if(fDebug>2){
771 Printf("%s:%d Jet tracks: %d %d",(char*)__FILE__,__LINE__,nJT,nJetTracks);
772 if(nJT != nJetTracks) Printf("%s:%d Mismatch Jet Tracks: %d %d",(char*)__FILE__,__LINE__,nJT,nJetTracks);
773 }
774
775 fh2NtracksJet ->Fill(JetPt,fTrackList->GetEntries());
776 fProNtracksJet ->Fill(JetPt,fTrackList->GetEntries());
777
778 for (Int_t j =0; j< fTrackList->GetEntries(); j++){
779 if(fTrackType==kTrackUndef)continue;
780 Float_t TrackEta=-99.0; Float_t TrackPt=-99.0; Float_t TrackPhi=-99.0;
781 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;
782 if(fTrackType==kTrackAOD){
783 AliAODTrack *trackaod = dynamic_cast<AliAODTrack*>(fTrackList->At(j));
784 if(!trackaod)continue;
785 TrackEta = trackaod->Eta();
786 TrackPhi = trackaod->Phi();
787 TrackPt = trackaod->Pt();
788 }//if kTrackAOD
789 else if(fTrackType==kTrackAODMC){
790 AliAODMCParticle* trackmc = dynamic_cast<AliAODMCParticle*>(fTrackList->At(j));
791 if(!trackmc)continue;
792 TrackEta = trackmc->Eta();
793 TrackPhi = trackmc->Phi();
794 TrackPt = trackmc->Pt();
795 }//if kTrackAODMC
796 else if(fTrackType==kTrackKine){
797 AliVParticle* trackkine = dynamic_cast<AliVParticle*>(fTrackList->At(j));
798 if(!trackkine)continue;
799 TrackEta = trackkine->Eta();
800 TrackPhi = trackkine->Phi();
801 TrackPt = trackkine->Pt();
802 }//kTrackKine
803 if(JetPt)FF = TrackPt/JetPt;
804 DelEta = TMath::Abs(JetEta - TrackEta);
805 DelPhi = TMath::Abs(JetPhi - TrackPhi);
448afeb9 806 if(DelPhi>TMath::Pi())DelPhi = TMath::Abs(DelPhi-TMath::TwoPi());
59a5935a 807 DelR = TMath::Sqrt(DelEta*DelEta + DelPhi*DelPhi);
808 AreaJ = TMath::Pi()*DelR*DelR;
809 fh2EtaTrack ->Fill(JetPt,TrackEta);
810 fh2PhiTrack ->Fill(JetPt,TrackPhi);
811 fh2PtTrack ->Fill(JetPt,TrackPt);
812 fh2FF ->Fill(JetPt,FF);
813 fh2DelEta ->Fill(JetPt,DelEta);
814 fh2DelPhi ->Fill(JetPt,DelPhi);
815 fh2DelR ->Fill(JetPt,DelR);
816 }//track loop
817 fTrackList->Clear();
818 }//iJet loop
819}//FillJetProperties
820//_________________________________________________________________________________//
821
822void AliAnalysisTaskJetProperties::FillJetShape(TList *jetlist){
823 //filling up the histograms
824 if(fDebug > 1) printf("AliAnalysisTaskJetProperties::FillJetShape() \n");
825 Float_t JetEta; Float_t JetPhi; Float_t JetPt;
826 AliAODJet *jet = dynamic_cast<AliAODJet*>(jetlist->At(0));//Leading jet
827 if(jet){
828 JetEta = jet->Eta();
829 JetPhi = jet->Phi();
830 JetPt = jet->Pt();
831 fh1PtLeadingJet->Fill(JetPt);
efabdd32 832 Float_t NchSumA[50] = {0.};
833 Float_t PtSumA[50] = {0.};
834 Float_t delRPtSum80pc = 0;
835 Float_t delRNtrkSum80pc = 0;
836 Float_t PtSumDiffShape[10] = {0.0};
837 Float_t PtSumIntShape[10] = {0.0};
838 Int_t kNbinsR = 10;
839
59a5935a 840 Int_t nJT = GetListOfJetTracks(fTrackList,jet);
841 Int_t nJetTracks = 0;
842 if(nJT>=0) nJetTracks = fTrackList->GetEntries();
843 if(fDebug>2){
844 Printf("%s:%d Jet tracks: %d %d",(char*)__FILE__,__LINE__,nJT,nJetTracks);
845 if(nJT != nJetTracks) Printf("%s:%d Mismatch Jet Tracks: %d %d",(char*)__FILE__,__LINE__,nJT,nJetTracks);
846 }
847 fh2NtracksLeadingJet->Fill(JetPt,nJetTracks);
848 fProNtracksLeadingJet->Fill(JetPt,nJetTracks);
849 Int_t *index = new Int_t [nJetTracks];//dynamic array containing index
850 Float_t *delRA = new Float_t [nJetTracks];//dynamic array containing delR
851 Float_t *delEtaA = new Float_t [nJetTracks];//dynamic array containing delEta
852 Float_t *delPhiA = new Float_t [nJetTracks];//dynamic array containing delPhi
853 Float_t *trackPtA = new Float_t [nJetTracks];//dynamic array containing pt-track
854 Float_t *trackEtaA = new Float_t [nJetTracks];//dynamic array containing eta-track
855 Float_t *trackPhiA = new Float_t [nJetTracks];//dynamic array containing phi-track
856 for(Int_t ii=0; ii<nJetTracks; ii++){
857 index[ii] = 0;
858 delRA[ii] = 0.;
859 delEtaA[ii] = 0.;
860 delPhiA[ii] = 0.;
861 trackPtA[ii] = 0.;
862 trackEtaA[ii] = 0.;
863 trackPhiA[ii] = 0.;
864 }//ii loop
865
866 for (Int_t j =0; j< nJetTracks; j++){
867 if(fTrackType==kTrackUndef)continue;
868
869 Float_t TrackEta=-99.0; Float_t TrackPt=-99.0; Float_t TrackPhi=-99.0;
870 Float_t DelEta=-99.0; Float_t DelPhi=-99.0; Float_t DelR=-99.0; Float_t AreaJ=-99.0;
871
872 if(fTrackType==kTrackAOD){
873 AliAODTrack *trackaod = dynamic_cast<AliAODTrack*>(fTrackList->At(j));
874 if(!trackaod)continue;
875 TrackEta = trackaod->Eta();
876 TrackPhi = trackaod->Phi();
877 TrackPt = trackaod->Pt();
878 }//if kTrackAOD
879 else if(fTrackType==kTrackAODMC){
880 AliAODMCParticle* trackmc = dynamic_cast<AliAODMCParticle*>(fTrackList->At(j));
881 if(!trackmc)continue;
882 TrackEta = trackmc->Eta();
883 TrackPhi = trackmc->Phi();
884 TrackPt = trackmc->Pt();
885 }//if kTrackAODMC
886 else if(fTrackType==kTrackKine){
887 AliVParticle* trackkine = dynamic_cast<AliVParticle*>(fTrackList->At(j));
888 if(!trackkine)continue;
889 TrackEta = trackkine->Eta();
890 TrackPhi = trackkine->Phi();
891 TrackPt = trackkine->Pt();
892 }//if kTrackKine
893
894 DelEta = TMath::Abs(JetEta - TrackEta);
895 DelPhi = TMath::Abs(JetPhi - TrackPhi);
efabdd32 896 if(DelPhi>TMath::Pi())DelPhi = TMath::Abs(DelPhi-TMath::TwoPi());
59a5935a 897 DelR = TMath::Sqrt(DelEta*DelEta + DelPhi*DelPhi);
898 AreaJ = TMath::Pi()*DelR*DelR;
899
900 fh2AreaCh ->Fill(JetPt,AreaJ);
901 fProAreaCh->Fill(JetPt,AreaJ);
902 delRA[j] = DelR;
903 delEtaA[j] = DelEta;
904 delPhiA[j] = DelPhi;
905 trackPtA[j] = TrackPt;
906 trackEtaA[j] = TrackEta;
907 trackPhiA[j] = TrackPhi;
908
efabdd32 909 //calculating diff and integrated jet shapes
910 Float_t kDeltaR = 0.1;
911 Float_t RMin = kDeltaR/2.0;
912 Float_t RMax = kDeltaR/2.0;
913 Float_t tmpR = 0.05;
914 for(Int_t ii1=0; ii1<kNbinsR;ii1++){
915 if((DelR > (tmpR-RMin)) && (DelR <=(tmpR+RMax)))PtSumDiffShape[ii1]+= TrackPt;
916 if(DelR>0.0 && DelR <=(tmpR+RMax))PtSumIntShape[ii1]+= TrackPt;
917 tmpR += 0.1;
918 }//ii1 loop
919
59a5935a 920 for(Int_t ibin=1; ibin<=50; ibin++){
921 Float_t xlow = 0.02*(ibin-1);
922 Float_t xup = 0.02*ibin;
923 if( xlow <= DelR && DelR < xup){
924 NchSumA[ibin-1]++;
925 PtSumA[ibin-1]+= TrackPt;
926 }//if loop
927 }//for ibin loop
928 }//track loop
929 fTrackList->Clear();
efabdd32 930
931 //---------------------
932 Float_t tmp1R = 0.05;
933 for(Int_t jj1=0; jj1<kNbinsR;jj1++){
934 if(JetPt>20 && JetPt<=100){
935 fProDiffJetShape->Fill(tmp1R,PtSumDiffShape[jj1]/JetPt);
936 fProIntJetShape ->Fill(tmp1R,PtSumIntShape[jj1]/JetPt);
937 }
938 Float_t jetPtMin0=20.0; Float_t jetPtMax0=30.0;
939 for(Int_t k=0; k<5; k++){
940 if(k==0){jetPtMin0=20.0;jetPtMax0=30.0;}
941 if(k==1){jetPtMin0=30.0;jetPtMax0=40.0;}
942 if(k==2){jetPtMin0=40.0;jetPtMax0=60.0;}
943 if(k==3){jetPtMin0=60.0;jetPtMax0=80.0;}
944 if(k==4){jetPtMin0=80.0;jetPtMax0=100.0;}
945 if(JetPt>jetPtMin0 && JetPt<=jetPtMax0){
946 fProDiffJetShapeA[k]->Fill(tmp1R,PtSumDiffShape[jj1]/JetPt);
947 fProIntJetShapeA[k] ->Fill(tmp1R,PtSumIntShape[jj1]/JetPt);
948 }//if
949 }//k loop
950 tmp1R +=0.1;
951 }//jj1 loop
952 //----------------------//
59a5935a 953 Float_t PtSum = 0;
954 Int_t NtrkSum = 0;
955 Bool_t iflagPtSum = kFALSE;
956 Bool_t iflagNtrkSum = kFALSE;
957 TMath::Sort(nJetTracks,delRA,index,0);
958 for(Int_t ii=0; ii<nJetTracks; ii++){
959 NtrkSum ++;
960 PtSum += trackPtA[index[ii]];
961 /*
962 cout << index[ii] << "\t" <<
963 delR[ii] << "\t" <<
964 delEta[ii]<< "\t" <<
965 delPhi[ii]<< "\t" <<
966 trackPt[ii]<< "\t" <<
967 trackEta[ii]<< "\t" <<
968 trackPhi[ii]<< "\t DelR " <<
969 delR[index[ii]] << endl;
970 */
971 if(!iflagNtrkSum){
972 if((Float_t)NtrkSum/(Float_t)nJetTracks > 0.79){
973 delRNtrkSum80pc = delRA[index[ii]];
974 iflagNtrkSum = kTRUE;
975 }//if NtrkSum
976 }//if iflag
977 if(!iflagPtSum){
978 if(PtSum/JetPt >= 0.8000){
979 delRPtSum80pc = delRA[index[ii]];
980 iflagPtSum = kTRUE;
981 }//if PtSum
982 }//if iflag
983 }//track loop 2nd
984 delete [] index;
985 delete [] delRA;
986 delete [] delEtaA;
987 delete [] delPhiA;
988 delete [] trackPtA;
989 delete [] trackEtaA;
990 delete [] trackPhiA;
991 fh2DelR80pcNch ->Fill(JetPt,delRNtrkSum80pc);
992 fProDelR80pcNch->Fill(JetPt,delRNtrkSum80pc);
993 fh2DelR80pcPt ->Fill(JetPt,delRPtSum80pc);
994 fProDelR80pcPt ->Fill(JetPt,delRPtSum80pc);
995 for(Int_t ibin=0; ibin<50; ibin++){
996 Float_t iR = 0.02*ibin + 0.01;
997 fh3PtDelRNchSum->Fill(JetPt,iR,NchSumA[ibin]);
998 fh3PtDelRPtSum ->Fill(JetPt,iR,PtSumA[ibin]);
999 Float_t jetPtMin=20.0; Float_t jetPtMax=30.0;
1000 for(Int_t k=0; k<5; k++){
1001 if(k==0){jetPtMin=20.0;jetPtMax=30.0;}
1002 if(k==1){jetPtMin=30.0;jetPtMax=40.0;}
1003 if(k==2){jetPtMin=40.0;jetPtMax=60.0;}
1004 if(k==3){jetPtMin=60.0;jetPtMax=80.0;}
1005 if(k==4){jetPtMin=80.0;jetPtMax=100.0;}
1006 if(JetPt>jetPtMin && JetPt<jetPtMax){
1007 fProDelRNchSum[k]->Fill(iR,NchSumA[ibin]);
1008 fProDelRPtSum[k]->Fill(iR,PtSumA[ibin]);
1009 }//if
1010 }//k loop
1011 }//ibin loop
1012 }//if(jet)
1013}//FillJetShape()
1014//_________________________________________________________________________________//
1015
1016