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