o First Version of TRDnSigma implementation (Xianguo) o still requires some catching...
[u/mrichter/AliRoot.git] / PWGJE / AliAnalysisTaskJetSpectrum2.cxx
CommitLineData
3b7ffecf 1// **************************************
3170a3f8 2// used for the correction of determiantion of reconstructed jet spectra
3b7ffecf 3// Compares input (gen) and output (rec) jets
4// *******************************************
5
6
7/**************************************************************************
8 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
9 * *
10 * Author: The ALICE Off-line Project. *
11 * Contributors are mentioned in the code where appropriate. *
12 * *
13 * Permission to use, copy, modify and distribute this software and its *
14 * documentation strictly for non-commercial purposes is hereby granted *
15 * without fee, provided that the above copyright notice appears in all *
16 * copies and that both the copyright notice and this permission notice *
17 * appear in the supporting documentation. The authors make no claims *
18 * about the suitability of this software for any purpose. It is *
19 * provided "as is" without express or implied warranty. *
20 **************************************************************************/
3b7ffecf 21
22#include <TROOT.h>
23#include <TRandom.h>
742ee86c 24#include <TRandom3.h>
3b7ffecf 25#include <TSystem.h>
26#include <TInterpreter.h>
27#include <TChain.h>
28#include <TFile.h>
29#include <TKey.h>
30#include <TH1F.h>
31#include <TH2F.h>
32#include <TH3F.h>
33#include <TProfile.h>
2eded560 34#include <TProfile2D.h>
3b7ffecf 35#include <TList.h>
36#include <TLorentzVector.h>
37#include <TClonesArray.h>
37eb26ea 38#include <TRefArray.h>
3b7ffecf 39#include "TDatabasePDG.h"
40
41#include "AliAnalysisTaskJetSpectrum2.h"
42#include "AliAnalysisManager.h"
43#include "AliJetFinder.h"
b99018c6 44#include "AliTHn.h"
3b7ffecf 45#include "AliJetHeader.h"
46#include "AliJetReader.h"
47#include "AliJetReaderHeader.h"
48#include "AliUA1JetHeaderV1.h"
49#include "AliESDEvent.h"
50#include "AliAODEvent.h"
51#include "AliAODHandler.h"
52#include "AliAODTrack.h"
53#include "AliAODJet.h"
c2785065 54#include "AliAODJetEventBackground.h"
3b7ffecf 55#include "AliAODMCParticle.h"
56#include "AliMCEventHandler.h"
57#include "AliMCEvent.h"
58#include "AliStack.h"
59#include "AliGenPythiaEventHeader.h"
60#include "AliJetKineReaderHeader.h"
61#include "AliGenCocktailEventHeader.h"
62#include "AliInputEventHandler.h"
152a0744 63#include "AliCFContainer.h"
3b7ffecf 64
65#include "AliAnalysisHelperJetTasks.h"
66
67ClassImp(AliAnalysisTaskJetSpectrum2)
68
37eb26ea 69AliAnalysisTaskJetSpectrum2::AliAnalysisTaskJetSpectrum2():
adf3920d 70AliAnalysisTaskSE(),
a31b8a87 71 fJetHeaderRec(0x0),
72 fJetHeaderGen(0x0),
73 fAODIn(0x0),
74 fAODOut(0x0),
75 fAODExtension(0x0),
b99018c6 76 fhnJetContainer(0x0),
a31b8a87 77 fhnCorrelation(0x0),
90a006c2 78 fhnEvent(0x0),
c8eabe24 79 f1PtScale(0x0),
a31b8a87 80 fBranchRec("jets"),
81 fBranchGen(""),
37eb26ea 82 fBranchBkgRec(""),
83 fBranchBkgGen(""),
a31b8a87 84 fNonStdFile(""),
742ee86c 85 fRandomizer(0x0),
a31b8a87 86 fUseAODJetInput(kFALSE),
87 fUseAODTrackInput(kFALSE),
88 fUseAODMCInput(kFALSE),
89 fUseGlobalSelection(kFALSE),
90 fUseExternalWeightOnly(kFALSE),
91 fLimitGenJetEta(kFALSE),
b99018c6 92 fDoMatching(kFALSE),
a31b8a87 93 fNMatchJets(5),
7b822bc5 94 fNRPBins(3),
4f6d1cc1 95 fTRP(0),
96 fDebug(0),
65c43de8 97 fJetTriggerExcludeMask(AliAODJet::kHighTrackPtTriggered),
d7396b04 98 fJetTriggerBestMask(AliAODJet::kHighTrackPtBest),
a31b8a87 99 fFilterMask(0),
f2dd0695 100 fEventSelectionMask(0),
3d1dba29 101 fNTrigger(0),
102 fTriggerBit(0x0),
db745d29 103 fNAcceptance(0),
85b5b73e 104 fNBinsLeadingTrackPt(10),
105 fNBinsMult(20),
a31b8a87 106 fAnalysisType(0),
3b7ffecf 107 fTrackTypeRec(kTrackUndef),
108 fTrackTypeGen(kTrackUndef),
f4132e7d 109 fEventClass(0),
2eded560 110 fRPMethod(0),
3b7ffecf 111 fAvgTrials(1),
112 fExternalWeight(1),
a31b8a87 113 fJetRecEtaWindow(0.5),
114 fTrackRecEtaWindow(0.5),
9280dfa6 115 fMinJetPt(0),
a31b8a87 116 fMinTrackPt(0.15),
117 fDeltaPhiWindow(90./180.*TMath::Pi()),
69a0d320 118 fAcceptancePhiMin(0x0),
119 fAcceptancePhiMax(0x0),
120 fAcceptanceEtaMin(0x0),
121 fAcceptanceEtaMax(0x0),
742ee86c 122 fCentrality(100),
7b822bc5 123 fRPAngle(0),
37eb26ea 124 fMultRec(0),
125 fMultGen(0),
3d1dba29 126 fTriggerName(0x0),
3b7ffecf 127 fh1Xsec(0x0),
128 fh1Trials(0x0),
432abb2b 129 fh1AvgTrials(0x0),
3b7ffecf 130 fh1PtHard(0x0),
131 fh1PtHardNoW(0x0),
132 fh1PtHardTrials(0x0),
57ca1193 133 fh1ZVtx(0x0),
7b822bc5 134 fh1RP(0x0),
742ee86c 135 fh1Centrality(0x0),
c8eabe24 136 fh1TmpRho(0x0),
37eb26ea 137 fh2MultRec(0x0),
138 fh2MultGen(0x0),
742ee86c 139 fh2RPCentrality(0x0),
a31b8a87 140 fh2PtFGen(0x0),
e39a4f22 141 fh2deltaPt1Pt2(0x0),
9a83d4af 142 fh2RelPtFGen(0x0),
346e5a93 143 fh3RelPtFGenLeadTrkPt(0x0),
0c964c83 144 fHistList(0x0),
145 fh1EvtSelection(0),
146 fMaxVertexZ(100.),
147 fMinNcontributors(0),
148 fRejectPileup(0)
149
3b7ffecf 150{
c17f867b 151
a31b8a87 152 for(int ij = 0;ij <kJetTypes;++ij){
742ee86c 153 fFlagJetType[ij] = 1; // default = on
a31b8a87 154 fh1NJets[ij] = 0;
155 fh1SumPtTrack[ij] = 0;
156 fh1PtJetsIn[ij] = 0;
65c43de8 157 fh1PtJetsInRej[ij] = 0;
b01d3ab7 158 fh1PtJetsInBest[ij] = 0;
a31b8a87 159 fh1PtTracksIn[ij] = 0;
cb883243 160 fh1PtTracksInLow[ij] = 0;
a31b8a87 161 fh2NJetsPt[ij] = 0;
162 fh2NTracksPt[ij] = 0;
2eded560 163 fp2MultRPPhiTrackPt[ij] = 0;
164 fp2CentRPPhiTrackPt[ij] = 0;
90a006c2 165 fhnJetPt[ij] = 0;
17f4943e 166 fhnJetPtBest[ij] = 0;
167 fhnJetPtRej[ij] = 0;
90a006c2 168 fhnJetPtQA[ij] = 0;
169 fhnTrackPt[ij] = 0;
170 fhnTrackPtQA[ij] = 0;
37eb26ea 171 for(int i = 0;i <= kMaxJets;++i){
d7117cbd 172 fh2LTrackPtJetPt[ij][i] = 0;
a31b8a87 173 fh1PtIn[ij][i] = 0;
a31b8a87 174 }
175
176 fh1DijetMinv[ij] = 0;
177 fh2DijetDeltaPhiPt[ij] = 0;
178 fh2DijetAsymPt[ij] = 0;
179 fh2DijetPt2vsPt1[ij] = 0;
180 fh2DijetDifvsSum[ij] = 0;
181
182 }
3b7ffecf 183}
184
185AliAnalysisTaskJetSpectrum2::AliAnalysisTaskJetSpectrum2(const char* name):
186 AliAnalysisTaskSE(name),
187 fJetHeaderRec(0x0),
188 fJetHeaderGen(0x0),
a31b8a87 189 fAODIn(0x0),
190 fAODOut(0x0),
191 fAODExtension(0x0),
b99018c6 192 fhnJetContainer(0x0),
3b7ffecf 193 fhnCorrelation(0x0),
90a006c2 194 fhnEvent(0x0),
c8eabe24 195 f1PtScale(0x0),
3b7ffecf 196 fBranchRec("jets"),
197 fBranchGen(""),
37eb26ea 198 fBranchBkgRec(""),
199 fBranchBkgGen(""),
a31b8a87 200 fNonStdFile(""),
742ee86c 201 fRandomizer(0x0),
565584e8 202 fUseAODJetInput(kFALSE),
203 fUseAODTrackInput(kFALSE),
204 fUseAODMCInput(kFALSE),
b5cc0c6d 205 fUseGlobalSelection(kFALSE),
3b7ffecf 206 fUseExternalWeightOnly(kFALSE),
207 fLimitGenJetEta(kFALSE),
b99018c6 208 fDoMatching(kFALSE),
a31b8a87 209 fNMatchJets(5),
7b822bc5 210 fNRPBins(3),
4f6d1cc1 211 fTRP(0),
212 fDebug(0),
65c43de8 213 fJetTriggerExcludeMask(AliAODJet::kHighTrackPtTriggered),
d7396b04 214 fJetTriggerBestMask(AliAODJet::kHighTrackPtBest),
3b7ffecf 215 fFilterMask(0),
f2dd0695 216 fEventSelectionMask(0),
3d1dba29 217 fNTrigger(0),
218 fTriggerBit(0x0),
db745d29 219 fNAcceptance(0),
85b5b73e 220 fNBinsLeadingTrackPt(10),
221 fNBinsMult(20),
3b7ffecf 222 fAnalysisType(0),
223 fTrackTypeRec(kTrackUndef),
224 fTrackTypeGen(kTrackUndef),
f4132e7d 225 fEventClass(0),
2eded560 226 fRPMethod(0),
3b7ffecf 227 fAvgTrials(1),
228 fExternalWeight(1),
a31b8a87 229 fJetRecEtaWindow(0.5),
230 fTrackRecEtaWindow(0.5),
9280dfa6 231 fMinJetPt(0),
a31b8a87 232 fMinTrackPt(0.15),
233 fDeltaPhiWindow(90./180.*TMath::Pi()),
69a0d320 234 fAcceptancePhiMin(0x0),
235 fAcceptancePhiMax(0x0),
236 fAcceptanceEtaMin(0x0),
237 fAcceptanceEtaMax(0x0),
742ee86c 238 fCentrality(100),
7b822bc5 239 fRPAngle(0),
37eb26ea 240 fMultRec(0),
241 fMultGen(0),
3d1dba29 242 fTriggerName(0x0),
3b7ffecf 243 fh1Xsec(0x0),
244 fh1Trials(0x0),
432abb2b 245 fh1AvgTrials(0x0),
3b7ffecf 246 fh1PtHard(0x0),
247 fh1PtHardNoW(0x0),
248 fh1PtHardTrials(0x0),
57ca1193 249 fh1ZVtx(0x0),
7b822bc5 250 fh1RP(0x0),
742ee86c 251 fh1Centrality(0x0),
c8eabe24 252 fh1TmpRho(0x0),
37eb26ea 253 fh2MultRec(0x0),
254 fh2MultGen(0x0),
742ee86c 255 fh2RPCentrality(0x0),
a31b8a87 256 fh2PtFGen(0x0),
e39a4f22 257 fh2deltaPt1Pt2(0x0),
9a83d4af 258 fh2RelPtFGen(0x0),
346e5a93 259 fh3RelPtFGenLeadTrkPt(0x0),
0c964c83 260 fHistList(0x0),
261 fh1EvtSelection(0),
262 fMaxVertexZ(100.),
263 fMinNcontributors(0),
264 fRejectPileup(0)
3b7ffecf 265{
266
a31b8a87 267 for(int ij = 0;ij <kJetTypes;++ij){
742ee86c 268 fFlagJetType[ij] = 1; // default = on
a31b8a87 269 fh1NJets[ij] = 0;
270 fh1SumPtTrack[ij] = 0;
271 fh1PtJetsIn[ij] = 0;
65c43de8 272 fh1PtJetsInRej[ij] = 0;
b01d3ab7 273 fh1PtJetsInBest[ij] = 0;
a31b8a87 274 fh1PtTracksIn[ij] = 0;
cb883243 275 fh1PtTracksInLow[ij] = 0;
a31b8a87 276 fh2NJetsPt[ij] = 0;
277 fh2NTracksPt[ij] = 0;
2eded560 278 fp2MultRPPhiTrackPt[ij] = 0;
279 fp2CentRPPhiTrackPt[ij] = 0;
90a006c2 280 fhnJetPt[ij] = 0;
17f4943e 281 fhnJetPtBest[ij] = 0;
282 fhnJetPtRej[ij] = 0;
90a006c2 283 fhnJetPtQA[ij] = 0;
284 fhnTrackPt[ij] = 0;
285 fhnTrackPtQA[ij] = 0;
37eb26ea 286 for(int i = 0;i <= kMaxJets;++i){
d7117cbd 287 fh2LTrackPtJetPt[ij][i] = 0;
a31b8a87 288 fh1PtIn[ij][i] = 0;
a31b8a87 289 }
290
291 fh1DijetMinv[ij] = 0;
292 fh2DijetDeltaPhiPt[ij] = 0;
293 fh2DijetAsymPt[ij] = 0;
294 fh2DijetPt2vsPt1[ij] = 0;
295 fh2DijetDifvsSum[ij] = 0;
296 }
297
37eb26ea 298 DefineOutput(1, TList::Class());
3b7ffecf 299}
300
301
302
303Bool_t AliAnalysisTaskJetSpectrum2::Notify()
304{
a31b8a87 305
306
307
3b7ffecf 308 //
309 // Implemented Notify() to read the cross sections
310 // and number of trials from pyxsec.root
311 //
a31b8a87 312
313 // Fetch the aod also from the input in,
314 // have todo it in notify
315
316
317 fAODIn = dynamic_cast<AliAODEvent*>(InputEvent());
318 // assume that the AOD is in the general output...
319 fAODOut = AODEvent();
320
321 if(fNonStdFile.Length()!=0){
322 // case that we have an AOD extension we need can fetch the jets from the extended output
323 AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
e855f5c5 324 fAODExtension = (aodH?aodH->GetExtension(fNonStdFile.Data()):0);
a31b8a87 325 if(!fAODExtension){
326 if(fDebug>1)Printf("AODExtension found for %s",fNonStdFile.Data());
327 }
328 }
329
3b7ffecf 330
331 TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
332 Float_t xsection = 0;
333 Float_t ftrials = 1;
334
335 fAvgTrials = 1;
336 if(tree){
337 TFile *curfile = tree->GetCurrentFile();
338 if (!curfile) {
339 Error("Notify","No current file");
340 return kFALSE;
341 }
342 if(!fh1Xsec||!fh1Trials){
346cfe4a 343 Printf("%s%d No Histogram fh1Xsec",(char*)__FILE__,__LINE__);
3b7ffecf 344 return kFALSE;
345 }
346 AliAnalysisHelperJetTasks::PythiaInfoFromFile(curfile->GetName(),xsection,ftrials);
347 fh1Xsec->Fill("<#sigma>",xsection);
348 // construct a poor man average trials
349 Float_t nEntries = (Float_t)tree->GetTree()->GetEntries();
a221560c 350 if(ftrials>=nEntries && nEntries>0.)fAvgTrials = ftrials/nEntries;
432abb2b 351 fh1Trials->Fill("#sum{ntrials}",ftrials);
3b7ffecf 352 }
a31b8a87 353
d7117cbd 354 if(fDebug)Printf("Reading File %s",fInputHandler->GetTree()->GetCurrentFile()->GetName());
a31b8a87 355
3b7ffecf 356 return kTRUE;
357}
358
359void AliAnalysisTaskJetSpectrum2::UserCreateOutputObjects()
360{
361
742ee86c 362
3b7ffecf 363 // Connect the AOD
364
3b7ffecf 365 if (fDebug > 1) printf("AnalysisTaskJetSpectrum2::UserCreateOutputObjects() \n");
3b7ffecf 366 OpenFile(1);
37eb26ea 367 if(!fHistList)fHistList = new TList();
40440b28 368 PostData(1, fHistList); // post data in any case once
369
742ee86c 370 if(!fRandomizer)fRandomizer = new TRandom3(0);
371
5cb0f438 372 fHistList->SetOwner(kTRUE);
37eb26ea 373 Bool_t oldStatus = TH1::AddDirectoryStatus();
3b7ffecf 374 TH1::AddDirectory(kFALSE);
375
90a006c2 376
377
378 // event npsparse cent, mult
e5ee1c19 379 const Int_t nBinsSparse0 = 4;
380 const Int_t nBins0[nBinsSparse0] = { 100, 500,fNTrigger,125};
381 const Double_t xmin0[nBinsSparse0] = { 0, 0, -0.5,-2};
382 const Double_t xmax0[nBinsSparse0] = { 100,5000,fNTrigger-0.5,248};
90a006c2 383
384
e5ee1c19 385 fhnEvent = new THnSparseF("fhnEvent",";cent;mult:trigger;#rho",nBinsSparse0,
90a006c2 386 nBins0,xmin0,xmax0);
387 fHistList->Add(fhnEvent);
388
b99018c6 389 if(fDoMatching){
390 MakeJetContainer();
391 fHistList->Add(fhnCorrelation);
392 fHistList->Add(fhnJetContainer);
393 }
3b7ffecf 394 //
395 // Histogram
396
90a006c2 397
0c964c83 398 fh1EvtSelection = new TH1F("fh1EvtSelection", "Event Selection", 7, -0.5, 6.5);
399 fh1EvtSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
400 fh1EvtSelection->GetXaxis()->SetBinLabel(2,"event selection: rejected");
401 fh1EvtSelection->GetXaxis()->SetBinLabel(3,"event class: rejected");
402 fh1EvtSelection->GetXaxis()->SetBinLabel(4,"vertex Ncontr: rejected");
403 fh1EvtSelection->GetXaxis()->SetBinLabel(5,"vertex z: rejected");
404 fh1EvtSelection->GetXaxis()->SetBinLabel(6,"vertex type: rejected");
405 fh1EvtSelection->GetXaxis()->SetBinLabel(7,"pileup: rejected");
90a006c2 406
0c964c83 407 fHistList->Add(fh1EvtSelection);
408
66fd4b2b 409 const Int_t nBinPt = 150;
3b7ffecf 410 Double_t binLimitsPt[nBinPt+1];
411 for(Int_t iPt = 0;iPt <= nBinPt;iPt++){
412 if(iPt == 0){
ed5c9f3f 413 binLimitsPt[iPt] = -50.0;
3b7ffecf 414 }
415 else {// 1.0
66fd4b2b 416 binLimitsPt[iPt] = binLimitsPt[iPt-1] + 2.;
3b7ffecf 417 }
418 }
edfbe476 419 const Int_t nBinPhi = 90;
3b7ffecf 420 Double_t binLimitsPhi[nBinPhi+1];
421 for(Int_t iPhi = 0;iPhi<=nBinPhi;iPhi++){
422 if(iPhi==0){
37eb26ea 423 binLimitsPhi[iPhi] = 0;
3b7ffecf 424 }
425 else{
426 binLimitsPhi[iPhi] = binLimitsPhi[iPhi-1] + 1/(Float_t)nBinPhi * TMath::Pi()*2;
427 }
428 }
429
edfbe476 430 const Int_t nBinEta = 40;
431 Double_t binLimitsEta[nBinEta+1];
432 for(Int_t iEta = 0;iEta<=nBinEta;iEta++){
433 if(iEta==0){
434 binLimitsEta[iEta] = -2.0;
435 }
436 else{
cf049e94 437 binLimitsEta[iEta] = binLimitsEta[iEta-1] + 0.1;
edfbe476 438 }
439 }
440
a31b8a87 441
3b7ffecf 442 fh1Xsec = new TProfile("fh1Xsec","xsec from pyxsec.root",1,0,1);
443 fh1Xsec->GetXaxis()->SetBinLabel(1,"<#sigma>");
90a006c2 444 fHistList->Add(fh1Xsec);
3b7ffecf 445 fh1Trials = new TH1F("fh1Trials","trials root file",1,0,1);
446 fh1Trials->GetXaxis()->SetBinLabel(1,"#sum{ntrials}");
a31b8a87 447 fHistList->Add(fh1Trials);
432abb2b 448 fh1AvgTrials = new TH1F("fh1AvgTrials","trials root file",1,0,1);
449 fh1AvgTrials->GetXaxis()->SetBinLabel(1,"#sum{avg ntrials}");
450 fHistList->Add(fh1AvgTrials);
3b7ffecf 451 fh1PtHard = new TH1F("fh1PtHard","PYTHIA Pt hard;p_{T,hard}",nBinPt,binLimitsPt);
a31b8a87 452 fHistList->Add(fh1PtHard);
3b7ffecf 453 fh1PtHardNoW = new TH1F("fh1PtHardNoW","PYTHIA Pt hard no weight;p_{T,hard}",nBinPt,binLimitsPt);
a31b8a87 454 fHistList->Add(fh1PtHardNoW);
3b7ffecf 455 fh1PtHardTrials = new TH1F("fh1PtHardTrials","PYTHIA Pt hard weight with trials;p_{T,hard}",nBinPt,binLimitsPt);
a31b8a87 456 fHistList->Add(fh1PtHardTrials);
a31b8a87 457
57ca1193 458 fh1ZVtx = new TH1F("fh1ZVtx","z vtx;z_{vtx} (cm)",400,-20,20);
a31b8a87 459 fHistList->Add(fh1ZVtx);
7b822bc5 460
2eded560 461
9939142f 462 fh1RP = new TH1F("fh1RP","RP;#Psi",440, -1.*TMath::Pi(), 2.*TMath::Pi());
7b822bc5 463 fHistList->Add(fh1RP);
464
eb1739b7 465 fh1Centrality = new TH1F("fh1Centrality","cent;cent (%)",103,-1,102);
742ee86c 466 fHistList->Add(fh1Centrality);
467
3170a3f8 468 fh2MultRec = new TH2F("fh2MultRec","multiplicity rec;# tracks;# jetrefs",500,0,5000,500,0.,5000);
37eb26ea 469 fHistList->Add(fh2MultRec);
3170a3f8 470 fh2MultGen = new TH2F("fh2MultGen","multiplicity gen;# tracks;# jetrefs",400,0,5000,500,0.,5000);
37eb26ea 471 fHistList->Add(fh2MultGen);
a31b8a87 472
a31b8a87 473
742ee86c 474 fh2RPCentrality = new TH2F("fh2RPCentrality" ,"Reaction Plane Angle" , 20, 0.,100., 180, 0, TMath::Pi());
475 fHistList->Add(fh2RPCentrality);
a31b8a87 476
742ee86c 477 fh2PtFGen = new TH2F("fh2PtFGen",Form("%s vs. %s;p_{T,gen};p_{T,rec}",fBranchRec.Data(),fBranchGen.Data()),nBinPt,binLimitsPt,nBinPt,binLimitsPt);
478 fHistList->Add(fh2PtFGen);
a31b8a87 479
e39a4f22 480 // fh2deltaPt1Pt2(0x0),
481 fh2deltaPt1Pt2 = new TH3F("fh2deltaPt1Pt2",Form("%s vs. %s;delta;p_{T,gen};p_{T,rec}",fBranchRec.Data(),fBranchGen.Data()),nBinPt,binLimitsPt,nBinPt,binLimitsPt,nBinPt,binLimitsPt);
482 fHistList->Add(fh2deltaPt1Pt2);
483
742ee86c 484 fh2RelPtFGen = new TH2F("fh2RelPtFGen",";p_{T,gen};p_{T,rec}-p_{T,gen}/p_{T,Gen}",nBinPt,binLimitsPt,241,-2.41,2.41);
485 fHistList->Add(fh2RelPtFGen);
346e5a93 486
487 const Int_t nBinsLeadingTrackPt = 10;
488 const Double_t binArrayLeadingTrackPt[nBinsLeadingTrackPt+1] = {0.,1.,2.,3.,4.,5.,6.,8.,10.,12.,200.}; //store pT of leading track in jet
489
490 const Int_t nBinDeltaPt = 241;
491 Double_t binLimitsDeltaPt[nBinDeltaPt+1];
492 for(Int_t iDeltaPt = 0;iDeltaPt <= nBinDeltaPt;iDeltaPt++){
493 if(iDeltaPt == 0){
494 binLimitsDeltaPt[iDeltaPt] = -2.41;
495 }
496 else {
497 binLimitsDeltaPt[iDeltaPt] = binLimitsDeltaPt[iDeltaPt-1] + 0.02;
498 }
499 }
500
501 fh3RelPtFGenLeadTrkPt = new TH3F("fh3RelPtFGenLeadTrkPt",";p_{T,gen};p_{T,rec}-p_{T,gen}/p_{T,Gen};p_{T}^{leading track}",nBinPt,binLimitsPt,nBinDeltaPt,binLimitsDeltaPt,nBinsLeadingTrackPt,binArrayLeadingTrackPt);
502 fHistList->Add(fh3RelPtFGenLeadTrkPt);
503
742ee86c 504 for(int ij = 0;ij <kJetTypes;++ij){
505 TString cAdd = "";
506 TString cJetBranch = "";
507 if(ij==kJetRec){
508 cAdd = "Rec";
509 cJetBranch = fBranchRec.Data();
510 }
511 else if (ij==kJetGen){
512 cAdd = "Gen";
513 cJetBranch = fBranchGen.Data();
514 }
515 else if (ij==kJetRecFull){
516 cAdd = "RecFull";
517 cJetBranch = fBranchRec.Data();
518 }
519 else if (ij==kJetGenFull){
520 cAdd = "GenFull";
521 cJetBranch = fBranchGen.Data();
522 }
a31b8a87 523
742ee86c 524 if(cJetBranch.Length()==0)fFlagJetType[ij] = 0;
525 if(!fFlagJetType[ij])continue;
7b822bc5 526
742ee86c 527 fh1NJets[ij] =new TH1F(Form("fh1N%sJets",cAdd.Data()),Form("N %s jets",cAdd.Data()),50,-0.5,49.5);
528 fHistList->Add(fh1NJets[ij]);
529
530 fh1PtJetsIn[ij] = new TH1F(Form("fh1PtJets%sIn",cAdd.Data()),Form("%s jets p_T;p_{T} (GeV/c)",cAdd.Data()),nBinPt,binLimitsPt);
531 fHistList->Add(fh1PtJetsIn[ij]);
65c43de8 532
533 fh1PtJetsInRej[ij] = new TH1F(Form("fh1PtJets%sInRej",cAdd.Data()),Form("%s jets p_T;p_{T} (GeV/c)",cAdd.Data()),nBinPt,binLimitsPt);
534 fHistList->Add(fh1PtJetsInRej[ij]);
d7396b04 535
536 fh1PtJetsInBest[ij] = new TH1F(Form("fh1PtJets%sInBest",cAdd.Data()),Form("%s jets p_T;p_{T} (GeV/c)",cAdd.Data()),nBinPt,binLimitsPt);
537 fHistList->Add(fh1PtJetsInBest[ij]);
742ee86c 538
539 fh1PtTracksIn[ij] = new TH1F(Form("fh1PtTracks%sIn",cAdd.Data()),Form("%s track p_T;p_{T} (GeV/c)",cAdd.Data()),nBinPt,binLimitsPt);
540 fHistList->Add(fh1PtTracksIn[ij]);
541
2eded560 542 fh1PtTracksInLow[ij] = new TH1F(Form("fh1PtTracks%sInLow",cAdd.Data()),Form("%s track p_T;p_{T} (GeV/c)",cAdd.Data()),100,0.,5.);
742ee86c 543 fHistList->Add(fh1PtTracksInLow[ij]);
544
3170a3f8 545 fh1SumPtTrack[ij] = new TH1F(Form("fh1SumPtTrack%s",cAdd.Data()),Form("Sum %s track p_T;p_{T} (GeV/c)",cAdd.Data()),1000,0.,3000.);
742ee86c 546 fHistList->Add(fh1SumPtTrack[ij]);
547
548 fh2NJetsPt[ij] = new TH2F(Form("fh2N%sJetsPt",cAdd.Data()),Form("Number of %s jets above threshhold;p_{T,cut} (GeV/c);N_{jets}",cAdd.Data()),nBinPt,binLimitsPt,50,-0.5,49.5);
549 fHistList->Add(fh2NJetsPt[ij]);
550
3170a3f8 551 fh2NTracksPt[ij] = new TH2F(Form("fh2N%sTracksPt",cAdd.Data()),Form("Number of %s tracks above threshhold;p_{T,cut} (GeV/c);N_{tracks}",cAdd.Data()),nBinPt,binLimitsPt,1000,0.,5000);
742ee86c 552 fHistList->Add(fh2NTracksPt[ij]);
90a006c2 553
2eded560 554
555 fp2MultRPPhiTrackPt[ij] = new TProfile2D(Form("fp2MultRPPhiTrackPt%s",cAdd.Data()),"RP phi vs Number of tracks;# tracks;#Delta#phi_{RP}; <p_{T}>",20,0,4000,181,-1./180.*TMath::Pi(),TMath::Pi(),"S");
556 fHistList->Add(fp2MultRPPhiTrackPt[ij]);
557 fp2CentRPPhiTrackPt[ij] = new TProfile2D(Form("fp2CentRPPhiTrackPt%s",cAdd.Data()),"RP phi vs cent;# cent;#Delta#phi_{RP}; <p_{T}>",10,0,100,181,-1./180.*TMath::Pi(),TMath::Pi(),"S");
558 fHistList->Add(fp2CentRPPhiTrackPt[ij]);
559
556ccab2 560 // Bins: Jet number: pTJet, cent, mult, RP, Area, trigger, leading track pT total bins = 4.5M
db745d29 561 const Int_t nBinsSparse1 = 9;
152a0744 562 const Int_t nBinsArea = 10;
4f6d1cc1 563 Int_t nbinr5=120;
564 Int_t nbinr5min=-50;
565
566 if(cJetBranch.Contains("KT05")){
567 nbinr5=132;
568 nbinr5min=-80;}
569
570 Int_t nBins1[nBinsSparse1] = { kMaxJets+1,nbinr5, 10, fNBinsMult, fNRPBins, nBinsArea,fNTrigger,fNBinsLeadingTrackPt,fNAcceptance+1};
d01931c3 571 if(cJetBranch.Contains("RandomCone")){
572 nBins1[1] = 600;
573 nBins1[5] = 1;
574 }
2942f542 575 const Double_t xmin1[nBinsSparse1] = { -0.5, static_cast<Double_t>(nbinr5min), 0, 0, -0.5, 0., -0.5, 0., -0.5,};
db745d29 576 const Double_t xmax1[nBinsSparse1] = {kMaxJets+0.5,250,100,4000,fNRPBins-0.5,1.0,fNTrigger-0.5,200.,fNAcceptance+0.5};
4f6d1cc1 577
adf3920d 578 const Double_t binArrayArea[nBinsArea+1] = {xmin1[5],0.07,0.1,0.15,0.2,0.25,0.3,0.4,0.5,0.6,xmax1[5]};
346e5a93 579
7bee4353 580
db745d29 581 fhnJetPt[ij] = new THnSparseF(Form("fhnJetPt%s",cAdd.Data()),";jet number;p_{T,jet};cent;# tracks;RP;area;trigger;leading track p_{T};acceptance bin",nBinsSparse1,nBins1,xmin1,xmax1);
152a0744 582 fhnJetPt[ij]->SetBinEdges(5,binArrayArea);
85b5b73e 583 if(fNBinsLeadingTrackPt>1) fhnJetPt[ij]->SetBinEdges(7,binArrayLeadingTrackPt);
2eded560 584 fHistList->Add(fhnJetPt[ij]);
17f4943e 585
586
587 // Bins: pTJet, cent, trigger,
588 const Int_t nBinsSparse1b = 3;
589 Int_t nBins1b[nBinsSparse1b] = {120, 10,fNTrigger};
590 const Double_t xmin1b[nBinsSparse1b] = {-50, 0,-0.5};
591 const Double_t xmax1b[nBinsSparse1b] = {250,100,fNTrigger-0.5};
592
593 fhnJetPtBest[ij] = new THnSparseF(Form("fhnJetPtBest%s",cAdd.Data()),";p_{T,jet};cent;trigger",nBinsSparse1b,nBins1b,xmin1b,xmax1b);
594 fHistList->Add(fhnJetPtBest[ij]);
595
596 fhnJetPtRej[ij] = new THnSparseF(Form("fhnJetPtRej%s",cAdd.Data()),";p_{T,jet};cent;trigger",nBinsSparse1b,nBins1b,xmin1b,xmax1b);
597 fHistList->Add(fhnJetPtRej[ij]);
2eded560 598
adf3920d 599 // Bins: Jet number: pTJet, cent, eta, phi, Area, trigger, acceptance, signed pT leading
5662640e 600 const Int_t nBinsSparse2 = 9;
fedb5a47 601 Int_t nBins2[nBinsSparse2] = { kMaxJets+1, 60, 8, 6, 90, 100,fNTrigger,fNAcceptance+1,20};
d01931c3 602 if(cJetBranch.Contains("RandomCone")){
603 nBins2[5] = 1;
604 }
adf3920d 605 const Double_t xmin2[nBinsSparse2] = { -0.5, -50, 0,-0.9, 0, 0., -0.5, -0.5, -100};
606 const Double_t xmax2[nBinsSparse2] = {kMaxJets+0.5, 250, 80, 0.9, 2.*TMath::Pi(),1.0,fNTrigger-0.5,fNAcceptance+0.5, 100};
5662640e 607 fhnJetPtQA[ij] = new THnSparseF(Form("fhnJetPtQA%s",cAdd.Data()),";jet number;p_{T,jet};cent;#eta;#phi;area;trigger;acceptance bin;signed pt leading",nBinsSparse2,nBins2,xmin2,xmax2);
e39a4f22 608 // fhnJetPt[ij]->SetBinEdges(5,binArrayArea);
d01931c3 609 fHistList->Add(fhnJetPtQA[ij]);
610
5a7c81dc 611 // Bins:track number pTtrack, cent, mult, RP. charge total bins = 224 k
612 const Int_t nBinsSparse3 = 7;
85b5b73e 613 const Int_t nRPBinsSparse3 = 3;
5a7c81dc 614 const Int_t nBins3[nBinsSparse3] = { 2, 100, 10, 1, nRPBinsSparse3,fNTrigger,2};
615 const Double_t xmin3[nBinsSparse3] = { -0.5, 0, 0, 0, -0.5,-0.5,-1.5};
616 const Double_t xmax3[nBinsSparse3] = { 1.5, 200, 100, 4000,nRPBinsSparse3-0.5,fNTrigger-0.5,1.5};
d01931c3 617
adf3920d 618 // change the binning of the pT axis:
d01931c3 619 Double_t *xPt3 = new Double_t[nBins3[1]+1];
620 xPt3[0] = 0.;
621 for(int i = 1; i<=nBins3[1];i++){
622 if(xPt3[i-1]<2)xPt3[i] = xPt3[i-1] + 0.05; // 1 - 40
623 else if(xPt3[i-1]<4)xPt3[i] = xPt3[i-1] + 0.2; // 41 - 50
624 else if(xPt3[i-1]<10)xPt3[i] = xPt3[i-1] + 0.5; // 50 - 62
17f4943e 625 else if(xPt3[i-1]<15)xPt3[i] = xPt3[i-1] + 1.; // 62 - 67
626 else if(xPt3[i-1]<20)xPt3[i] = xPt3[i-1] + 2.; // 67 - 72
627 else if(xPt3[i-1]<60)xPt3[i] = xPt3[i-1] + 5; // 72 - 76
628 else xPt3[i] = xPt3[i-1] + 10; // 76 - 100 = 140
d01931c3 629 }
630
adf3920d 631 fhnTrackPt[ij] = new THnSparseF(Form("fhnTrackPt%s",cAdd.Data()),";track number;p_{T};cent;#tracks;RP;trigger",nBinsSparse3,nBins3,xmin3,xmax3);
632 fhnTrackPt[ij]->SetBinEdges(1,xPt3);
633 fHistList->Add(fhnTrackPt[ij]);
634 delete [] xPt3;
635
636 // Track QA bins track nr, pTrack, cent, eta, phi bins 5.4 M
5a7c81dc 637 const Int_t nBinsSparse4 = 6;
638 const Int_t nBins4[nBinsSparse4] = { 2, 50, 10, 20, 180,2};
639 const Double_t xmin4[nBinsSparse4] = { -0.5, 0, 0, -1.0, 0.,-1.5};
640 const Double_t xmax4[nBinsSparse4] = { 1.5,150, 100, 1.0,2.*TMath::Pi(),1.5};
adf3920d 641
642 // change the binning ot the pT axis:
643 Double_t *xPt4 = new Double_t[nBins4[1]+1];
644 xPt4[0] = 0.;
645 for(int i = 1; i<=nBins4[1];i++){
646 if(xPt4[i-1]<2)xPt4[i] = xPt4[i-1] + 0.1;
647 else if(xPt4[i-1]<10)xPt4[i] = xPt4[i-1] + 0.5;
648 else if(xPt4[i-1]<20)xPt4[i] = xPt4[i-1] + 1.;
649 else if(xPt4[i-1]<30)xPt4[i] = xPt4[i-1] + 2.5;
650 else xPt4[i] = xPt4[i-1] + 5.;
651 }
5a7c81dc 652 fhnTrackPtQA[ij] = new THnSparseF(Form("fhnTrackPtQA%s",cAdd.Data()),";track number;p_{T};cent;#eta;#phi;sign",nBinsSparse4,nBins4,xmin4,xmax4);
adf3920d 653 fhnTrackPtQA[ij]->SetBinEdges(1,xPt4);
654 fHistList->Add(fhnTrackPtQA[ij]);
655 delete [] xPt4;
656
657 for(int i = 0;i <= kMaxJets;++i){
658 fh1PtIn[ij][i] = new TH1F(Form("fh1Pt%sIn_j%d",cAdd.Data(),i),Form("%s p_T input ;p_{T}",cAdd.Data()),nBinPt,binLimitsPt);
659 fHistList->Add(fh1PtIn[ij][i]);
a31b8a87 660
661
adf3920d 662 if(!fh1TmpRho)fh1TmpRho = new TH1F("fh1TmpRho","tmp histo for jet shape",40,0.,2);
663 fh2LTrackPtJetPt[ij][i] = new TH2F(Form("fh2LTrackPtJetPt%s_j%d",cAdd.Data(),i),
664 Form("pt of leadin track within a jet vs jet %s;p_{T,lead in jet};p_{T.jet};",
665 cAdd.Data()),
666 200,0.,200.,nBinPt,binLimitsPt);
667 fHistList->Add(fh2LTrackPtJetPt[ij][i]);
668 }
669
a31b8a87 670
adf3920d 671 fh1DijetMinv[ij] = new TH1F(Form("fh1Dijet%sMinv",cAdd.Data()),"Dijet invariant mass;m_{JJ}",nBinPt,binLimitsPt);
672 fHistList->Add(fh1DijetMinv[ij]);
a31b8a87 673
adf3920d 674 fh2DijetDeltaPhiPt[ij] = new TH2F(Form("fh2Dijet%sDeltaPhiPt",cAdd.Data()),"Difference in the azimuthal angle;#Delta#phi;p_{T,2};Entries",180,0.,TMath::Pi(),nBinPt,binLimitsPt);
675 fHistList->Add(fh2DijetDeltaPhiPt[ij]);
a31b8a87 676
adf3920d 677 fh2DijetAsymPt[ij] = new TH2F(Form("fh2Dijet%sAsym",cAdd.Data()),"Pt asymmetry;#Deltap_{T}/(p_{T,1}+p_{T,2});p_{T,1};Entries",50,0.,1.,nBinPt,binLimitsPt);
678 fHistList->Add(fh2DijetAsymPt[ij]);
679
680 fh2DijetPt2vsPt1[ij] = new TH2F(Form("fh2Dijet%sPt2vsPt1",cAdd.Data()),"Pt2 versus Pt1;p_{T,1} (GeV/c);p_{T,2} (GeV/c)",250,0.,250.,250,0.,250.);
681 fHistList->Add(fh2DijetPt2vsPt1[ij]);
682 fh2DijetDifvsSum[ij] = new TH2F(Form("fh2Dijet%sDifvsSum",cAdd.Data()),"Pt difference vs Pt sum;p_{T,1}+p_{T,2} (GeV/c);#Deltap_{T} (GeV/c)",400,0.,400.,150,0.,150.);
683 fHistList->Add( fh2DijetDifvsSum[ij]);
684 }
3b7ffecf 685 // =========== Switch on Sumw2 for all histos ===========
686 for (Int_t i=0; i<fHistList->GetEntries(); ++i) {
687 TH1 *h1 = dynamic_cast<TH1*>(fHistList->At(i));
688 if (h1){
689 h1->Sumw2();
690 continue;
691 }
692 THnSparse *hn = dynamic_cast<THnSparse*>(fHistList->At(i));
693 if(hn)hn->Sumw2();
694 }
695 TH1::AddDirectory(oldStatus);
696}
697
698void AliAnalysisTaskJetSpectrum2::Init()
699{
700 //
701 // Initialization
702 //
703
3b7ffecf 704 if (fDebug > 1) printf("AnalysisTaskJetSpectrum2::Init() \n");
705
706}
707
a31b8a87 708void AliAnalysisTaskJetSpectrum2::UserExec(Option_t */*option*/){
709
b5cc0c6d 710
f2dd0695 711 Bool_t selected = kTRUE;
712
713 if(fUseGlobalSelection&&fEventSelectionMask==0){
714 selected = AliAnalysisHelperJetTasks::Selected();
715 }
716 else if(fUseGlobalSelection&&fEventSelectionMask>0){
45a11aef 717 selected = AliAnalysisHelperJetTasks::TestSelectInfo(fEventSelectionMask);
f2dd0695 718 }
719
0c964c83 720 if(!selected){
721 // no selection by the service task, we continue
722 if (fDebug > 1)Printf("Not selected %s:%d SelectInfo %d Class %d",(char*)__FILE__,__LINE__, AliAnalysisHelperJetTasks::Selected(),AliAnalysisHelperJetTasks::EventClass());
723 fh1EvtSelection->Fill(1.);
724 PostData(1, fHistList);
725 return;
726 }
727
f4132e7d 728 if(fEventClass>0){
729 selected = selected&&(AliAnalysisHelperJetTasks::EventClass()==fEventClass);
730 }
731
f2dd0695 732 if(!selected){
b5cc0c6d 733 // no selection by the service task, we continue
f4132e7d 734 if (fDebug > 1)Printf("Not selected %s:%d SelectInfo %d Class %d",(char*)__FILE__,__LINE__, AliAnalysisHelperJetTasks::Selected(),AliAnalysisHelperJetTasks::EventClass());
0c964c83 735 fh1EvtSelection->Fill(2.);
b5cc0c6d 736 PostData(1, fHistList);
737 return;
738 }
739
f2dd0695 740
a31b8a87 741 static AliAODEvent* aod = 0;
742
743 // take all other information from the aod we take the tracks from
744 if(!aod){
adf3920d 745 if(fUseAODTrackInput)aod = fAODIn;
746 else aod = fAODOut;
a31b8a87 747 }
748
749
3b7ffecf 750 //
751 // Execute analysis for current event
752 //
a31b8a87 753 if (fDebug > 1)printf("Analysing event # %5d\n", (Int_t) fEntry);
754 AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
755 if(!aodH){
756 Printf("%s:%d no output aodHandler found Jet",(char*)__FILE__,__LINE__);
757 return;
758 }
759
760 if (fDebug > 10)Printf("%s:%d",(char*)__FILE__,__LINE__);
761 TClonesArray *aodRecJets = 0;
762
763 if(fAODOut&&!aodRecJets){
764 aodRecJets = dynamic_cast<TClonesArray*>(fAODOut->FindListObject(fBranchRec.Data()));
765 }
766 if(fAODExtension&&!aodRecJets){
767 aodRecJets = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fBranchRec.Data()));
768 }
769 if(fAODIn&&!aodRecJets){
770 aodRecJets = dynamic_cast<TClonesArray*>(fAODIn->FindListObject(fBranchRec.Data()));
771 }
772
773
774
775 if(!aodRecJets){
d95fc15a 776 if(fDebug){
777
adf3920d 778 Printf("%s:%d no reconstructed Jet array with name %s in AOD",(char*)__FILE__,__LINE__,fBranchRec.Data());
779 if(fAODIn){
780 Printf("Input AOD >>>>");
781 fAODIn->Print();
782 }
783 if(fAODExtension){
784 Printf("AOD Extension >>>>");
785 fAODExtension->Print();
786 }
787 if(fAODOut){
788 Printf("Output AOD >>>>");
789 fAODOut->Print();
790 }
d95fc15a 791 }
adf3920d 792 return;
a31b8a87 793 }
794
795 TClonesArray *aodGenJets = 0;
796 if(fBranchGen.Length()>0){
797 if(fAODOut&&!aodGenJets){
798 aodGenJets = dynamic_cast<TClonesArray*>(fAODOut->FindListObject(fBranchGen.Data()));
799 }
800 if(fAODExtension&&!aodGenJets){
801 aodGenJets = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fBranchGen.Data()));
802 }
803 if(fAODIn&&!aodGenJets){
804 aodGenJets = dynamic_cast<TClonesArray*>(fAODIn->FindListObject(fBranchGen.Data()));
805 }
806
807 if(!aodGenJets){
808 Printf("%s:%d no generated Jet array with name %s in AOD",(char*)__FILE__,__LINE__,fBranchGen.Data());
3b7ffecf 809 return;
810 }
3b7ffecf 811 }
a31b8a87 812
37eb26ea 813 TClonesArray *aodBackRecJets = 0;
814 if(fBranchBkgRec.Length()>0){
815 if(fAODOut&&!aodBackRecJets){
816 aodBackRecJets = dynamic_cast<TClonesArray*>(fAODOut->FindListObject(fBranchBkgRec.Data()));
817 }
818 if(fAODExtension&&!aodBackRecJets){
819 aodBackRecJets = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fBranchBkgRec.Data()));
820 }
821 if(fAODIn&&!aodBackRecJets){
822 aodBackRecJets = dynamic_cast<TClonesArray*>(fAODIn->FindListObject(fBranchBkgRec.Data()));
823 }
824
825 if(!aodBackRecJets){
826 Printf("%s:%d no background rec Jet array with name %s in AOD",(char*)__FILE__,__LINE__,fBranchBkgRec.Data());
827 return;
828 }
829 }
830
831
832 TClonesArray *aodBackGenJets = 0;
833
834 if(fBranchBkgGen.Length()>0){
835 if(fAODOut&&!aodBackGenJets){
836 aodBackGenJets = dynamic_cast<TClonesArray*>(fAODOut->FindListObject(fBranchBkgGen.Data()));
837 }
838 if(fAODExtension&&!aodBackGenJets){
839 aodBackGenJets = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fBranchBkgGen.Data()));
840 }
841 if(fAODIn&&!aodBackGenJets){
842 aodBackGenJets = dynamic_cast<TClonesArray*>(fAODIn->FindListObject(fBranchBkgGen.Data()));
843 }
844
845 if(!aodBackGenJets){
846 Printf("%s:%d no background rec Jet array with name %s in AOD",(char*)__FILE__,__LINE__,fBranchBkgGen.Data());
847 return;
848 }
849 }
850
0c964c83 851 AliAODVertex* primVtx = aod->GetPrimaryVertex();
852 Int_t nTracksPrim = primVtx->GetNContributors();
853
854
855 if (fDebug > 1) Printf("%s:%d primary vertex selection: %d", (char*)__FILE__,__LINE__,nTracksPrim);
856 if(nTracksPrim < fMinNcontributors){
857 if (fDebug > 1) Printf("%s:%d primary vertex selection: event REJECTED...",(char*)__FILE__,__LINE__);
858 fh1EvtSelection->Fill(3.);
859 PostData(1, fHistList);
860 return;
861 }
862
863 if(TMath::Abs(primVtx->GetZ())>fMaxVertexZ){
864 if (fDebug > 1) Printf("%s:%d primary vertex z = %f: event REJECTED...",(char*)__FILE__,__LINE__,primVtx->GetZ());
865 fh1EvtSelection->Fill(4.);
866 PostData(1, fHistList);
867 return;
868 }
869
870 TString primVtxName(primVtx->GetName());
871
872 if(primVtxName.CompareTo("TPCVertex",TString::kIgnoreCase) == 1){
873 if (fDebug > 1) Printf("%s:%d primary vertex selection: TPC vertex, event REJECTED...",(char*)__FILE__,__LINE__);
874 fh1EvtSelection->Fill(5.);
875 PostData(1, fHistList);
876 return;
877 }
878
879 if(fRejectPileup && AliAnalysisHelperJetTasks::IsPileUp()){
880 if (fDebug > 1) Printf("%s:%d SPD pileup: event REJECTED...",(char*)__FILE__,__LINE__);
881 fh1EvtSelection->Fill(6.);
882 PostData(1, fHistList);
883 return;
884 }
885
886 if (fDebug > 1) Printf("%s:%d event ACCEPTED ...",(char*)__FILE__,__LINE__);
887 fh1EvtSelection->Fill(0.);
888
a31b8a87 889
890 // new Scheme
891 // first fill all the pure histograms, i.e. full jets
892 // in case of the correaltion limit to the n-leading jets
893
894 // reconstructed
895
896
897 // generated
898
899
900 // second fill the correlation histograms, here we limit to the n-th leading jet in case of the reconstructed
901
902
742ee86c 903
a31b8a87 904 TList genJetsList; // full acceptance
905 TList genJetsListCut; // acceptance cut
906 TList recJetsList; // full acceptance
907 TList recJetsListCut; // acceptance cut
908
909
910 GetListOfJets(&recJetsList,aodRecJets,0);
911 GetListOfJets(&recJetsListCut,aodRecJets,1);
912
913 if(fBranchGen.Length()>0){
914 GetListOfJets(&genJetsList,aodGenJets,0);
915 GetListOfJets(&genJetsListCut,aodGenJets,1);
916 }
917
918 Double_t eventW = 1;
919 Double_t ptHard = 0;
920 Double_t nTrials = 1; // Trials for MC trigger
432abb2b 921 fh1AvgTrials->Fill("#sum{avg ntrials}",fAvgTrials);
a31b8a87 922
742ee86c 923 // Getting some global properties
924 fCentrality = GetCentrality();
aa3ba8d2 925 if(fCentrality<=0)fCentrality = 0;
742ee86c 926 fh1Centrality->Fill(fCentrality);
927
928
a31b8a87 929 if((fAnalysisType&kAnaMCESD)==kAnaMCESD){
930 // this is the part we only use when we have MC information
931 AliMCEvent* mcEvent = MCEvent();
932 // AliStack *pStack = 0;
933 if(!mcEvent){
934 Printf("%s:%d no mcEvent",(char*)__FILE__,__LINE__);
3b7ffecf 935 return;
936 }
a31b8a87 937 AliGenPythiaEventHeader* pythiaGenHeader = AliAnalysisHelperJetTasks::GetPythiaEventHeader(mcEvent);
938 if(pythiaGenHeader){
939 nTrials = pythiaGenHeader->Trials();
940 ptHard = pythiaGenHeader->GetPtHard();
941 int iProcessType = pythiaGenHeader->ProcessType();
942 // 11 f+f -> f+f
943 // 12 f+barf -> f+barf
944 // 13 f+barf -> g+g
945 // 28 f+g -> f+g
946 // 53 g+g -> f+barf
947 // 68 g+g -> g+g
948 if (fDebug > 10)Printf("%d iProcessType %d",__LINE__, iProcessType);
949 if(fDebug>20)AliAnalysisHelperJetTasks::PrintStack(mcEvent);
950 }
951 }// (fAnalysisType&kMCESD)==kMCESD)
952
953
954 // we simply fetch the tracks/mc particles as a list of AliVParticles
955
956 TList recParticles;
957 TList genParticles;
958
959 Int_t nT = GetListOfTracks(&recParticles,fTrackTypeRec);
4f6d1cc1 960
a31b8a87 961 if(fDebug>2)Printf("%s:%d Selected Rec tracks: %d %d",(char*)__FILE__,__LINE__,nT,recParticles.GetEntries());
962 nT = GetListOfTracks(&genParticles,fTrackTypeGen);
963 if(fDebug>2)Printf("%s:%d Selected Gen tracks: %d %d",(char*)__FILE__,__LINE__,nT,genParticles.GetEntries());
964
d3a3f33d 965 // CalculateReactionPlaneAngle(&recParticles);
2eded560 966 fRPAngle = 0;
967
968 if(fRPMethod==0)fRPAngle = aod->GetHeader()->GetEventplane();
969 else if(fRPMethod==1||fRPMethod==2){
970 fRPAngle = aod->GetHeader()->GetQTheta(fRPMethod);
971 }
972 fh1RP->Fill(fRPAngle);
973 fh2RPCentrality->Fill(fCentrality,fRPAngle);
a31b8a87 974 // Event control and counting ...
975 // MC
976 fh1PtHard->Fill(ptHard,eventW);
977 fh1PtHardNoW->Fill(ptHard,1);
978 fh1PtHardTrials->Fill(ptHard,nTrials);
979
980 // Real
784b1747 981 if(aod->GetPrimaryVertex()){// No vtx for pure MC
982 fh1ZVtx->Fill(aod->GetPrimaryVertex()->GetZ());
983 }
a31b8a87 984
a2522483 985
37eb26ea 986 Int_t recMult1 = recParticles.GetEntries();
987 Int_t genMult1 = genParticles.GetEntries();
988
742ee86c 989 Int_t recMult2 = MultFromJetRefs(aodBackRecJets);
990 Int_t genMult2 = MultFromJetRefs(aodBackGenJets);
37eb26ea 991
992 fh2MultRec->Fill(recMult1,recMult2);
993 fh2MultGen->Fill(genMult1,genMult2);
994 fMultRec = recMult1;
995 if(fMultRec<=0)fMultRec = recMult2;
996 fMultGen = genMult1;
997 if(fMultGen<=0)fMultGen = genMult2;
a2522483 998
e5ee1c19 999 Double_t var0[4] = {0,};
90a006c2 1000 var0[0] = fCentrality;
1001 var0[1] = fMultRec;
3d1dba29 1002 for(int it=0;it<fNTrigger;it++){
1003 if(fInputHandler->IsEventSelected()&fTriggerBit[it]){
1004 var0[2] = it;
e5ee1c19 1005 var0[3] = GetRho(recJetsList);
3d1dba29 1006 fhnEvent->Fill(var0);
1007 }
1008 }
a31b8a87 1009 // the loops for rec and gen should be indentical... pass it to a separate
1010 // function ...
1011 // Jet Loop
1012 // Track Loop
1013 // Jet Jet Loop
1014 // Jet Track loop
1015
1016 FillJetHistos(recJetsListCut,recParticles,kJetRec);
a2522483 1017 FillJetHistos(recJetsList,recParticles,kJetRecFull);
a31b8a87 1018 FillTrackHistos(recParticles,kJetRec);
1019
1020 FillJetHistos(genJetsListCut,genParticles,kJetGen);
a2522483 1021 FillJetHistos(genJetsList,genParticles,kJetGenFull);
a31b8a87 1022 FillTrackHistos(genParticles,kJetGen);
1023
1024 // Here follows the jet matching part
1025 // delegate to a separated method?
1026
b99018c6 1027 if(fDoMatching){
1028 FillMatchHistos(recJetsList,genJetsList);
1029 }
1030
a31b8a87 1031 if (fDebug > 10)Printf("%s:%d",(char*)__FILE__,__LINE__);
1032 PostData(1, fHistList);
90a006c2 1033}
a31b8a87 1034
1035void AliAnalysisTaskJetSpectrum2::FillJetHistos(TList &jetsList,TList &particlesList,Int_t iType){
1036
1037 if(iType>=kJetTypes){
1038 return;
3b7ffecf 1039 }
742ee86c 1040 if(!fFlagJetType[iType])return;
a31b8a87 1041
37eb26ea 1042 Int_t refMult = fMultRec;
1043 if(iType==kJetGen||iType==kJetGenFull){
1044 refMult = fMultGen;
1045 }
1046
a31b8a87 1047 Int_t nJets = jetsList.GetEntries();
1048 fh1NJets[iType]->Fill(nJets);
1049
1050 if(nJets<=0)return;
3b7ffecf 1051
a31b8a87 1052 Float_t pT1 = 0;
1053 Float_t phi0 = 0;
1054 Float_t phi1 = 0;
1055 Int_t ij0 = -1;
1056 Int_t ij1 = -1;
1057
db745d29 1058 Double_t var1[9] = {0,}; // jet number;p_{T,jet};cent;# tracks;RP;area;trigger;leadingTrackPt
17f4943e 1059 Double_t var1b[3] = {0,}; // p_{T,jet};cent;trigger;
7bee4353 1060
90a006c2 1061 var1[2] = fCentrality;
17f4943e 1062 var1b[1] = fCentrality;
90a006c2 1063 var1[3] = refMult;
1064
17f4943e 1065
1066
5662640e 1067 Double_t var2[9] = {0,}; // jet number;p_{T,jet};cent;#eta;#phi;area;trigger
90a006c2 1068 var2[2] = fCentrality;
a31b8a87 1069
1070 for(int ij = 0;ij < nJets;ij++){
1071 AliAODJet *jet = (AliAODJet*)jetsList.At(ij);
1072 Float_t ptJet = jet->Pt();
17f4943e 1073
1074
4a828dd2 1075 if(ptJet<0.150)ptJet = jet->GetPtSubtracted(0);
17f4943e 1076
1077 var1b[0] = ptJet;
d7396b04 1078 if(jet->Trigger()&fJetTriggerBestMask){
1079 fh1PtJetsInBest[iType]->Fill(ptJet);
17f4943e 1080 for(int it = 0;it <fNTrigger;it++){
1081 if(fInputHandler->IsEventSelected()&fTriggerBit[it]){
1082 var1b[2] = it;
1083 fhnJetPtBest[iType]->Fill(var1b);
1084 }
1085 }
d7396b04 1086 }
65c43de8 1087 if(jet->Trigger()&fJetTriggerExcludeMask){
1088 fh1PtJetsInRej[iType]->Fill(ptJet);
17f4943e 1089 for(int it = 0;it <fNTrigger;it++){
1090 if(fInputHandler->IsEventSelected()&fTriggerBit[it]){
1091 var1b[2] = it;
1092 fhnJetPtRej[iType]->Fill(var1b);
1093 }
1094 }
65c43de8 1095 continue;
1096 }
a31b8a87 1097 fh1PtJetsIn[iType]->Fill(ptJet);
a31b8a87 1098
b5d90baf 1099 // find the dijets assume sorting and acceptance cut...
a31b8a87 1100 if(ij==0){
1101 ij0 = ij;
a31b8a87 1102 phi0 = jet->Phi();
1103 if(phi0<0)phi0+=TMath::Pi()*2.;
1104 }
1105 else if(ptJet>pT1){
1106 // take only the backward hemisphere??
1107 phi1 = jet->Phi();
b5d90baf 1108 if(phi1<0)phi1+=TMath::Pi()*2.;
a31b8a87 1109 Float_t dPhi = phi1 - phi0;
1110 if(dPhi>TMath::Pi())dPhi = dPhi - 2.*TMath::Pi();
1111 if(dPhi<(-1.*TMath::Pi()))dPhi = dPhi + 2.*TMath::Pi();
1112 if(TMath::Abs(TMath::Pi()-dPhi)<fDeltaPhiWindow){
1113 ij1 = ij;
1114 pT1 = ptJet;
1115 }
1116 }
1117 // fill jet histos for kmax jets
adf3920d 1118 Float_t phiJet = jet->Phi();
1119 Float_t etaJet = jet->Eta();
1120 if(phiJet<0)phiJet+=TMath::Pi()*2.;
1121 fh1TmpRho->Reset();
1122 if(ij<kMaxJets)fh1PtIn[iType][ij]->Fill(ptJet);
1123
1124 fh1PtIn[iType][kMaxJets]->Fill(ptJet);
1125 // fill leading jets...
1126 AliVParticle *leadTrack = LeadingTrackFromJetRefs(jet);
4f6d1cc1 1127 //AliVParticle *leadTrack = LeadingTrackInCone(jet,&particlesList);
adf3920d 1128 Int_t phiBin = GetPhiBin(phiJet-fRPAngle);
1129 Double_t ptLead = jet->GetPtLeading(); //pT of leading jet
4f6d1cc1 1130 AliVParticle *tt = (AliVParticle*)particlesList.At(0);
1131 Float_t ttphi=tt->Phi();
1132 if(ttphi<0)ttphi+=TMath::Pi()*2.;
1133 Float_t ttpt=tt->Pt();
1134 Int_t phiBintt = GetPhiBin(ttphi-fRPAngle);
fedb5a47 1135 Double_t dphitrigjet=RelativePhi(ttphi,phiJet);
4f6d1cc1 1136 if(fTRP==1){
fedb5a47 1137 if(TMath::Abs(dphitrigjet)<TMath::Pi()-0.6) continue;
1138 var1[1] = ptJet;
1139 var1[4] = phiBintt;
1140 var1[5] = jet->EffectiveAreaCharged();
1141 var1[7] = ttpt;
1142 var1[8] = CheckAcceptance(phiJet,etaJet);
1143 }
4f6d1cc1 1144
fedb5a47 1145 if(fTRP==0){
1146 var1[1] = ptJet;
1147 var1[4] = phiBin;
1148 var1[5] = jet->EffectiveAreaCharged();
1149 var1[7] = ptLead;
1150 var1[8] = CheckAcceptance(phiJet,etaJet);
1151 }
1152
1153 //jet number;p_{T,jet};cent;#eta;#phi;area;trigger;acceptance bin;signed pt leading
adf3920d 1154 var2[1] = ptJet;
1155 var2[3] = etaJet;
1156 var2[4] = phiJet;
1157 var2[5] = jet->EffectiveAreaCharged();
1158 var2[7] = var1[8];
1159 var2[8] = (leadTrack?leadTrack->Charge()*ptLead:0);//pT of leading jet x charge
1160
1161 if(ij<kMaxJets){
1162 fh2LTrackPtJetPt[iType][ij]->Fill(jet->GetPtLeading(),ptJet);
1163 var1[0] = ij;
1164 var2[0] = ij;
3d1dba29 1165 for(int it = 0;it <fNTrigger;it++){
1166 if(fInputHandler->IsEventSelected()&fTriggerBit[it]){
1167 var1[6] = it;
adf3920d 1168 var2[6] = it;
3d1dba29 1169 fhnJetPt[iType]->Fill(var1);
1d478650 1170 fhnJetPtQA[iType]->Fill(var2);
3d1dba29 1171 }
1172 }
adf3920d 1173 }
1174 var1[0] = kMaxJets;// fill for all jets
1175 var2[0] = kMaxJets;// fill for all jets
1176 for(int it = 0;it <fNTrigger;it++){
1177 if(fInputHandler->IsEventSelected()&fTriggerBit[it]){
1178 var1[6] = it;
94c3d5e4 1179 var2[6] = it;
adf3920d 1180 fhnJetPt[iType]->Fill(var1);
1181 fhnJetPtQA[iType]->Fill(var2);
1182 }
1183 }
3d1dba29 1184
adf3920d 1185 fh2LTrackPtJetPt[iType][kMaxJets]->Fill(jet->GetPtLeading(),ptJet);
a31b8a87 1186
adf3920d 1187 if(particlesList.GetSize()&&ij<kMaxJets){
1188 // Particles... correlated with jets...
1189 for(int it = 0;it<particlesList.GetEntries();++it){
1190 AliVParticle *part = (AliVParticle*)particlesList.At(it);
1191 Float_t deltaR = jet->DeltaR(part);
1192 if(ptJet>0)fh1TmpRho->Fill(deltaR,part->Pt()/ptJet);
1193 }
1194 // fill the jet shapes
1195 }// if we have particles
a31b8a87 1196 }// Jet Loop
1197
1198
1199 // do something with dijets...
1200 if(ij0>=0&&ij1>0){
1201 AliAODJet *jet0 = (AliAODJet*)jetsList.At(ij0);
1202 Double_t ptJet0 = jet0->Pt();
1203 Double_t phiJet0 = jet0->Phi();
1204 if(phiJet0<0)phiJet0+=TMath::Pi()*2.;
1205
1206 AliAODJet *jet1 = (AliAODJet*)jetsList.At(ij1);
1207 Double_t ptJet1 = jet1->Pt();
1208 Double_t phiJet1 = jet1->Phi();
1209 if(phiJet1<0)phiJet1+=TMath::Pi()*2.;
1210
1211 Float_t deltaPhi = phiJet0 - phiJet1;
1212 if(deltaPhi>TMath::Pi())deltaPhi = deltaPhi - 2.*TMath::Pi();
1213 if(deltaPhi<(-1.*TMath::Pi()))deltaPhi = deltaPhi + 2.*TMath::Pi();
1214 deltaPhi = TMath::Abs(deltaPhi);
b5d90baf 1215 fh2DijetDeltaPhiPt[iType]->Fill(deltaPhi,ptJet1);
a31b8a87 1216
1217 Float_t asym = 9999;
1218 if((ptJet0+ptJet1)>0)asym = (ptJet0-ptJet1)/(ptJet0+ptJet1);
adf3920d 1219 fh2DijetAsymPt[iType]->Fill(asym,ptJet0);
1220 fh2DijetPt2vsPt1[iType]->Fill(ptJet0,ptJet1);
1221 fh2DijetDifvsSum[iType]->Fill(ptJet0+ptJet1,ptJet0-ptJet1);
1222 Float_t minv = 2.*(jet0->P()*jet1->P()-
1223 jet0->Px()*jet1->Px()-
1224 jet0->Py()*jet1->Py()-
1225 jet0->Pz()*jet1->Pz()); // assume mass == 0;
1226 if(minv<0)minv=0; // prevent numerical instabilities
1227 minv = TMath::Sqrt(minv);
1228 fh1DijetMinv[iType]->Fill(minv);
a31b8a87 1229 }
1230
1231
1232
1233 // count down the jets above thrueshold
1234 Int_t nOver = nJets;
1235 if(nOver>0){
1236 TIterator *jetIter = jetsList.MakeIterator();
1237 AliAODJet *tmpJet = (AliAODJet*)(jetIter->Next());
a31b8a87 1238 if(tmpJet){
e855f5c5 1239 Float_t pt = tmpJet->Pt();
a31b8a87 1240 for(int i = 1;i <= fh2NJetsPt[iType]->GetNbinsX();i++){
1241 Float_t ptCut = fh2NJetsPt[iType]->GetXaxis()->GetBinCenter(i);
1242 while(pt<ptCut&&tmpJet){
1243 nOver--;
1244 tmpJet = (AliAODJet*)(jetIter->Next());
1245 if(tmpJet){
1246 pt = tmpJet->Pt();
1247 }
1248 }
1249 if(nOver<=0)break;
1250 fh2NJetsPt[iType]->Fill(ptCut,nOver);
1251 }
1252 }
1253 delete jetIter;
1254 }
1255}
1256
1257void AliAnalysisTaskJetSpectrum2::FillTrackHistos(TList &particlesList,int iType){
1258
90a006c2 1259 if(fFlagJetType[iType]<=0)return;
7b822bc5 1260 Int_t refMult = fMultRec;
1261 if(iType==kJetGen||iType==kJetGenFull){
1262 refMult = fMultGen;
742ee86c 1263
7b822bc5 1264 }
1265
90a006c2 1266 //
5a7c81dc 1267 Double_t var3[7]; // track number;p_{T};cent;#tracks;RP;cahrge
90a006c2 1268 var3[2] = fCentrality;
1269 var3[3] = refMult;
5a7c81dc 1270 Double_t var4[6]; // track number;p_{T};cent;#eta;#phi;sign
90a006c2 1271 var4[2] = fCentrality;
a31b8a87 1272 Int_t nTrackOver = particlesList.GetSize();
1273 // do the same for tracks and jets
1274 if(nTrackOver>0){
1275 TIterator *trackIter = particlesList.MakeIterator();
1276 AliVParticle *tmpTrack = (AliVParticle*)(trackIter->Next());
1277 Float_t pt = tmpTrack->Pt();
1278 for(int i = 1;i <= fh2NTracksPt[iType]->GetNbinsX();i++){
1279 Float_t ptCut = fh2NTracksPt[iType]->GetXaxis()->GetBinCenter(i);
1280 while(pt<ptCut&&tmpTrack){
1281 nTrackOver--;
1282 tmpTrack = (AliVParticle*)(trackIter->Next());
1283 if(tmpTrack){
1284 pt = tmpTrack->Pt();
1285 }
1286 }
1287 if(nTrackOver<=0)break;
1288 fh2NTracksPt[iType]->Fill(ptCut,nTrackOver);
1289 }
1290
1291
1292 trackIter->Reset();
1293 AliVParticle *leading = (AliVParticle*)particlesList.At(0);
1294 Float_t sumPt = 0;
1295
1296 while((tmpTrack = (AliVParticle*)(trackIter->Next()))){
1297 Float_t tmpPt = tmpTrack->Pt();
1298 fh1PtTracksIn[iType]->Fill(tmpPt);
cb883243 1299 fh1PtTracksInLow[iType]->Fill(tmpPt);
742ee86c 1300
a31b8a87 1301 sumPt += tmpPt;
1302 Float_t tmpPhi = tmpTrack->Phi();
1303 if(tmpPhi<0)tmpPhi+=TMath::Pi()*2.;
2eded560 1304
1305
1306 Float_t phiRP = tmpPhi-fRPAngle;
1307 if(phiRP>TMath::Pi())phiRP -= TMath::Pi();
1308 if(phiRP<0)phiRP += TMath::Pi();
1309 if(phiRP<0)phiRP += TMath::Pi();
1310 const float allPhi = -1./180.*TMath::Pi();
1311
dae7dd67 1312 if(tmpPt<100){
1313 fp2MultRPPhiTrackPt[iType]->Fill(refMult,phiRP,tmpPt);
1314 fp2MultRPPhiTrackPt[iType]->Fill(refMult,allPhi,tmpPt);
1315
1316 fp2CentRPPhiTrackPt[iType]->Fill(fCentrality,phiRP,tmpPt);
1317 fp2CentRPPhiTrackPt[iType]->Fill(fCentrality,allPhi,tmpPt);
1318 }
742ee86c 1319 Int_t phiBin = GetPhiBin(tmpPhi-fRPAngle);
90a006c2 1320 var3[0] = 1;
1321 var3[1] = tmpPt;
1322 var3[4] = phiBin;
5a7c81dc 1323 var3[6] = tmpTrack->Charge();
90a006c2 1324
1325 var4[0] = 1;
1326 var4[1] = tmpPt;
1327 var4[3] = tmpTrack->Eta();
1328 var4[4] = tmpPhi;
5662640e 1329 var4[5] = tmpTrack->Charge();
90a006c2 1330
17f4943e 1331
1332 for(int it = 0;it <fNTrigger;it++){
1333 if(fInputHandler->IsEventSelected()&fTriggerBit[it]){
1334 var3[0] = 1;
1335 var3[5] = it;
1336 fhnTrackPt[iType]->Fill(var3);
1337 if(tmpTrack==leading){
1338 var3[0] = 0;
1339 fhnTrackPt[iType]->Fill(var3);
1340 }
1341 }
1342 }
1343
1344
1345
5a7c81dc 1346
1347 fhnTrackPtQA[iType]->Fill(var4);
1348 if(tmpTrack==leading){
1349 var4[0] = 0;
1617846c 1350 fhnTrackPtQA[iType]->Fill(var4);
5a7c81dc 1351 continue;
1617846c 1352 }
17f4943e 1353
1354
90a006c2 1355
1617846c 1356
1357
a31b8a87 1358 }
1359 fh1SumPtTrack[iType]->Fill(sumPt);
1360 delete trackIter;
1361 }
1362
1363}
1364
1365
1366void AliAnalysisTaskJetSpectrum2::FillMatchHistos(TList &recJetsList,TList &genJetsList){
1367
1368
1369 // Fill al the matching histos
1370 // It is important that the acceptances for the mathing are as large as possible
1371 // to avoid false matches on the edge of acceptance
1372 // therefore we add some extra matching jets as overhead
1373
1374 static TArrayI aGenIndex(recJetsList.GetEntries());
1375 if(aGenIndex.GetSize()<recJetsList.GetEntries())aGenIndex.Set(recJetsList.GetEntries());
3b7ffecf 1376
a31b8a87 1377 static TArrayI aRecIndex(genJetsList.GetEntries());
1378 if(aRecIndex.GetSize()<genJetsList.GetEntries())aRecIndex.Set(genJetsList.GetEntries());
3b7ffecf 1379
a31b8a87 1380 if(fDebug){
1381 Printf("New Gens List %d rec index Array %d",genJetsList.GetEntries(),aRecIndex.GetSize());
1382 Printf("New Rec List %d gen indey Array %d",recJetsList.GetEntries(),aGenIndex.GetSize());
1383 }
1384 AliAnalysisHelperJetTasks::GetClosestJets(&genJetsList,TMath::Min((Int_t)fNMatchJets,(Int_t)genJetsList.GetEntries()),
1385 &recJetsList,TMath::Min((Int_t)fNMatchJets,(Int_t)recJetsList.GetEntries()),
1386 aGenIndex,aRecIndex,fDebug);
1387
1388 if(fDebug){
1389 for(int i = 0;i< aGenIndex.GetSize();++i){
1390 if(aGenIndex[i]>=0)Printf("iGenFound: %d -> %d",i,aGenIndex[i]);
1391 }
1392 for(int i = 0;i< aRecIndex.GetSize();++i){
1393 if(aRecIndex[i]>=0)Printf("iRecFound: %d -> %d",i,aRecIndex[i]);
1394 }
1395 }
1396
152a0744 1397 Double_t container[8];
1398 Double_t containerRec[4];
1399 Double_t containerGen[4];
a31b8a87 1400
1401 // loop over generated jets
1402 // consider the
1403 for(int ig = 0;ig < genJetsList.GetEntries();++ig){
1404 AliAODJet *genJet = (AliAODJet*)genJetsList.At(ig);
1405 Double_t ptGen = genJet->Pt();
1406 Double_t phiGen = genJet->Phi();
1407 if(phiGen<0)phiGen+=TMath::Pi()*2.;
1408 Double_t etaGen = genJet->Eta();
152a0744 1409 containerGen[0] = ptGen;
1410 containerGen[1] = etaGen;
1411 containerGen[2] = phiGen;
1412 containerGen[3] = genJet->GetPtLeading();
1413
1414 fhnJetContainer->Fill(containerGen,kStep0); //all generated jets
1415
346e5a93 1416 if(JetSelected(genJet))
1417 fhnJetContainer->Fill(containerGen,kStep1); // all generated jets in eta window
1418
152a0744 1419 Int_t ir = aRecIndex[ig];
1420 if(ir>=0&&ir<recJetsList.GetEntries()){
346e5a93 1421 AliAODJet* recJet = (AliAODJet*)recJetsList.At(ir);
152a0744 1422
346e5a93 1423 fhnJetContainer->Fill(containerGen,kStep2); // all generated jets with reconstructed partner
152a0744 1424
346e5a93 1425 if(JetSelected(genJet)){
152a0744 1426 fhnJetContainer->Fill(containerGen,kStep3); // all generated jets in eta window with reconstructed partner
346e5a93 1427
152a0744 1428 if(JetSelected(recJet)) {
1429 fhnJetContainer->Fill(containerGen,kStep4); // all generated jets in eta window with reconstructed partner in eta window
346e5a93 1430
152a0744 1431 }
346e5a93 1432 containerGen[3] = recJet->GetPtLeading();
1433 fhnJetContainer->Fill(containerGen,kStep5); // all generated jets in eta window with reconstructed partner with leading track on reconstructed level
a31b8a87 1434 }
1435 }
1436 }// loop over generated jets used for matching...
1437
1438
1439
1440 // loop over reconstructed jets
1441 for(int ir = 0;ir < recJetsList.GetEntries();++ir){
1442 AliAODJet *recJet = (AliAODJet*)recJetsList.At(ir);
1443 Double_t etaRec = recJet->Eta();
1444 Double_t ptRec = recJet->Pt();
1445 Double_t phiRec = recJet->Phi();
1446 if(phiRec<0)phiRec+=TMath::Pi()*2.;
a31b8a87 1447
152a0744 1448 containerRec[0] = ptRec;
1449 containerRec[1] = etaRec;
1450 containerRec[2] = phiRec;
1451 containerRec[3] = recJet->GetPtLeading();
1452
a31b8a87 1453 container[0] = ptRec;
1454 container[1] = etaRec;
1455 container[2] = phiRec;
152a0744 1456 container[3] = recJet->GetPtLeading();
a31b8a87 1457
a31b8a87 1458 if (fDebug > 10)Printf("%s:%d",(char*)__FILE__,__LINE__);
1459
1460 if(JetSelected(recJet)){
152a0744 1461 fhnJetContainer->Fill(containerRec,kStep7); //all rec jets in eta window
a31b8a87 1462 // Fill Correlation
1463 Int_t ig = aGenIndex[ir];
1464 if(ig>=0 && ig<genJetsList.GetEntries()){
a31b8a87 1465 if (fDebug > 10)Printf("%s:%d ig = %d ir = %d",(char*)__FILE__,__LINE__,ig,ir);
1466 AliAODJet *genJet = (AliAODJet*)genJetsList.At(ig);
1467 Double_t ptGen = genJet->Pt();
1468 Double_t phiGen = genJet->Phi();
1469 if(phiGen<0)phiGen+=TMath::Pi()*2.;
1470 Double_t etaGen = genJet->Eta();
152a0744 1471
1472 containerGen[0] = ptGen;
1473 containerGen[1] = etaGen;
1474 containerGen[2] = phiGen;
1475 containerGen[3] = genJet->GetPtLeading();
1476
1477 container[4] = ptGen;
1478 container[5] = etaGen;
1479 container[6] = phiGen;
1480 container[7] = genJet->GetPtLeading();
1481
1482 fhnJetContainer->Fill(containerGen,kStep6); // all rec jets in eta window with generated partner
1483
1484 fhnCorrelation->Fill(container);
e39a4f22 1485
1486 Float_t delta = (ptRec-ptGen)/ptGen;
1487 fh2RelPtFGen->Fill(ptGen,delta);
1488 fh3RelPtFGenLeadTrkPt->Fill(ptGen,delta,recJet->GetPtLeading());
1489 fh2PtFGen->Fill(ptGen,ptRec);
1490
1491 fh2deltaPt1Pt2->Fill(ptRec-ptGen,ptRec,ptGen);
1492
a31b8a87 1493 }
a31b8a87 1494 }// loop over reconstructed jets
1495 }
1496 if (fDebug > 10)Printf("%s:%d",(char*)__FILE__,__LINE__);
1497}
1498
1499
3b7ffecf 1500void AliAnalysisTaskJetSpectrum2::MakeJetContainer(){
1501 //
1502 // Create the particle container for the correction framework manager and
1503 // link it
1504 //
152a0744 1505 const Int_t kNvar = 4 ; //number of variables on the grid:pt,eta, phi
d95fc15a 1506 const Double_t kPtmin = 0.0, kPtmax = 250.; // we do not want to have empty bins at the beginning...
152a0744 1507 const Double_t kEtamin = -1.0, kEtamax = 1.0;
3b7ffecf 1508 const Double_t kPhimin = 0., kPhimax = 2. * TMath::Pi();
152a0744 1509 const Double_t kPtLeadingTrackPtMin = 0., kPtLeadingTrackPtMax = 200.;
3b7ffecf 1510
1511 // can we neglect migration in eta and phi?
1512 // phi should be no problem since we cover full phi and are phi symmetric
1513 // eta migration is more difficult due to needed acceptance correction
1514 // in limited eta range
1515
3b7ffecf 1516 //arrays for the number of bins in each dimension
1517 Int_t iBin[kNvar];
d95fc15a 1518 iBin[0] = 125; //bins in pt
152a0744 1519 iBin[1] = 4; //bins in eta
1520 iBin[2] = 1; // bins in phi
1521 iBin[3] = 10; // bins in leading track Pt
3b7ffecf 1522
1523 //arrays for lower bounds :
1524 Double_t* binEdges[kNvar];
1525 for(Int_t ivar = 0; ivar < kNvar; ivar++)
1526 binEdges[ivar] = new Double_t[iBin[ivar] + 1];
1527
1528 //values for bin lower bounds
1529 // for(Int_t i=0; i<=iBin[0]; i++) binEdges[0][i]=(Double_t)TMath::Power(10,TMath::Log10(kPtmin) + (TMath::Log10(kPtmax)-TMath::Log10(kPtmin))/iBin[0]*(Double_t)i);
a923bd34 1530 for(Int_t i=0; i<=iBin[0]; i++) binEdges[0][i]=(Double_t)kPtmin + (kPtmax-kPtmin)/(Double_t)iBin[0]*(Double_t)i;
3b7ffecf 1531 for(Int_t i=0; i<=iBin[1]; i++) binEdges[1][i]=(Double_t)kEtamin + (kEtamax-kEtamin)/iBin[1]*(Double_t)i;
1532 for(Int_t i=0; i<=iBin[2]; i++) binEdges[2][i]=(Double_t)kPhimin + (kPhimax-kPhimin)/iBin[2]*(Double_t)i;
152a0744 1533 binEdges[3][0]= kPtLeadingTrackPtMin;
1534 binEdges[3][1]= 1.;
1535 binEdges[3][2]= 2.;
1536 binEdges[3][3]= 3.;
1537 binEdges[3][4]= 4.;
1538 binEdges[3][5]= 5.;
1539 binEdges[3][6]= 6.;
1540 binEdges[3][7]= 8.;
1541 binEdges[3][8]= 10.;
1542 binEdges[3][9]= 12.;
1543 binEdges[3][10]= kPtLeadingTrackPtMax;
1544
1545 fhnJetContainer = new AliCFContainer(Form("fhnJetContainer"),Form("AliCFContainer jet info"),kMaxStep,kNvar,iBin);
b99018c6 1546 for (int k=0; k<kNvar; k++) {
1547 fhnJetContainer->SetBinLimits(k,binEdges[k]);
3b7ffecf 1548 }
b99018c6 1549
3b7ffecf 1550 //create correlation matrix for unfolding
1551 Int_t thnDim[2*kNvar];
1552 for (int k=0; k<kNvar; k++) {
1553 //first half : reconstructed
1554 //second half : MC
1555 thnDim[k] = iBin[k];
1556 thnDim[k+kNvar] = iBin[k];
1557 }
1558
152a0744 1559 fhnCorrelation = new THnSparseF("fhnCorrelation","THnSparseF with correlations",2*kNvar,thnDim);
3b7ffecf 1560 for (int k=0; k<kNvar; k++) {
152a0744 1561 fhnCorrelation->SetBinEdges(k,binEdges[k]);
1562 fhnCorrelation->SetBinEdges(k+kNvar,binEdges[k]);
3b7ffecf 1563 }
742ee86c 1564
1565 for(Int_t ivar = 0; ivar < kNvar; ivar++)
1566 delete [] binEdges[ivar];
1567
1568
3b7ffecf 1569}
1570
1571void AliAnalysisTaskJetSpectrum2::Terminate(Option_t */*option*/)
1572{
adf3920d 1573 // Terminate analysis
1574 //
1575 if (fDebug > 1) printf("AnalysisJetSpectrum2: Terminate() \n");
3b7ffecf 1576}
1577
1578
1579Int_t AliAnalysisTaskJetSpectrum2::GetListOfTracks(TList *list,Int_t type){
1580
565584e8 1581 if(fDebug>2)Printf("%s:%d Selecting tracks with %d",(char*)__FILE__,__LINE__,type);
3b7ffecf 1582
1583 Int_t iCount = 0;
565584e8 1584 if(type==kTrackAOD){
3b7ffecf 1585 AliAODEvent *aod = 0;
565584e8 1586 if(fUseAODTrackInput)aod = dynamic_cast<AliAODEvent*>(InputEvent());
1587 else aod = AODEvent();
3b7ffecf 1588 if(!aod){
1589 return iCount;
1590 }
1591 for(int it = 0;it < aod->GetNumberOfTracks();++it){
1592 AliAODTrack *tr = aod->GetTrack(it);
1593 if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask)))continue;
a31b8a87 1594 if(TMath::Abs(tr->Eta())>fTrackRecEtaWindow)continue;
1595 if(tr->Pt()<fMinTrackPt)continue;
38ecb6a5 1596 if(fDebug>0){
1597 if(tr->Pt()>20){
1598 Printf("High pT track found in Event %d with p_T, %E",(int)Entry(),tr->Pt());
fd5db301 1599 Printf("%s read event, %d",fInputHandler->GetTree()->GetCurrentFile()->GetName(),(Int_t)fInputHandler->GetTree()->GetReadEntry());
38ecb6a5 1600 tr->Print();
1601 // tr->Dump();
1602 AliESDEvent *fESD = dynamic_cast<AliESDEvent*> (InputEvent());
1603 if(fESD){
a2522483 1604 AliESDtrack *esdTr = (AliESDtrack*)fESD->GetTrack(TMath::Abs(tr->GetID()+1));
38ecb6a5 1605 esdTr->Print("");
bb3d13aa 1606 // esdTr->Dump();
38ecb6a5 1607 }
1608 }
1609 }
3b7ffecf 1610 list->Add(tr);
1611 iCount++;
1612 }
1613 }
1614 else if (type == kTrackKineAll||type == kTrackKineCharged){
1615 AliMCEvent* mcEvent = MCEvent();
1616 if(!mcEvent)return iCount;
1617 // we want to have alivpartilces so use get track
1618 for(int it = 0;it < mcEvent->GetNumberOfTracks();++it){
1619 if(!mcEvent->IsPhysicalPrimary(it))continue;
1620 AliMCParticle* part = (AliMCParticle*)mcEvent->GetTrack(it);
a31b8a87 1621 if(part->Pt()<fMinTrackPt)continue;
3b7ffecf 1622 if(type == kTrackKineAll){
1623 list->Add(part);
1624 iCount++;
1625 }
1626 else if(type == kTrackKineCharged){
1627 if(part->Particle()->GetPDG()->Charge()==0)continue;
1628 list->Add(part);
1629 iCount++;
1630 }
1631 }
1632 }
565584e8 1633 else if (type == kTrackAODMCCharged || type == kTrackAODMCAll || type == kTrackAODMCChargedAcceptance) {
1634 AliAODEvent *aod = 0;
5301f754 1635 if(fUseAODMCInput)aod = dynamic_cast<AliAODEvent*>(InputEvent());
565584e8 1636 else aod = AODEvent();
5010a3f7 1637 if(!aod)return iCount;
3b7ffecf 1638 TClonesArray *tca = dynamic_cast<TClonesArray*>(aod->FindListObject(AliAODMCParticle::StdBranchName()));
1639 if(!tca)return iCount;
1640 for(int it = 0;it < tca->GetEntriesFast();++it){
1641 AliAODMCParticle *part = dynamic_cast<AliAODMCParticle*>(tca->At(it));
e855f5c5 1642 if(!part)continue;
a31b8a87 1643 if(part->Pt()<fMinTrackPt)continue;
3b7ffecf 1644 if(!part->IsPhysicalPrimary())continue;
c4631cdb 1645 if(type == kTrackAODMCAll){
3b7ffecf 1646 list->Add(part);
1647 iCount++;
1648 }
565584e8 1649 else if (type == kTrackAODMCCharged || type == kTrackAODMCChargedAcceptance ){
3b7ffecf 1650 if(part->Charge()==0)continue;
565584e8 1651 if(kTrackAODMCCharged){
1652 list->Add(part);
1653 }
1654 else {
a31b8a87 1655 if(TMath::Abs(part->Eta())>fTrackRecEtaWindow)continue;
565584e8 1656 list->Add(part);
1657 }
3b7ffecf 1658 iCount++;
1659 }
1660 }
1661 }// AODMCparticle
cc0649e4 1662 list->Sort();
c4631cdb 1663 return iCount;
3b7ffecf 1664
1665}
a31b8a87 1666
1667
742ee86c 1668Float_t AliAnalysisTaskJetSpectrum2::GetCentrality(){
adf3920d 1669 AliAODEvent *aod = 0;
1670 if(fUseAODTrackInput)aod = dynamic_cast<AliAODEvent*>(InputEvent());
1671 else aod = AODEvent();
1672 if(!aod){
1673 return 101;
1674 }
1675 return aod->GetHeader()->GetCentrality();
742ee86c 1676}
1677
1678
a31b8a87 1679
1680Bool_t AliAnalysisTaskJetSpectrum2::JetSelected(AliAODJet *jet){
1681 Bool_t selected = false;
1682
1683 if(!jet)return selected;
1684
1685 if(fabs(jet->Eta())<fJetRecEtaWindow&&jet->Pt()>fMinJetPt){
1686 selected = kTRUE;
1687 }
1688 return selected;
1689
1690}
1691
1692Int_t AliAnalysisTaskJetSpectrum2::GetListOfJets(TList *list,TClonesArray* jarray,Int_t type){
1693
1694 if(fDebug>2)Printf("%s:%d Selecting jets with cuts %d",(char*)__FILE__,__LINE__,type);
1695 Int_t iCount = 0;
1696
1697 if(!jarray){
1698 Printf("%s:%d no Jet array",(char*)__FILE__,__LINE__);
1699 return iCount;
1700 }
1701
1702
1703 for(int ij=0;ij<jarray->GetEntries();ij++){
1704 AliAODJet* jet = (AliAODJet*)jarray->At(ij);
1705 if(!jet)continue;
1706 if(type==0){
1707 // No cut at all, main purpose here is sorting
1708 list->Add(jet);
1709 iCount++;
1710 }
1711 else if(type == 1){
1712 // eta cut
1713 if(JetSelected(jet)){
1714 list->Add(jet);
1715 iCount++;
1716 }
1717 }
1718 }
1719
1720 list->Sort();
1721 return iCount;
1722
1723}
1724
1725
37eb26ea 1726Int_t AliAnalysisTaskJetSpectrum2::MultFromJetRefs(TClonesArray* jets){
1727 if(!jets)return 0;
1728
1729 Int_t refMult = 0;
1730 for(int ij = 0;ij < jets->GetEntries();++ij){
1731 AliAODJet* jet = (AliAODJet*)jets->At(ij);
1732 if(!jet)continue;
1733 TRefArray *refs = jet->GetRefTracks();
1734 if(!refs)continue;
1735 refMult += refs->GetEntries();
1736 }
1737 return refMult;
1738
1739}
d7117cbd 1740
1741
1742AliVParticle *AliAnalysisTaskJetSpectrum2::LeadingTrackFromJetRefs(AliAODJet* jet){
1743 if(!jet)return 0;
1744 TRefArray *refs = jet->GetRefTracks();
1745 if(!refs) return 0;
1746 AliVParticle *leading = 0;
1747 Float_t fMaxPt = 0;
1748 for(int i = 0;i<refs->GetEntries();i++){
1749 AliVParticle *tmp = dynamic_cast<AliVParticle*>(refs->At(i));
1750 if(!tmp)continue;
1751 if(tmp->Pt()>fMaxPt){
1752 leading = tmp;
1753 fMaxPt = tmp->Pt();
1754 }
1755 }
1756 return leading;
1757}
1758
1759
1760AliVParticle *AliAnalysisTaskJetSpectrum2::LeadingTrackInCone(AliAODJet* jet,TList *list,Float_t r){
1761 if(!jet)return 0;
1762 if(!list) return 0;
1763 AliVParticle *leading = 0;
1764 Float_t fMaxPt = 0;
1765 for(int i = 0;i<list->GetEntries();i++){
1766 AliVParticle *tmp = (AliVParticle*)(list->At(i));
1767 if(!tmp)continue;
1768 if(jet->DeltaR(tmp)>r)continue;
1769 if(tmp->Pt()>fMaxPt){
1770 leading = tmp;
1771 fMaxPt = tmp->Pt();
1772 }
1773 }
1774 return leading;
1775}
742ee86c 1776
4f6d1cc1 1777
1778Double_t AliAnalysisTaskJetSpectrum2::RelativePhi(Double_t mphi,Double_t vphi){
1779
1780 if (vphi < -1*TMath::Pi()) vphi += (2*TMath::Pi());
1781 else if (vphi > TMath::Pi()) vphi -= (2*TMath::Pi());
1782 if (mphi < -1*TMath::Pi()) mphi += (2*TMath::Pi());
1783 else if (mphi > TMath::Pi()) mphi -= (2*TMath::Pi());
1784 double dphi = mphi-vphi;
1785 if (dphi < -1*TMath::Pi()) dphi += (2*TMath::Pi());
1786 else if (dphi > TMath::Pi()) dphi -= (2*TMath::Pi());
1787 return dphi;//dphi in [-Pi, Pi]
1788}
1789
1790
1791
1792
1793
1794
742ee86c 1795Int_t AliAnalysisTaskJetSpectrum2::GetPhiBin(Double_t phi)
1796{
adf3920d 1797 Int_t phibin=-1;
1798 if(!(TMath::Abs(phi)<=2*TMath::Pi())){AliError("phi w.r.t. RP out of defined range");return -1;}
1799 Double_t phiwrtrp=TMath::ACos(TMath::Abs(TMath::Cos(phi)));
1800 phibin=Int_t(fNRPBins*phiwrtrp/(0.5*TMath::Pi()));
1801 if(phibin<0||phibin>=fNRPBins){AliError("Phi Bin not defined");}
1802 return phibin;
742ee86c 1803}
1804
3d1dba29 1805void AliAnalysisTaskJetSpectrum2::SetNTrigger(Int_t n){
db745d29 1806 //
1807 // set number of trigger bins
1808 //
3d1dba29 1809 if(n>0){
adf3920d 1810 fNTrigger = n;
1811 delete [] fTriggerName;
1812 fTriggerName = new TString [fNTrigger];
1813 delete [] fTriggerBit;fTriggerBit = 0;
1814 fTriggerBit = new UInt_t [fNTrigger];
3d1dba29 1815 }
1816 else{
1817 fNTrigger = 0;
1818 }
1819}
1820
db745d29 1821
3d1dba29 1822void AliAnalysisTaskJetSpectrum2::SetTrigger(Int_t i,UInt_t it,const char* c){
db745d29 1823 //
1824 // set trigger bin
1825 //
3d1dba29 1826 if(i<fNTrigger){
1827 fTriggerBit[i] = it;
1828 fTriggerName[i] = c;
1829 }
1830}
1831
db745d29 1832
1833void AliAnalysisTaskJetSpectrum2::SetNAcceptance(Int_t n){
1834 //
1835 // Set number of acceptance bins
1836 //
1837
1838 if(n>0){
adf3920d 1839 fNAcceptance = n;
1840 delete [] fAcceptancePhiMin;
1841 delete [] fAcceptancePhiMax;
1842 delete [] fAcceptanceEtaMin;
1843 delete [] fAcceptanceEtaMax;
1844
1845 fAcceptancePhiMin = new Float_t[fNAcceptance];
1846 fAcceptancePhiMax = new Float_t[fNAcceptance];
1847 fAcceptanceEtaMin = new Float_t[fNAcceptance];
1848 fAcceptanceEtaMax = new Float_t[fNAcceptance];
db745d29 1849 }
1850 else{
1851 fNTrigger = 0;
1852 }
1853}
1854
1855void AliAnalysisTaskJetSpectrum2::SetAcceptance(Int_t i,Float_t phiMin,Float_t phiMax,Float_t etaMin,Float_t etaMax){
1856 //
1857 // Set acceptance bins
1858 //
1859
1860 if(i<fNAcceptance){
24cc3482 1861 Printf("%s:%d Setting acceptance %d phi %3.3f-%3.3f eta %3.3f-%3.3f",(char*)__FILE__,__LINE__,i,phiMin,phiMax,etaMin,etaMax);
1862
1863 fAcceptancePhiMin[i] = phiMin;
1864 fAcceptancePhiMax[i] = phiMax;
1865 fAcceptanceEtaMin[i] = etaMin;
1866 fAcceptanceEtaMax[i] = etaMax;
db745d29 1867 }
1868 else{
1869 fNTrigger = 0;
1870 }
1871}
1872
1873
1874Int_t AliAnalysisTaskJetSpectrum2::CheckAcceptance(Float_t phi,Float_t eta){
1875
1876 //
1877 // return aceptnace bin do no use ovelapping bins
1878 //
1879
24cc3482 1880 for(int i = 0;i<fNAcceptance;i++){
db745d29 1881 if(phi<fAcceptancePhiMin[i])continue;
1882 if(phi>fAcceptancePhiMax[i])continue;
1883 if(eta<fAcceptanceEtaMin[i])continue;
1884 if(eta>fAcceptanceEtaMax[i])continue;
1885 return i;
1886 }
1887 // catch the rest
1888 return fNAcceptance;
1889}
1890
e5ee1c19 1891Float_t AliAnalysisTaskJetSpectrum2::GetRho(TList &list){
db745d29 1892
e5ee1c19 1893 // invert the correction
1894 AliAODJet *jet = (AliAODJet*)list.At(0); // highest pt jet
1895 if(!jet)return -1;
1896 if(jet->EffectiveAreaCharged()<=0)return -1;
1897 Float_t rho = jet->ChargedBgEnergy()/jet->EffectiveAreaCharged();
1898 return rho;
1899
1900}
db745d29 1901
1902
1903
3d1dba29 1904AliAnalysisTaskJetSpectrum2::~AliAnalysisTaskJetSpectrum2(){
1905 //
1906 delete [] fTriggerBit;
1907 delete [] fTriggerName;
db745d29 1908
1909 delete [] fAcceptancePhiMin;
1910 delete [] fAcceptancePhiMax;
1911 delete [] fAcceptanceEtaMin;
1912 delete [] fAcceptanceEtaMax;
1913
1914
1915
3d1dba29 1916}