]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGJE/AliAnalysisTaskJetCore.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Updates in the drawing macros
[u/mrichter/AliRoot.git] / PWGJE / AliAnalysisTaskJetCore.cxx
CommitLineData
a9e585a7 1
568f8fa2 2// ******************************************
3// This task computes several jet observables like
4// the fraction of energy in inner and outer coronnas,
da93bb11 5// jet-track correlations,triggered jet shapes and
568f8fa2 6// correlation strength distribution of particles inside jets.
7// Author: lcunquei@cern.ch
8// *******************************************
9
10
11/**************************************************************************
12 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
13 * *
14 * Author: The ALICE Off-line Project. *
15 * Contributors are mentioned in the code where appropriate. *
16 * *
17 * Permission to use, copy, modify and distribute this software and its *
18 * documentation strictly for non-commercial purposes is hereby granted *
19 * without fee, provided that the above copyright notice appears in all *
20 * copies and that both the copyright notice and this permission notice *
21 * appear in the supporting documentation. The authors make no claims *
22 * about the suitability of this software for any purpose. It is *
23 * provided "as is" without express or implied warranty. *
24 **************************************************************************/
25
26
75bf77e3 27#include "TChain.h"
28#include "TTree.h"
29#include "TMath.h"
30#include "TH1F.h"
31#include "TH2F.h"
32#include "TH3F.h"
33#include "THnSparse.h"
34#include "TCanvas.h"
35
36#include "AliLog.h"
37
38#include "AliAnalysisTask.h"
39#include "AliAnalysisManager.h"
40
41#include "AliVEvent.h"
42#include "AliESDEvent.h"
43#include "AliESDInputHandler.h"
44#include "AliCentrality.h"
45#include "AliAnalysisHelperJetTasks.h"
46#include "AliInputEventHandler.h"
47#include "AliAODJetEventBackground.h"
85b5b73e 48#include "AliAODMCParticle.h"
75bf77e3 49#include "AliAnalysisTaskFastEmbedding.h"
75bf77e3 50#include "AliAODEvent.h"
ea693273 51#include "AliAODHandler.h"
75bf77e3 52#include "AliAODJet.h"
53
54#include "AliAnalysisTaskJetCore.h"
55
56ClassImp(AliAnalysisTaskJetCore)
57
58AliAnalysisTaskJetCore::AliAnalysisTaskJetCore() :
59AliAnalysisTaskSE(),
60fESD(0x0),
5bd732d4 61fAODIn(0x0),
62fAODOut(0x0),
ea693273 63fAODExtension(0x0),
75bf77e3 64fBackgroundBranch(""),
ea693273 65fNonStdFile(""),
75bf77e3 66fIsPbPb(kTRUE),
67fOfflineTrgMask(AliVEvent::kAny),
68fMinContribVtx(1),
a9e585a7 69fVtxZMin(-10.),
70fVtxZMax(10.),
75bf77e3 71fEvtClassMin(0),
72fEvtClassMax(4),
8b47ec90 73fFilterMask(0),
873306cf 74fFilterMaskBestPt(16),
ae73f4c9 75fFilterType(2),
ea693273 76fRadioFrac(0.2),
77fMinDist(0.1),
75bf77e3 78fCentMin(0.),
79fCentMax(100.),
80fNInputTracksMin(0),
81fNInputTracksMax(-1),
ea693273 82fAngStructCloseTracks(0),
da93bb11 83fCheckMethods(0),
529e2916 84fDoEventMixing(0),
447d2a5b 85fFlagPhiBkg(0),
ca1009f1 86fFlagEtaBkg(0),
16a4c8e8 87fFlagJetHadron(0),
d5e7475c 88fFlagRandom(0),
ae1e07c1 89fFlagOnlyRecoil(0),
e22029ea 90fFlagOnlyHardest(0),
85b5b73e 91fTrackTypeRec(kTrackUndef),
d5e7475c 92fRPAngle(0),
85b5b73e 93fNRPBins(100),
75bf77e3 94fJetEtaMin(-.5),
95fJetEtaMax(.5),
447d2a5b 96fNevents(0),
97fTindex(0),
98fTrigBufferIndex(0),
99fCountAgain(0),
75bf77e3 100fJetPtMin(20.),
101fJetTriggerExcludeMask(AliAODJet::kHighTrackPtTriggered),
102fJetPtFractionMin(0.5),
103fNMatchJets(4),
104fMatchMaxDist(0.8),
105fKeepJets(kFALSE),
75bf77e3 106fkNbranches(2),
107fkEvtClasses(12),
108fOutputList(0x0),
109fbEvent(kTRUE),
110fHistEvtSelection(0x0),
da93bb11 111fhnDeltaR(0x0),
529e2916 112fhnMixedEvents(0x0),
75bf77e3 113fh2JetCoreMethod1C10(0x0),
114fh2JetCoreMethod2C10(0x0),
ea693273 115fh2JetCoreMethod1C20(0x0),
116fh2JetCoreMethod2C20(0x0),
75bf77e3 117fh2JetCoreMethod1C30(0x0),
118fh2JetCoreMethod2C30(0x0),
75bf77e3 119fh2JetCoreMethod1C60(0x0),
120fh2JetCoreMethod2C60(0x0),
d5e7475c 121fh3JetTrackC3060(0x0),
e03c221a 122fh3JetTrackC20(0x0),
ea693273 123fh2AngStructpt1C10(0x0),
124fh2AngStructpt2C10(0x0),
125fh2AngStructpt3C10(0x0),
126fh2AngStructpt4C10(0x0),
127fh2AngStructpt1C20(0x0),
128fh2AngStructpt2C20(0x0),
129fh2AngStructpt3C20(0x0),
130fh2AngStructpt4C20(0x0),
131fh2AngStructpt1C30(0x0),
132fh2AngStructpt2C30(0x0),
133fh2AngStructpt3C30(0x0),
134fh2AngStructpt4C30(0x0),
135fh2AngStructpt1C60(0x0),
136fh2AngStructpt2C60(0x0),
137fh2AngStructpt3C60(0x0),
da93bb11 138fh2AngStructpt4C60(0x0),
8205e054 139fh2Ntriggers(0x0),
16a4c8e8 140fh2Ntriggers2(0x0),
ae1e07c1 141fh3JetDensity(0x0),
142fh3JetDensityA4(0x0),
d5e7475c 143fh2RPJets(0x0),
85b5b73e 144fh2RPT(0x0),
145fh3spectriggeredC20RP(0x0),
3353f803 146fh3spectriggeredC20(0x0),
d5e7475c 147fh3spectriggeredC3060(0x0)
148
da93bb11 149
75bf77e3 150{
151 // default Constructor
152
529e2916 153
154 // Trigger buffer.
155 for(Int_t i=0; i<10; i++) {
447d2a5b 156 for(Int_t j=0; j<6; j++) {
529e2916 157 fTrigBuffer[i][j]=0;
158 }
159 }
160
161
162
163
164
75bf77e3 165 fJetBranchName[0] = "";
166 fJetBranchName[1] = "";
167
168 fListJets[0] = new TList;
169 fListJets[1] = new TList;
170}
171
172AliAnalysisTaskJetCore::AliAnalysisTaskJetCore(const char *name) :
173AliAnalysisTaskSE(name),
174fESD(0x0),
5bd732d4 175fAODIn(0x0),
176fAODOut(0x0),
ea693273 177fAODExtension(0x0),
75bf77e3 178fBackgroundBranch(""),
ea693273 179fNonStdFile(""),
75bf77e3 180fIsPbPb(kTRUE),
181fOfflineTrgMask(AliVEvent::kAny),
182fMinContribVtx(1),
a9e585a7 183fVtxZMin(-10.),
184fVtxZMax(10.),
75bf77e3 185fEvtClassMin(0),
186fEvtClassMax(4),
8b47ec90 187fFilterMask(0),
873306cf 188fFilterMaskBestPt(16),
ae73f4c9 189fFilterType(2),
ea693273 190fRadioFrac(0.2),
191fMinDist(0.1),
75bf77e3 192fCentMin(0.),
193fCentMax(100.),
194fNInputTracksMin(0),
195fNInputTracksMax(-1),
ea693273 196fAngStructCloseTracks(0),
da93bb11 197fCheckMethods(0),
529e2916 198fDoEventMixing(0),
447d2a5b 199fFlagPhiBkg(0),
ca1009f1 200fFlagEtaBkg(0),
16a4c8e8 201fFlagJetHadron(0),
d5e7475c 202fFlagRandom(0),
ae1e07c1 203fFlagOnlyRecoil(0),
e22029ea 204fFlagOnlyHardest(0),
85b5b73e 205fTrackTypeRec(kTrackUndef),
d5e7475c 206fRPAngle(0),
85b5b73e 207fNRPBins(100),
75bf77e3 208fJetEtaMin(-.5),
209fJetEtaMax(.5),
447d2a5b 210fNevents(0),
211fTindex(0),
212fTrigBufferIndex(0),
213fCountAgain(0),
75bf77e3 214fJetPtMin(20.),
215fJetTriggerExcludeMask(AliAODJet::kHighTrackPtTriggered),
216fJetPtFractionMin(0.5),
217fNMatchJets(4),
218fMatchMaxDist(0.8),
219fKeepJets(kFALSE),
75bf77e3 220fkNbranches(2),
221fkEvtClasses(12),
222fOutputList(0x0),
223fbEvent(kTRUE),
224fHistEvtSelection(0x0),
da93bb11 225fhnDeltaR(0x0),
529e2916 226fhnMixedEvents(0x0),
75bf77e3 227fh2JetCoreMethod1C10(0x0),
228fh2JetCoreMethod2C10(0x0),
ea693273 229fh2JetCoreMethod1C20(0x0),
230fh2JetCoreMethod2C20(0x0),
75bf77e3 231fh2JetCoreMethod1C30(0x0),
232fh2JetCoreMethod2C30(0x0),
75bf77e3 233fh2JetCoreMethod1C60(0x0),
234fh2JetCoreMethod2C60(0x0),
d5e7475c 235fh3JetTrackC3060(0x0),
e03c221a 236fh3JetTrackC20(0x0),
ea693273 237fh2AngStructpt1C10(0x0),
238fh2AngStructpt2C10(0x0),
239fh2AngStructpt3C10(0x0),
240fh2AngStructpt4C10(0x0),
241fh2AngStructpt1C20(0x0),
242fh2AngStructpt2C20(0x0),
243fh2AngStructpt3C20(0x0),
244fh2AngStructpt4C20(0x0),
245fh2AngStructpt1C30(0x0),
246fh2AngStructpt2C30(0x0),
247fh2AngStructpt3C30(0x0),
248fh2AngStructpt4C30(0x0),
249fh2AngStructpt1C60(0x0),
250fh2AngStructpt2C60(0x0),
251fh2AngStructpt3C60(0x0),
da93bb11 252fh2AngStructpt4C60(0x0),
8205e054 253fh2Ntriggers(0x0),
16a4c8e8 254fh2Ntriggers2(0x0),
ae1e07c1 255fh3JetDensity(0x0),
256fh3JetDensityA4(0x0),
d5e7475c 257fh2RPJets(0x0),
85b5b73e 258fh2RPT(0x0),
259fh3spectriggeredC20RP(0x0),
3353f803 260fh3spectriggeredC20(0x0),
d5e7475c 261fh3spectriggeredC3060(0x0)
262
75bf77e3 263 {
264 // Constructor
265
529e2916 266
267 for(Int_t i=0; i<10; i++) {
447d2a5b 268 for(Int_t j=0; j<6; j++) {
529e2916 269 fTrigBuffer[i][j]=0;
270 }
271 }
272
273
274
75bf77e3 275 fJetBranchName[0] = "";
276 fJetBranchName[1] = "";
277
278 fListJets[0] = new TList;
279 fListJets[1] = new TList;
280
281 DefineOutput(1, TList::Class());
282}
283
284AliAnalysisTaskJetCore::~AliAnalysisTaskJetCore()
285{
286 delete fListJets[0];
287 delete fListJets[1];
288}
289
290void AliAnalysisTaskJetCore::SetBranchNames(const TString &branch1, const TString &branch2)
291{
292 fJetBranchName[0] = branch1;
293 fJetBranchName[1] = branch2;
294}
295
296void AliAnalysisTaskJetCore::Init()
297{
298
299 // check for jet branches
300 if(!strlen(fJetBranchName[0].Data()) || !strlen(fJetBranchName[1].Data())){
301 AliError("Jet branch name not set.");
302 }
303
304}
305
306void AliAnalysisTaskJetCore::UserCreateOutputObjects()
307{
308 // Create histograms
309 // Called once
310 OpenFile(1);
311 if(!fOutputList) fOutputList = new TList;
312 fOutputList->SetOwner(kTRUE);
313
314 Bool_t oldStatus = TH1::AddDirectoryStatus();
315 TH1::AddDirectory(kFALSE);
316
317
318 fHistEvtSelection = new TH1I("fHistEvtSelection", "event selection", 6, -0.5, 5.5);
319 fHistEvtSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
320 fHistEvtSelection->GetXaxis()->SetBinLabel(2,"events IN");
321 fHistEvtSelection->GetXaxis()->SetBinLabel(3,"event selection (rejected)");
322 fHistEvtSelection->GetXaxis()->SetBinLabel(4,"vertex cut (rejected)");
323 fHistEvtSelection->GetXaxis()->SetBinLabel(5,"centrality (rejected)");
324 fHistEvtSelection->GetXaxis()->SetBinLabel(6,"multiplicity (rejected)");
325
529e2916 326 UInt_t entries = 0; // bit coded, see GetDimParams() below
327 entries = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<6 |1<<7;
da93bb11 328 fhnDeltaR = NewTHnSparseF("fhnDeltaR", entries);
c42c6838 329
20dcc500 330 //change binning in pTtrack
331 Double_t *xPt3 = new Double_t[10];
c42c6838 332 xPt3[0] = 0.;
20dcc500 333 for(int i = 1; i<=9;i++){
4e90d948 334 if(xPt3[i-1]<2)xPt3[i] = xPt3[i-1] + 0.4; // 1 - 5
335 else if(xPt3[i-1]<11)xPt3[i] = xPt3[i-1] + 3; // 5 - 12
20dcc500 336 else xPt3[i] = xPt3[i-1] + 150.; // 18
c42c6838 337 }
338 fhnDeltaR->SetBinEdges(2,xPt3);
339 delete [] xPt3;
340
20dcc500 341 //change binning in HTI
342 Double_t *xPt4 = new Double_t[14];
343 xPt4[0] = 0.;
344 for(int i = 1; i<=13;i++){
345 if(xPt4[i-1]<10)xPt4[i] = xPt4[i-1] + 1; // 1 - 10
346 else if(xPt4[i-1]<20)xPt4[i] = xPt4[i-1] + 5; // 10 - 12
347 else xPt4[i] = xPt4[i-1] + 30.; // 13
348 }
349 fhnDeltaR->SetBinEdges(6,xPt4);
350 delete [] xPt4;
351
352
353
c42c6838 354
355
356
529e2916 357
c42c6838 358 if(fDoEventMixing){
529e2916 359 UInt_t cifras = 0; // bit coded, see GetDimParams() below
447d2a5b 360 cifras = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<7;
c42c6838 361 fhnMixedEvents = NewTHnSparseF("fhnMixedEvents", cifras);}
ea693273 362
da93bb11 363 if(fCheckMethods){
75bf77e3 364 fh2JetCoreMethod1C10 = new TH2F("JetCoreMethod1C10","",150, 0., 150.,100, 0., 1.5);
365 fh2JetCoreMethod2C10 = new TH2F("JetCoreMethod2C10","",150, 0., 150.,100, 0., 1.5);
ea693273 366 fh2JetCoreMethod1C20 = new TH2F("JetCoreMethod1C20","",150, 0., 150.,100, 0., 1.5);
367 fh2JetCoreMethod2C20 = new TH2F("JetCoreMethod2C20","",150, 0., 150.,100, 0., 1.5);
75bf77e3 368 fh2JetCoreMethod1C30 = new TH2F("JetCoreMethod1C30","",150, 0., 150.,100, 0., 1.5);
369 fh2JetCoreMethod2C30 = new TH2F("JetCoreMethod2C30","",150, 0., 150.,100, 0., 1.5);
75bf77e3 370 fh2JetCoreMethod1C60 = new TH2F("JetCoreMethod1C60","",150, 0., 150.,100, 0., 1.5);
da93bb11 371 fh2JetCoreMethod2C60 = new TH2F("JetCoreMethod2C60","",150, 0., 150.,100, 0., 1.5);}
372
ae73f4c9 373 fh3JetTrackC3060=new TH3F("JetTrackC3060","",50,0,50,150,0.,150.,35,0.,3.5);
374 fh3JetTrackC20=new TH3F("JetTrackC20","",50,0,50,150,0.,150.,35,0.,3.5);
da93bb11 375 if(fAngStructCloseTracks>0){
a9e585a7 376 fh2AngStructpt1C10 = new TH2F("Ang struct pt1 C10","",15,0.,1.5,150,0.,10.);
377 fh2AngStructpt2C10 = new TH2F("Ang struct pt2 C10","",15,0.,1.5,150,0.,10.);
378 fh2AngStructpt3C10 = new TH2F("Ang struct pt3 C10","",15,0.,1.5,150,0.,10.);
379 fh2AngStructpt4C10 = new TH2F("Ang struct pt4 C10","",15,0.,1.5,150,0.,10.);
380 fh2AngStructpt1C20 = new TH2F("Ang struct pt1 C20","",15,0.,1.5,150,0.,10.);
381 fh2AngStructpt2C20 = new TH2F("Ang struct pt2 C20","",15,0.,1.5,150,0.,10.);
382 fh2AngStructpt3C20 = new TH2F("Ang struct pt3 C20","",15,0.,1.5,150,0.,10.);
383 fh2AngStructpt4C20 = new TH2F("Ang struct pt4 C20","",15,0.,1.5,150,0.,10.);
384 fh2AngStructpt1C30 = new TH2F("Ang struct pt1 C30","",15,0.,1.5,150,0.,10.);
385 fh2AngStructpt2C30 = new TH2F("Ang struct pt2 C30","",15,0.,1.5,150,0.,10.);
386 fh2AngStructpt3C30 = new TH2F("Ang struct pt3 C30","",15,0.,1.5,150,0.,10.);
387 fh2AngStructpt4C30 = new TH2F("Ang struct pt4 C30","",15,0.,1.5,150,0.,10.);
388 fh2AngStructpt1C60 = new TH2F("Ang struct pt1 C60","",15,0.,1.5,150,0.,10.);
389 fh2AngStructpt2C60 = new TH2F("Ang struct pt2 C60","",15,0.,1.5,150,0.,10.);
390 fh2AngStructpt3C60 = new TH2F("Ang struct pt3 C60","",15,0.,1.5,150,0.,10.);
391 fh2AngStructpt4C60 = new TH2F("Ang struct pt4 C60","",15,0.,1.5,150,0.,10.); }
20dcc500 392
ba143e7f 393
394
8205e054 395 fh2Ntriggers=new TH2F("# of triggers","",10,0.,100.,50,0.,50.);
85b5b73e 396 fh2Ntriggers2=new TH2F("# of triggers2","",50,0.,50.,50,0.,50.);
16a4c8e8 397
ae1e07c1 398 fh3JetDensity=new TH3F("Jet density vs mutliplicity A>0.4","",100,0.,4000.,100,0.,5.,10,0.,50.);
399 fh3JetDensityA4=new TH3F("Jet density vs multiplicity A>0.4","",100,0.,4000.,100,0.,5.,10,0.,50.);
ae73f4c9 400 fh2RPJets=new TH2F("RPJet","",fNRPBins,0.,fNRPBins,150,0.,150.);
401 fh2RPT=new TH2F("RPTrigger","",fNRPBins,0.,fNRPBins,150,0.,150.);
85b5b73e 402 fh3spectriggeredC20RP = new TH3F("Triggered spectrumC20RP","",3,0.,3.,140,-80.,200.,10,0.,50.);
ae1e07c1 403 fh3spectriggeredC20 = new TH3F("Triggered spectrumC20","",100,0.,1.,140,-80.,200.,10,0.,50.);
404 fh3spectriggeredC3060 = new TH3F("Triggered spectrumC3060","",100,0.,1.,140,-80.,200.,10,0.,50.);
3353f803 405
d5e7475c 406
407
75bf77e3 408 fOutputList->Add(fHistEvtSelection);
a9e585a7 409
da93bb11 410 fOutputList->Add(fhnDeltaR);
411
529e2916 412 fOutputList->Add(fhnMixedEvents);
da93bb11 413
414
415
416 if(fCheckMethods){
529e2916 417
75bf77e3 418 fOutputList->Add(fh2JetCoreMethod1C10);
419 fOutputList->Add(fh2JetCoreMethod2C10);
ea693273 420 fOutputList->Add(fh2JetCoreMethod1C20);
421 fOutputList->Add(fh2JetCoreMethod2C20);
75bf77e3 422 fOutputList->Add(fh2JetCoreMethod1C30);
423 fOutputList->Add(fh2JetCoreMethod2C30);
75bf77e3 424 fOutputList->Add(fh2JetCoreMethod1C60);
da93bb11 425 fOutputList->Add(fh2JetCoreMethod2C60);}
a9e585a7 426
d5e7475c 427 fOutputList->Add(fh3JetTrackC3060);
e03c221a 428 fOutputList->Add(fh3JetTrackC20);
85b5b73e 429
4e90d948 430
75bf77e3 431
a9e585a7 432
433
434 if(fAngStructCloseTracks>0){
ea693273 435 fOutputList->Add(fh2AngStructpt1C10);
436 fOutputList->Add(fh2AngStructpt2C10);
437 fOutputList->Add(fh2AngStructpt3C10);
438 fOutputList->Add(fh2AngStructpt4C10);
439 fOutputList->Add(fh2AngStructpt1C20);
440 fOutputList->Add(fh2AngStructpt2C20);
441 fOutputList->Add(fh2AngStructpt3C20);
442 fOutputList->Add(fh2AngStructpt4C20);
443 fOutputList->Add(fh2AngStructpt1C30);
444 fOutputList->Add(fh2AngStructpt2C30);
445 fOutputList->Add(fh2AngStructpt3C30);
446 fOutputList->Add(fh2AngStructpt4C30);
447 fOutputList->Add(fh2AngStructpt1C60);
448 fOutputList->Add(fh2AngStructpt2C60);
449 fOutputList->Add(fh2AngStructpt3C60);
a9e585a7 450 fOutputList->Add(fh2AngStructpt4C60);}
20dcc500 451
452
d5e7475c 453
5bd732d4 454
ba143e7f 455
8205e054 456 fOutputList->Add(fh2Ntriggers);
16a4c8e8 457 fOutputList->Add(fh2Ntriggers2);
ae1e07c1 458 fOutputList->Add(fh3JetDensity);
459 fOutputList->Add(fh3JetDensityA4);
d5e7475c 460 fOutputList->Add(fh2RPJets);
85b5b73e 461 fOutputList->Add(fh2RPT);
ae73f4c9 462 fOutputList->Add(fh3spectriggeredC20RP);
3353f803 463 fOutputList->Add(fh3spectriggeredC20);
464 fOutputList->Add(fh3spectriggeredC3060);
465
d5e7475c 466
75bf77e3 467 // =========== Switch on Sumw2 for all histos ===========
468 for (Int_t i=0; i<fOutputList->GetEntries(); ++i) {
469 TH1 *h1 = dynamic_cast<TH1*>(fOutputList->At(i));
470 if (h1){
471 h1->Sumw2();
472 continue;
473 }
ea693273 474 THnSparse *hn = dynamic_cast<THnSparse*>(fOutputList->At(i));
475 if (hn){
476 hn->Sumw2();
477 }
75bf77e3 478 }
479 TH1::AddDirectory(oldStatus);
480
481 PostData(1, fOutputList);
482}
483
484void AliAnalysisTaskJetCore::UserExec(Option_t *)
485{
486
487
488 if(!strlen(fJetBranchName[0].Data()) || !strlen(fJetBranchName[1].Data())){
489 AliError("Jet branch name not set.");
490 return;
491 }
492
493 fESD=dynamic_cast<AliESDEvent*>(InputEvent());
494 if (!fESD) {
495 AliError("ESD not available");
5bd732d4 496 fAODIn = dynamic_cast<AliAODEvent*>(InputEvent());
497 }
498 fAODOut = dynamic_cast<AliAODEvent*>(AODEvent());
499
500 static AliAODEvent* aod = 0;
501 // take all other information from the aod we take the tracks from
502 if(!aod){
503 if(!fESD)aod = fAODIn;
504 else aod = fAODOut;}
505
506
ea693273 507
508 if(fNonStdFile.Length()!=0){
509 // case that we have an AOD extension we need can fetch the jets from the extended output
510 AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
511 fAODExtension = (aodH?aodH->GetExtension(fNonStdFile.Data()):0);
512 if(!fAODExtension){
513 if(fDebug>1)Printf("AODExtension found for %s",fNonStdFile.Data());
514 }}
515
516
75bf77e3 517 // -- event selection --
518 fHistEvtSelection->Fill(1); // number of events before event selection
519
520 // physics selection
521 AliInputEventHandler* inputHandler = (AliInputEventHandler*)
522 ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
e22029ea 523 cout<<inputHandler->IsEventSelected()<<" "<<fOfflineTrgMask<<endl;
75bf77e3 524 if(!(inputHandler->IsEventSelected() & fOfflineTrgMask)){
525 if(fDebug) Printf(" Trigger Selection: event REJECTED ... ");
526 fHistEvtSelection->Fill(2);
527 PostData(1, fOutputList);
528 return;
529 }
530
531 // vertex selection
5bd732d4 532 if(!aod){
52721bee 533 if(fDebug) Printf("%s:%d No AOD",(char*)__FILE__,__LINE__);
534 fHistEvtSelection->Fill(3);
535 PostData(1, fOutputList);
536 }
5bd732d4 537 AliAODVertex* primVtx = aod->GetPrimaryVertex();
52721bee 538
539 if(!primVtx){
540 if(fDebug) Printf("%s:%d No primVtx",(char*)__FILE__,__LINE__);
541 fHistEvtSelection->Fill(3);
d27719bf 542 PostData(1, fOutputList);
543 return;
52721bee 544 }
545
75bf77e3 546 Int_t nTracksPrim = primVtx->GetNContributors();
547 if ((nTracksPrim < fMinContribVtx) ||
548 (primVtx->GetZ() < fVtxZMin) ||
549 (primVtx->GetZ() > fVtxZMax) ){
550 if(fDebug) Printf("%s:%d primary vertex z = %f: event REJECTED...",(char*)__FILE__,__LINE__,primVtx->GetZ());
551 fHistEvtSelection->Fill(3);
552 PostData(1, fOutputList);
553 return;
554 }
555
556 // event class selection (from jet helper task)
557 Int_t eventClass = AliAnalysisHelperJetTasks::EventClass();
558 if(fDebug) Printf("Event class %d", eventClass);
559 if (eventClass < fEvtClassMin || eventClass > fEvtClassMax){
560 fHistEvtSelection->Fill(4);
561 PostData(1, fOutputList);
562 return;
563 }
564
565 // centrality selection
566 AliCentrality *cent = 0x0;
da93bb11 567 Double_t centValue = 0.;
8205e054 568 if(fIsPbPb){
20dcc500 569 if(fESD) {cent = fESD->GetCentrality();
570 if(cent) centValue = cent->GetCentralityPercentile("V0M");}
5bd732d4 571 else centValue=aod->GetHeader()->GetCentrality();
20dcc500 572
75bf77e3 573 if(fDebug) printf("centrality: %f\n", centValue);
3353f803 574 if (centValue < fCentMin || centValue > fCentMax){
575 fHistEvtSelection->Fill(4);
576 PostData(1, fOutputList);
577 return;
8205e054 578 }}
75bf77e3 579
580
568f8fa2 581 fHistEvtSelection->Fill(0);
582 // accepted events
75bf77e3 583 // -- end event selection --
ea693273 584
75bf77e3 585 // get background
586 AliAODJetEventBackground* externalBackground = 0;
5bd732d4 587 if(fAODOut&&!externalBackground&&fBackgroundBranch.Length()){
588 externalBackground = (AliAODJetEventBackground*)(fAODOut->FindListObject(fBackgroundBranch.Data()));
75bf77e3 589 if(!externalBackground)Printf("%s:%d Background branch not found %s",(char*)__FILE__,__LINE__,fBackgroundBranch.Data());;
590 }
ea693273 591 if(fAODExtension&&!externalBackground&&fBackgroundBranch.Length()){
592 externalBackground = (AliAODJetEventBackground*)(fAODExtension->GetAOD()->FindListObject(fBackgroundBranch.Data()));
593 if(!externalBackground)Printf("%s:%d Background branch not found %s",(char*)__FILE__,__LINE__,fBackgroundBranch.Data());;
594 }
5bd732d4 595
596 if(fAODIn&&!externalBackground&&fBackgroundBranch.Length()){
597 externalBackground = (AliAODJetEventBackground*)(fAODIn->FindListObject(fBackgroundBranch.Data()));
598 if(!externalBackground)Printf("%s:%d Background branch not found %s",(char*)__FILE__,__LINE__,fBackgroundBranch.Data());;
599 }
ea693273 600
75bf77e3 601 Float_t rho = 0;
75bf77e3 602
d5e7475c 603 if(fFlagRandom==0){
604 if(externalBackground)rho = externalBackground->GetBackground(0);}
605 if(fFlagRandom==1){
606 if(externalBackground)rho = externalBackground->GetBackground(2);}
75bf77e3 607
608 // fetch jets
609 TClonesArray *aodJets[2];
ea693273 610 aodJets[0]=0;
5bd732d4 611 if(fAODOut&&!aodJets[0]){
612 aodJets[0] = dynamic_cast<TClonesArray*>(fAODOut->FindListObject(fJetBranchName[0].Data()));
613 aodJets[1] = dynamic_cast<TClonesArray*>(fAODOut->FindListObject(fJetBranchName[1].Data())); }
ea693273 614 if(fAODExtension && !aodJets[0]){
529e2916 615 aodJets[0] = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fJetBranchName[0].Data()));
616 aodJets[1] = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fJetBranchName[1].Data())); }
5bd732d4 617 if(fAODIn&&!aodJets[0]){
618 aodJets[0] = dynamic_cast<TClonesArray*>(fAODIn->FindListObject(fJetBranchName[0].Data()));
619 aodJets[1] = dynamic_cast<TClonesArray*>(fAODIn->FindListObject(fJetBranchName[1].Data())); }
620
ea693273 621
529e2916 622 //Double_t ptsub[aodJets[0]->GetEntriesFast()];
623 //Int_t inord[aodJets[0]->GetEntriesFast()];
624 //for(Int_t n=0;n<aodJets[0]->GetEntriesFast();n++){
625 // ptsub[n]=0;
626 // inord[n]=0;}
a9e585a7 627
ea693273 628 TList ParticleList;
629 Int_t nT = GetListOfTracks(&ParticleList);
630 for (Int_t iJetType = 0; iJetType < 2; iJetType++) {
75bf77e3 631 fListJets[iJetType]->Clear();
632 if (!aodJets[iJetType]) continue;
75bf77e3 633 if(fDebug) Printf("%s: %d jets",fJetBranchName[iJetType].Data(),aodJets[iJetType]->GetEntriesFast());
75bf77e3 634 for (Int_t iJet = 0; iJet < aodJets[iJetType]->GetEntriesFast(); iJet++) {
635 AliAODJet *jet = dynamic_cast<AliAODJet*>((*aodJets[iJetType])[iJet]);
529e2916 636 if (jet) fListJets[iJetType]->Add(jet);
637 // if(iJetType==0){
638 // ptsub[iJet]=jet->Pt()-rho*jet->EffectiveAreaCharged();}
639 }}
75bf77e3 640
641 Double_t etabig=0;
642 Double_t ptbig=0;
643 Double_t areabig=0;
644 Double_t phibig=0.;
645 Double_t etasmall=0;
646 Double_t ptsmall=0;
647 Double_t areasmall=0;
75bf77e3 648 Double_t phismall=0.;
a9e585a7 649
da93bb11 650
d5e7475c 651
529e2916 652 Int_t iCount=0;
653 Int_t trigJet=-1;
654 Int_t trigBBTrack=-1;
655 Int_t trigInTrack=-1;
d5e7475c 656 fRPAngle = aod->GetHeader()->GetEventplane();
85b5b73e 657
e22029ea 658 // AliVParticle *partback = (AliVParticle*)ParticleList.At(nT);
659 //if(!partback){
660 //PostData(1, fOutputList);
661 //return;}
662
663
664 for(Int_t tt=0;tt<ParticleList.GetEntries();tt++){
665
666 if(fFlagOnlyHardest!=0){if(tt!=nT) continue;}
667 AliVParticle *partback = (AliVParticle*)ParticleList.At(tt);
8205e054 668 fh2Ntriggers->Fill(centValue,partback->Pt());
ae73f4c9 669 Int_t phiBinT = GetPhiBin(RelativePhi(partback->Phi(),fRPAngle));
85b5b73e 670 if(centValue<20.) fh2RPT->Fill(phiBinT,partback->Pt());
16a4c8e8 671
8205e054 672 Double_t accep=2.*TMath::Pi()*1.8;
4e90d948 673 Int_t injet4=0;
674 Int_t injet=0;
75bf77e3 675 for(Int_t i=0; i<fListJets[0]->GetEntries(); ++i){
676 AliAODJet* jetbig = (AliAODJet*)(fListJets[0]->At(i));
677 etabig = jetbig->Eta();
678 phibig = jetbig->Phi();
679 ptbig = jetbig->Pt();
680 if(ptbig==0) continue;
ae73f4c9 681 Int_t phiBin = GetPhiBin(RelativePhi(phibig,fRPAngle));
75bf77e3 682 areabig = jetbig->EffectiveAreaCharged();
ea693273 683 Double_t ptcorr=ptbig-rho*areabig;
75bf77e3 684 if((etabig<fJetEtaMin)||(etabig>fJetEtaMax)) continue;
ae1e07c1 685 if(areabig>=0.07) injet=injet+1;
8205e054 686 if(areabig>=0.4) injet4=injet4+1;
ca1009f1 687 Double_t dphi=RelativePhi(partback->Phi(),phibig);
688
ae1e07c1 689 if(fFlagEtaBkg==1){
ca1009f1 690 Double_t etadif= partback->Eta()-etabig;
691 if(TMath::Abs(etadif)<=0.5){
85b5b73e 692
ca1009f1 693 if(centValue<20.) fh3JetTrackC20->Fill(partback->Pt(),ptcorr,TMath::Abs(dphi));
694 if(centValue>30. && centValue<60.) fh3JetTrackC3060->Fill(partback->Pt(),ptcorr,TMath::Abs(dphi));}}
ca1009f1 695 if(fFlagEtaBkg==0){
e03c221a 696 if(centValue<20.) fh3JetTrackC20->Fill(partback->Pt(),ptcorr,TMath::Abs(dphi));
ca1009f1 697 if(centValue>30. && centValue<60.) fh3JetTrackC3060->Fill(partback->Pt(),ptcorr,TMath::Abs(dphi));}
698
699
16a4c8e8 700 if(fFlagJetHadron==0){
447d2a5b 701 if(fFlagPhiBkg!=0) if((TMath::Abs(dphi)<TMath::Pi()/2.-0.1)||(TMath::Abs(dphi)>TMath::Pi()/2.+0.1)) continue;
16a4c8e8 702 if(fFlagPhiBkg==0) if(TMath::Abs(dphi)<TMath::Pi()-0.6) continue;}
703
704 if(fFlagJetHadron!=0) if(TMath::Abs(dphi)>0.4) continue;
705
706
d5e7475c 707 if(centValue<20.) fh2RPJets->Fill(phiBin, ptcorr);
431db4a1 708
75bf77e3 709 Double_t dismin=100.;
710 Double_t ptmax=-10.;
711 Int_t index1=-1;
712 Int_t index2=-1;
e22029ea 713
85b5b73e 714 if(centValue<20.) fh3spectriggeredC20RP->Fill(phiBin,ptcorr,partback->Pt());
3353f803 715 if(centValue<20.) fh3spectriggeredC20->Fill(jetbig->EffectiveAreaCharged(),ptcorr,partback->Pt());
716 if(centValue>30. && centValue<60.) fh3spectriggeredC3060->Fill(jetbig->EffectiveAreaCharged(),ptcorr,partback->Pt());
717
718 if(ptcorr<=0) continue;
16a4c8e8 719
720
721 AliAODTrack* leadtrack=0;
da93bb11 722 Int_t ippt=0;
16a4c8e8 723 Double_t ppt=-10;
724 if(fFlagJetHadron==0){
529e2916 725 TRefArray *genTrackList = jetbig->GetRefTracks();
726 Int_t nTracksGenJet = genTrackList->GetEntriesFast();
727 AliAODTrack* genTrack;
da93bb11 728 for(Int_t ir=0; ir<nTracksGenJet; ++ir){
729 genTrack = (AliAODTrack*)(genTrackList->At(ir));
730 if(genTrack->Pt()>ppt){ppt=genTrack->Pt();
16a4c8e8 731 ippt=ir;}}
732 leadtrack=(AliAODTrack*)(genTrackList->At(ippt));
733 if(!leadtrack) continue;}
734
735 AliVParticle* leadtrackb=0;
736 if(fFlagJetHadron!=0){
737 Int_t nTb = GetHardestTrackBackToJet(jetbig);
738 leadtrackb = (AliVParticle*)ParticleList.At(nTb);
739 if(!leadtrackb) continue;
740 }
741
742
743
744
d5e7475c 745
529e2916 746 //store one trigger info
447d2a5b 747 if(iCount==0){
529e2916 748 trigJet=i;
8205e054 749 trigBBTrack=nT;
529e2916 750 trigInTrack=ippt;
751 iCount=iCount+1;}
752
d5e7475c 753
529e2916 754 if(fCheckMethods){
755 for(Int_t j=0; j<fListJets[1]->GetEntries(); ++j){
75bf77e3 756 AliAODJet* jetsmall = (AliAODJet*)(fListJets[1]->At(j));
757 etasmall = jetsmall->Eta();
758 phismall = jetsmall->Phi();
759 ptsmall = jetsmall->Pt();
760 areasmall = jetsmall->EffectiveAreaCharged();
ea693273 761 Double_t tmpDeltaR=(phismall-phibig)*(phismall-phibig)+(etasmall-etabig)*(etasmall-etabig);
762 tmpDeltaR=TMath::Sqrt(tmpDeltaR);
763 //Fraction in the jet core
764 if((ptsmall>ptmax)&&(tmpDeltaR<=fRadioFrac)){ptmax=ptsmall;
75bf77e3 765 index2=j;}
ea693273 766 if(tmpDeltaR<=dismin){ dismin=tmpDeltaR;
767 index1=j;}} //en of loop over R=0.2 jets
75bf77e3 768 //method1:most concentric jet=core
769 if(dismin<fMinDist){ AliAODJet* jetmethod1 = (AliAODJet*)(fListJets[1]->At(index1));
ea693273 770 if(centValue<10) fh2JetCoreMethod1C10->Fill(ptcorr,jetmethod1->Pt()/ptbig);
771 if((centValue>20)&&(centValue<40)) fh2JetCoreMethod1C20->Fill(ptcorr,jetmethod1->Pt()/ptbig);
772 if((centValue>30)&&(centValue<60)) fh2JetCoreMethod1C30->Fill(ptcorr,jetmethod1->Pt()/ptbig);
773 if(centValue>60) fh2JetCoreMethod1C60->Fill(ptcorr,jetmethod1->Pt()/ptbig); }
75bf77e3 774 //method2:hardest contained jet=core
775 if(index2!=-1){
776 AliAODJet* jetmethod2 = (AliAODJet*)(fListJets[1]->At(index2));
ea693273 777 if(centValue<10) fh2JetCoreMethod2C10->Fill(ptcorr,jetmethod2->Pt()/ptbig);
778 if((centValue>20)&&(centValue<40)) fh2JetCoreMethod2C20->Fill(ptcorr,jetmethod2->Pt()/ptbig);
779 if((centValue>30)&&(centValue<60)) fh2JetCoreMethod2C30->Fill(ptcorr,jetmethod2->Pt()/ptbig);
da93bb11 780 if(centValue>60) fh2JetCoreMethod2C60->Fill(ptcorr,jetmethod2->Pt()/ptbig); }}
d5e7475c 781
5bd732d4 782
85b5b73e 783
e22029ea 784 // for(int it = 0;it<ParticleList.GetEntries();++it){
785 // AliVParticle *part = (AliVParticle*)ParticleList.At(it);
786 // Double_t deltaR = jetbig->DeltaR(part);
787 // Double_t deltaEta = etabig-part->Eta();
788 // if(centValue<20.) fh2Ntriggers2->Fill(partback->Pt(),part->Pt());
789 // if(fDoEventMixing==0 && fFlagOnlyRecoil==0){
790 // Double_t deltaPhi=phibig-part->Phi();
791 // if(deltaPhi<-0.5*TMath::Pi()) deltaPhi+=2.*TMath::Pi();
792 // if(deltaPhi>3./2.*TMath::Pi()) deltaPhi-=2.*TMath::Pi();
793 // Double_t pTcont=0;
794 // if(fFlagJetHadron==0) pTcont=leadtrack->Pt();
795 // if(fFlagJetHadron!=0) pTcont=leadtrackb->Pt();
796 // Double_t jetEntries[8] = {centValue,ptcorr,part->Pt(),deltaR,deltaEta,deltaPhi,pTcont,partback->Pt()};
797 // fhnDeltaR->Fill(jetEntries);}
798
799
800 // }
85b5b73e 801
d5e7475c 802
d5e7475c 803 //end of track loop, we only do it if EM is switched off
804
ba143e7f 805
806
807
808
809
810
811
20dcc500 812
813 }
ae1e07c1 814 if(injet>0) fh3JetDensity->Fill(ParticleList.GetEntries(),injet/accep,partback->Pt());
815 if(injet4>0)fh3JetDensityA4->Fill(ParticleList.GetEntries(),injet4/accep,partback->Pt());
529e2916 816 //end of jet loop
817
e22029ea 818 }
529e2916 819
820
447d2a5b 821 if(fDoEventMixing>0){
529e2916 822 //check before if the trigger exists
823 // fTrigBuffer[i][0] = zvtx
824 // fTrigBuffer[i][1] = phi
825 // fTrigBuffer[i][2] = eta
826 // fTrigBuffer[i][3] = pt_jet
827 // fTrigBuffer[i][4] = pt_trig
447d2a5b 828 // fTrigBuffer[i][5]= centrality
829 if(fTindex==10) fTindex=0;
529e2916 830 if(fTrigBuffer[fTindex][3]>0){
831 if (TMath::Abs(fTrigBuffer[fTindex][0]-primVtx->GetZ()<2.)){
447d2a5b 832 if (TMath::Abs(fTrigBuffer[fTindex][5]-centValue<5)){
529e2916 833
834 for(int it = 0;it<nT;++it){
20dcc500 835 AliVParticle *part = (AliVParticle*)ParticleList.At(it);
836 Double_t DPhi = fTrigBuffer[fTindex][1] - part->Phi();
529e2916 837 Double_t DEta = fTrigBuffer[fTindex][2] - part->Eta();
838 Double_t DR=TMath::Sqrt(DPhi*DPhi+DEta*DEta);
839 if(DPhi<-0.5*TMath::Pi()) DPhi+=2.*TMath::Pi();
840 if(DPhi>3./2.*TMath::Pi()) DPhi-=2.*TMath::Pi();
447d2a5b 841 Double_t triggerEntries[7] = {centValue,fTrigBuffer[fTindex][3],part->Pt(),DR,DEta,DPhi,fTrigBuffer[fTindex][4]};
529e2916 842 fhnMixedEvents->Fill(triggerEntries);
843 }
844 fNevents=fNevents+1;
447d2a5b 845 if(fNevents==10) fTindex=fTindex+1;
529e2916 846 }}}
ae1e07c1 847
447d2a5b 848 if(fTindex==10&&fNevents==10) fCountAgain=0;
529e2916 849
850 // Copy the triggers from the current event into the buffer.
851 //again, only if the trigger exists:
447d2a5b 852 if(fCountAgain==0){
529e2916 853 if(trigJet>-1){
447d2a5b 854 AliAODJet* jetT = (AliAODJet*)(fListJets[0]->At(trigJet)); AliVParticle *partT = (AliVParticle*)ParticleList.At(trigBBTrack);
529e2916 855 fTrigBuffer[fTrigBufferIndex][0] = primVtx->GetZ();
856 fTrigBuffer[fTrigBufferIndex][1] = jetT->Phi();
857 fTrigBuffer[fTrigBufferIndex][2] = jetT->Eta();
858 fTrigBuffer[fTrigBufferIndex][3] = jetT->Pt()-rho*jetT->EffectiveAreaCharged();
859 fTrigBuffer[fTrigBufferIndex][4] = partT->Pt();
447d2a5b 860 fTrigBuffer[fTrigBufferIndex][5] = centValue;
529e2916 861 fTrigBufferIndex++;
447d2a5b 862 if(fTrigBufferIndex==9) {fTrigBufferIndex=0;
863 fCountAgain=1;}
529e2916 864 }
447d2a5b 865 }
529e2916 866
447d2a5b 867 }
529e2916 868
869
da93bb11 870
ea693273 871 //////////////////ANGULAR STRUCTURE//////////////////////////////////////
872
873 //tracks up to R=0.8 distant from the jet axis
529e2916 874 // if(fAngStructCloseTracks==1){
875 // TList CloseTrackList;
876 // Int_t nn=GetListOfTracksCloseToJet(&CloseTrackList,jetbig);
877 // Double_t difR=0.04;
878 // for(Int_t l=0;l<15;l++){
879 // Double_t rr=l*0.1+0.1;
880 // for(int it = 0;it<nn;++it){
881 // AliVParticle *part1 = (AliVParticle*)CloseTrackList.At(it);
882 // for(int itu=it+1;itu<CloseTrackList.GetEntries();itu++){
883 // AliVParticle *part2 = (AliVParticle*)CloseTrackList.At(itu);
884 // Double_t ptm=part1->Pt();
885 // Double_t ptn=part2->Pt();
886 // Double_t Rnm = (part1->Eta()-part2->Eta())*(part1->Eta()-part2->Eta())+(part1->Phi()-part2->Phi())*(part1->Phi()-part2->Phi());
887 // Rnm=TMath::Sqrt(Rnm);
888 // Double_t deltag=(1./(TMath::Sqrt(2*TMath::Pi())*difR))*TMath::Exp(-1.*(rr-Rnm)*(rr-Rnm)/(2.*difR*difR));
889 // Double_t stepf=0.5*(1.+TMath::Erf((rr-Rnm)/(TMath::Sqrt(2.)*difR)));
890 // if((ptcorr<85.) && (ptcorr>=70.)){up1[l]=up1[l]+ptm*ptn*Rnm*Rnm*deltag;
891 // down1[l]=down1[l]+ptm*ptn*Rnm*Rnm*stepf;}
892 // if((ptcorr<100.) && (ptcorr>=85.)){up2[l]=up2[l]+ptm*ptn*Rnm*Rnm*deltag;
893 // down2[l]=down2[l]+ptm*ptn*Rnm*Rnm*stepf;}
894 // if((ptcorr<120.) && (ptcorr>=100.)){up3[l]=up3[l]+ptm*ptn*Rnm*Rnm*deltag;
895 // down3[l]=down3[l]+ptm*ptn*Rnm*Rnm*stepf;}
896 // if((ptcorr<140.) && (ptcorr>=120.)){up4[l]=up4[l]+ptm*ptn*Rnm*Rnm*deltag;
897 // down4[l]=down4[l]+ptm*ptn*Rnm*Rnm*stepf;}}}}
898 // }
ea693273 899
529e2916 900 // //only jet constituents
901 // if(fAngStructCloseTracks==2){
ea693273 902
529e2916 903 // Double_t difR=0.04;
904 // for(Int_t l=0;l<15;l++){
905 // Double_t rr=l*0.1+0.1;
ea693273 906
907
529e2916 908 // AliAODTrack* part1;
909 // AliAODTrack* part2;
da93bb11 910
529e2916 911 // TRefArray *genTrackListb = jetbig->GetRefTracks();
912 // Int_t nTracksGenJetb = genTrackListb->GetEntriesFast();
da93bb11 913
914
915
529e2916 916 // for(Int_t it=0; it<nTracksGenJetb; ++it){
917 // part1 = (AliAODTrack*)(genTrackListb->At(it));
918 // for(Int_t itu=0; itu<nTracksGenJetb; ++itu){
919 // part2 = (AliAODTrack*)(genTrackListb->At(itu));
920 // Double_t ptm=part1->Pt();
921 // Double_t ptn=part2->Pt();
922 // Double_t Rnm = (part1->Eta()-part2->Eta())*(part1->Eta()-part2->Eta())+(part1->Phi()-part2->Phi())*(part1->Phi()-part2->Phi());
923 // Rnm=TMath::Sqrt(Rnm);
924 // Double_t deltag=(1./(TMath::Sqrt(2*TMath::Pi())*difR))*TMath::Exp(-1.*(rr-Rnm)*(rr-Rnm)/(2.*difR*difR));
925 // Double_t stepf=0.5*(1.+TMath::Erf((rr-Rnm)/(TMath::Sqrt(2.)*difR)));
926 // if((ptcorr<85.) && (ptcorr>=70.)){up1[l]=up1[l]+ptm*ptn*Rnm*Rnm*deltag;
927 // down1[l]=down1[l]+ptm*ptn*Rnm*Rnm*stepf;}
928 // if((ptcorr<100.) && (ptcorr>=85.)){up2[l]=up2[l]+ptm*ptn*Rnm*Rnm*deltag;
929 // down2[l]=down2[l]+ptm*ptn*Rnm*Rnm*stepf;}
930 // if((ptcorr<120.) && (ptcorr>=100.)){up3[l]=up3[l]+ptm*ptn*Rnm*Rnm*deltag;
931 // down3[l]=down3[l]+ptm*ptn*Rnm*Rnm*stepf;}
932 // if((ptcorr<140.) && (ptcorr>=120.)){up4[l]=up4[l]+ptm*ptn*Rnm*Rnm*deltag;
933 // down4[l]=down4[l]+ptm*ptn*Rnm*Rnm*stepf;}}}}}
934 // }
935 // //end loop over R=0.4 jets
936 // if(fAngStructCloseTracks>0){
937 // for(Int_t l=0;l<15;l++){
938 // Double_t rr=l*0.1+0.1;
939 // if(down1[l]!=0){
940 // if(centValue<10.)fh2AngStructpt1C10->Fill(rr,rr*up1[l]/down1[l]);
941 // if(centValue>20. && centValue<40.) fh2AngStructpt1C20->Fill(rr,rr*up1[l]/down1[l]);
942 // if(centValue>30. && centValue<60.) fh2AngStructpt1C30->Fill(rr,rr*up1[l]/down1[l]);
943 // if(centValue>60.) fh2AngStructpt1C60->Fill(rr,rr*up1[l]/down1[l]);}
944 // if(down2[l]!=0){
945 // if(centValue<10.) fh2AngStructpt2C10->Fill(rr,rr*up2[l]/down2[l]);
946 // if(centValue>20. && centValue<40.) fh2AngStructpt2C20->Fill(rr,rr*up2[l]/down2[l]);
947 // if(centValue>30. && centValue<60.) fh2AngStructpt2C30->Fill(rr,rr*up2[l]/down2[l]);
948 // if(centValue>60.) fh2AngStructpt2C60->Fill(rr,rr*up2[l]/down2[l]);}
949 // if(down3[l]!=0){
950 // if(centValue<10.) fh2AngStructpt3C10->Fill(rr,rr*up3[l]/down3[l]);
951 // if(centValue>20. && centValue<40.) fh2AngStructpt3C20->Fill(rr,rr*up3[l]/down3[l]);
952 // if(centValue>30. && centValue<60.) fh2AngStructpt3C30->Fill(rr,rr*up3[l]/down3[l]);
953 // if(centValue>60.) fh2AngStructpt3C60->Fill(rr,rr*up3[l]/down3[l]);}
954 // if(down4[l]!=0){
955 // if(centValue<10.) fh2AngStructpt4C10->Fill(rr,rr*up4[l]/down4[l]);
956 // if(centValue>20. && centValue<40.) fh2AngStructpt4C20->Fill(rr,rr*up4[l]/down4[l]);
957 // if(centValue>30. && centValue<60.) fh2AngStructpt4C30->Fill(rr,rr*up4[l]/down4[l]);
958 // if(centValue>60.) fh2AngStructpt4C60->Fill(rr,rr*up4[l]/down4[l]);}}}
ea693273 959
960
961
962
963
75bf77e3 964
965
966 PostData(1, fOutputList);
da93bb11 967}
75bf77e3 968
969void AliAnalysisTaskJetCore::Terminate(const Option_t *)
970{
971 // Draw result to the screen
972 // Called once at the end of the query
973
974 if (!GetOutputData(1))
975 return;
976}
977
ea693273 978
979
5bd732d4 980
ea693273 981
982
983Int_t AliAnalysisTaskJetCore::GetListOfTracks(TList *list){
984
85b5b73e 985 Int_t iCount = 0;
986 AliAODEvent *aod = 0;
987
ba143e7f 988 if(!fESD)aod = fAODIn;
5bd732d4 989 else aod = fAODOut;
85b5b73e 990 Int_t index=-1;
8205e054 991 Double_t ptmax=-10;
85b5b73e 992
993
994
995 for(int it = 0;it < aod->GetNumberOfTracks();++it){
5bd732d4 996 AliAODTrack *tr = aod->GetTrack(it);
ae73f4c9 997 Bool_t bGood = false;
998 if(fFilterType == 0)bGood = true;
999 else if(fFilterType == 1)bGood = tr->IsHybridTPCConstrainedGlobal();
1000 else if(fFilterType == 2)bGood = tr->IsHybridGlobalConstrainedGlobal();
ea693273 1001 if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask)))continue;
ae73f4c9 1002 if(bGood==false) continue;
ea693273 1003 if(TMath::Abs(tr->Eta())>0.9)continue;
1004 if(tr->Pt()<0.15)continue;
1005 list->Add(tr);
ea693273 1006 iCount++;
873306cf 1007 if(fFilterMaskBestPt>0){// only set the trigger track index for good quality tracks
1008 if(tr->TestFilterBit(fFilterMaskBestPt)){
1009 if(tr->Pt()>ptmax){
1010 ptmax=tr->Pt();
1011 index=iCount-1;
1012 }
1013 }
1014 }
1015 else{
1016 if(tr->Pt()>ptmax){
1017 ptmax=tr->Pt();
1018 index=iCount-1;
1019 }
1020 }
85b5b73e 1021 }
ea693273 1022
529e2916 1023
85b5b73e 1024 // else if (type == kTrackAODMCCharged) {
1025 // TClonesArray *tca = dynamic_cast<TClonesArray*>(aod->FindListObject(AliAODMCParticle::StdBranchName()));
1026 // if(!tca)return iCount;
1027 // for(int it = 0;it < tca->GetEntriesFast();++it){
1028 // AliAODMCParticle *part = dynamic_cast<AliAODMCParticle*>(tca->At(it));
1029 // if(!part)continue;
1030 // if(part->Pt()<0.15)continue;
1031 // if(!part->IsPhysicalPrimary())continue;
1032 // if(part->Charge()==0)continue;
1033 // if(TMath::Abs(part->Eta())>0.9)continue;
1034 // list->Add(part);
1035 // iCount++;
1036 // if(part->Pt()>ptmax){ ptmax=part->Pt();
1037 // index=iCount-1;}}}
1038 return index;
ea693273 1039
1040}
1041
da93bb11 1042 Int_t AliAnalysisTaskJetCore::GetHardestTrackBackToJet(AliAODJet *jetbig){
5bd732d4 1043
1044 AliAODEvent *aod = 0;
1045 if(!fESD)aod = fAODIn;
1046 else aod = fAODOut;
da93bb11 1047 Int_t index=-1;
1048 Double_t ptmax=-10;
1049 Double_t dphi=0;
1050 Double_t dif=0;
1051 Int_t iCount=0;
5bd732d4 1052 for(int it = 0;it < aod->GetNumberOfTracks();++it){
1053 AliAODTrack *tr = aod->GetTrack(it);
da93bb11 1054 if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask)))continue;
1055 if(TMath::Abs(tr->Eta())>0.9)continue;
1056 if(tr->Pt()<0.15)continue;
1057 iCount=iCount+1;
529e2916 1058 dphi=RelativePhi(tr->Phi(),jetbig->Phi());
16a4c8e8 1059 if(TMath::Abs(dphi)<TMath::Pi()-0.6) continue;
da93bb11 1060 if(tr->Pt()>ptmax){ ptmax=tr->Pt();
b6844e9e 1061 index=iCount-1;
1062 dif=dphi; }}
da93bb11 1063
1064 return index;
1065
1066 }
1067
1068
1069
1070
1071
1072
1073
1074
1075
ea693273 1076 Int_t AliAnalysisTaskJetCore::GetListOfTracksCloseToJet(TList *list,AliAODJet *jetbig){
1077
1078 Int_t iCount = 0;
5bd732d4 1079 AliAODEvent *aod = 0;
1080 if(!fESD)aod = fAODIn;
1081 else aod = fAODOut;
8b47ec90 1082
16a4c8e8 1083 for(int it = 0;it < aod->GetNumberOfTracks();++it){
5bd732d4 1084 AliAODTrack *tr = aod->GetTrack(it);
ea693273 1085 if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask)))continue;
1086 if(TMath::Abs(tr->Eta())>0.9)continue;
1087 if(tr->Pt()<0.15)continue;
1088 Double_t disR=jetbig->DeltaR(tr);
1089 if(disR>0.8) continue;
1090 list->Add(tr);
1091 //cout<<fAOD->GetNumberOfTracks()<<" "<<tr->Pt()<<endl;
1092 iCount++;
1093 }
1094
1095 list->Sort();
1096 return iCount;
1097
1098}
1099
1100
1101
1102
1103
1104
1105
1106
1107
1108
1109
75bf77e3 1110Int_t AliAnalysisTaskJetCore::GetNInputTracks()
1111{
1112
1113 Int_t nInputTracks = 0;
5bd732d4 1114 AliAODEvent *aod = 0;
1115 if(!fESD)aod = fAODIn;
1116 else aod = fAODOut;
75bf77e3 1117 TString jbname(fJetBranchName[1]);
1118 //needs complete event, use jets without background subtraction
1119 for(Int_t i=1; i<=3; ++i){
1120 if(jbname.Contains(Form("B%d",i))) jbname.ReplaceAll(Form("B%d",i),"B0");
1121 }
1122 // use only HI event
1123 if(jbname.Contains("AODextraonly")) jbname.ReplaceAll("AODextraonly","AOD");
1124 if(jbname.Contains("AODextra")) jbname.ReplaceAll("AODextra","AOD");
1125
1126 if(fDebug) Printf("Multiplicity from jet branch %s", jbname.Data());
5bd732d4 1127 TClonesArray *tmpAODjets = dynamic_cast<TClonesArray*>(aod->FindListObject(jbname.Data()));
75bf77e3 1128 if(!tmpAODjets){
1129 Printf("Jet branch %s not found", jbname.Data());
1130 Printf("AliAnalysisTaskJetCore::GetNInputTracks FAILED");
1131 return -1;
1132 }
1133
1134 for (Int_t iJet=0; iJet<tmpAODjets->GetEntriesFast(); iJet++){
1135 AliAODJet *jet = dynamic_cast<AliAODJet*>((*tmpAODjets)[iJet]);
1136 if(!jet) continue;
1137 TRefArray *trackList = jet->GetRefTracks();
1138 Int_t nTracks = trackList->GetEntriesFast();
1139 nInputTracks += nTracks;
1140 if(fDebug) Printf("#jet%d: %d tracks", iJet, nTracks);
1141 }
1142 if(fDebug) Printf("---> input tracks: %d", nInputTracks);
1143
1144 return nInputTracks;
1145}
1146
1147
1148
ea693273 1149Double_t AliAnalysisTaskJetCore::RelativePhi(Double_t mphi,Double_t vphi){
1150
1151 if (vphi < -1*TMath::Pi()) vphi += (2*TMath::Pi());
1152 else if (vphi > TMath::Pi()) vphi -= (2*TMath::Pi());
1153 if (mphi < -1*TMath::Pi()) mphi += (2*TMath::Pi());
1154 else if (mphi > TMath::Pi()) mphi -= (2*TMath::Pi());
1155 double dphi = mphi-vphi;
1156 if (dphi < -1*TMath::Pi()) dphi += (2*TMath::Pi());
1157 else if (dphi > TMath::Pi()) dphi -= (2*TMath::Pi());
1158 return dphi;//dphi in [-Pi, Pi]
1159}
1160
d5e7475c 1161Int_t AliAnalysisTaskJetCore::GetPhiBin(Double_t phi)
1162{
1163 Int_t phibin=-1;
1164 if(!(TMath::Abs(phi)<=2*TMath::Pi())){AliError("phi w.r.t. RP out of defined range");return -1;}
1165 Double_t phiwrtrp=TMath::ACos(TMath::Abs(TMath::Cos(phi)));
1166 phibin=Int_t(fNRPBins*phiwrtrp/(0.5*TMath::Pi()));
1167 if(phibin<0||phibin>=fNRPBins){AliError("Phi Bin not defined");}
1168 return phibin;
1169}
1170
1171
75bf77e3 1172
1173
da93bb11 1174THnSparse* AliAnalysisTaskJetCore::NewTHnSparseF(const char* name, UInt_t entries)
1175{
1176 // generate new THnSparseF, axes are defined in GetDimParams()
1177
1178 Int_t count = 0;
1179 UInt_t tmp = entries;
1180 while(tmp!=0){
1181 count++;
1182 tmp = tmp &~ -tmp; // clear lowest bit
1183 }
1184
1185 TString hnTitle(name);
1186 const Int_t dim = count;
1187 Int_t nbins[dim];
1188 Double_t xmin[dim];
1189 Double_t xmax[dim];
1190
1191 Int_t i=0;
1192 Int_t c=0;
1193 while(c<dim && i<32){
1194 if(entries&(1<<i)){
1195
1196 TString label("");
1197 GetDimParams(i, label, nbins[c], xmin[c], xmax[c]);
1198 hnTitle += Form(";%s",label.Data());
1199 c++;
1200 }
1201
1202 i++;
1203 }
1204 hnTitle += ";";
1205
1206 return new THnSparseF(name, hnTitle.Data(), dim, nbins, xmin, xmax);
1207}
1208
1209void AliAnalysisTaskJetCore::GetDimParams(Int_t iEntry, TString &label, Int_t &nbins, Double_t &xmin, Double_t &xmax)
1210{
1211 // stores label and binning of axis for THnSparse
1212
1213 const Double_t pi = TMath::Pi();
1214
1215 switch(iEntry){
1216
1217 case 0:
1218 label = "V0 centrality (%)";
1219
1220 nbins = 10;
1221 xmin = 0.;
1222 xmax = 100.;
1223 break;
1224
1225
1226 case 1:
1227 label = "corrected jet pt";
c42c6838 1228 nbins = 20;
da93bb11 1229 xmin = 0.;
1230 xmax = 200.;
1231 break;
1232
1233
1234 case 2:
1235 label = "track pT";
1236
ba143e7f 1237 nbins = 9;
da93bb11 1238 xmin = 0.;
5bd732d4 1239 xmax = 150;
da93bb11 1240 break;
1241
1242
529e2916 1243 case 3:
da93bb11 1244 label = "deltaR";
1245 nbins = 15;
1246 xmin = 0.;
1247 xmax = 1.5;
1248 break;
529e2916 1249
1250
1251
da93bb11 1252 case 4:
1253 label = "deltaEta";
ba143e7f 1254 nbins = 8;
1255 xmin = -1.6;
1256 xmax = 1.6;
da93bb11 1257 break;
1258
1259
5bd732d4 1260 case 5:
da93bb11 1261 label = "deltaPhi";
1262 nbins = 90;
1263 xmin = -0.5*pi;
1264 xmax = 1.5*pi;
1265 break;
1266
1267
529e2916 1268
1269 case 6:
da93bb11 1270 label = "leading track";
ba143e7f 1271 nbins = 13;
da93bb11 1272 xmin = 0;
20dcc500 1273 xmax = 50;
da93bb11 1274 break;
1275
529e2916 1276 case 7:
da93bb11 1277
1278 label = "trigger track";
ba143e7f 1279 nbins =10;
da93bb11 1280 xmin = 0;
1281 xmax = 50;
1282 break;
529e2916 1283
1284
1285
1286
1287
1288
1289
1290
da93bb11 1291 }
1292
1293}
1294