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