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