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