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