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