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