]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGJE/AliAnalysisTaskJetResponseV2.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
allow for jet by jet comparisio
[u/mrichter/AliRoot.git] / PWGJE / AliAnalysisTaskJetResponseV2.cxx
CommitLineData
7a88ef30 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16//
d5ec3aef 17// task compares jets in two or three branches,
7a88ef30 18// written for analysis of jet embedding in HI events
19//
20
21
3c6a60f7 22#include "TChain.h"
23#include "TTree.h"
24#include "TMath.h"
25#include "TH1F.h"
8c483a6b 26#include "TH1I.h"
3c6a60f7 27#include "TH2F.h"
28#include "TH3F.h"
29#include "THnSparse.h"
30#include "TCanvas.h"
31
32#include "AliLog.h"
33
34#include "AliAnalysisTask.h"
35#include "AliAnalysisManager.h"
36
37#include "AliVEvent.h"
38#include "AliESDEvent.h"
39#include "AliESDInputHandler.h"
40#include "AliCentrality.h"
41#include "AliAnalysisHelperJetTasks.h"
42#include "AliInputEventHandler.h"
d314de7b 43#include "AliAODJetEventBackground.h"
31b9d515 44#include "AliAnalysisTaskFastEmbedding.h"
3c6a60f7 45
46#include "AliAODEvent.h"
47#include "AliAODJet.h"
48
49#include "AliAnalysisTaskJetResponseV2.h"
50
c64cb1f6 51using std::cout;
52using std::endl;
53
3c6a60f7 54ClassImp(AliAnalysisTaskJetResponseV2)
55
56AliAnalysisTaskJetResponseV2::AliAnalysisTaskJetResponseV2() :
31b9d515 57AliAnalysisTaskSE(),
d5ec3aef 58 fESD(0x0),
59 fAOD(0x0),
60 fBackgroundBranch(""),
61 fIsPbPb(kTRUE),
62 fOfflineTrgMask(AliVEvent::kAny),
63 fMinContribVtx(1),
64 fVtxZMin(-8.),
65 fVtxZMax(8.),
66 fEvtClassMin(0),
67 fEvtClassMax(4),
68 fCentMin(0.),
69 fCentMax(100.),
70 fNInputTracksMin(0),
71 fNInputTracksMax(-1),
72 fJetEtaMin(-.5),
73 fJetEtaMax(.5),
74 fJetPtMin(20.),
75 fJetTriggerExcludeMask(AliAODJet::kHighTrackPtTriggered),
76 fJetPtFractionMin(0.5),
77 fNMatchJets(4),
78 fMatchMaxDist(0.8),
79 fKeepJets(kFALSE),
80 fkNbranches(2),
81 fkEvtClasses(12),
82 fOutputList(0x0),
83 fbEvent(kTRUE),
84 fbJetsMismatch1(kTRUE),
85 fbJetsMismatch2(kTRUE),
86 fbJetsRp(kTRUE),
87 fbJetsDeltaPt(kTRUE),
88 fbJetsEta(kTRUE),
89 fbJetsPhi(kTRUE),
90 fbJetsArea(kTRUE),
91 fbJets3Branches(kFALSE),
92 fbJetsBeforeCut1(kTRUE),
93 fbJetsBeforeCut2(kTRUE),
94 fHistEvtSelection(0x0),
95 fHistJetSelection(0x0),
96 fh2JetSelection(0x0),
97 fhnEvent(0x0),
98 fhnJetsMismatch1(0x0),
99 fhnJetsMismatch2(0x0),
100 fhnJetsRp(0x0),
101 fhnJetsDeltaPt(0x0),
102 fhnJetsEta(0x0),
103 fhnJetsPhi(0x0),
104 fhnJetsArea(0x0),
105 fhnJets3Branches(0x0),
106 fhnJetsBeforeCut1(0x0),
107 fhnJetsBeforeCut2(0x0)
3c6a60f7 108{
d5ec3aef 109 // default Constructor
3c6a60f7 110
d5ec3aef 111 fJetBranchName[0] = "";
112 fJetBranchName[1] = "";
113 fJetBranchName[2] = "";
3c6a60f7 114
d5ec3aef 115 fListJets[0] = new TList;
116 fListJets[1] = new TList;
117 fListJets[2] = new TList;
3c6a60f7 118}
119
120AliAnalysisTaskJetResponseV2::AliAnalysisTaskJetResponseV2(const char *name) :
d5ec3aef 121 AliAnalysisTaskSE(name),
122 fESD(0x0),
123 fAOD(0x0),
124 fBackgroundBranch(""),
125 fIsPbPb(kTRUE),
126 fOfflineTrgMask(AliVEvent::kAny),
127 fMinContribVtx(1),
128 fVtxZMin(-8.),
129 fVtxZMax(8.),
130 fEvtClassMin(0),
131 fEvtClassMax(4),
132 fCentMin(0.),
133 fCentMax(100.),
134 fNInputTracksMin(0),
135 fNInputTracksMax(-1),
136 fJetEtaMin(-.5),
137 fJetEtaMax(.5),
138 fJetPtMin(20.),
139 fJetTriggerExcludeMask(AliAODJet::kHighTrackPtTriggered),
140 fJetPtFractionMin(0.5),
141 fNMatchJets(4),
142 fMatchMaxDist(0.8),
143 fKeepJets(kFALSE),
144 fkNbranches(2),
145 fkEvtClasses(12),
146 fOutputList(0x0),
147 fbEvent(kTRUE),
148 fbJetsMismatch1(kTRUE),
149 fbJetsMismatch2(kTRUE),
150 fbJetsRp(kTRUE),
151 fbJetsDeltaPt(kTRUE),
152 fbJetsEta(kTRUE),
153 fbJetsPhi(kTRUE),
154 fbJetsArea(kTRUE),
155 fbJets3Branches(kFALSE),
156 fbJetsBeforeCut1(kTRUE),
157 fbJetsBeforeCut2(kTRUE),
158 fHistEvtSelection(0x0),
159 fHistJetSelection(0x0),
160 fh2JetSelection(0x0),
161 fhnEvent(0x0),
162 fhnJetsMismatch1(0x0),
163 fhnJetsMismatch2(0x0),
164 fhnJetsRp(0x0),
165 fhnJetsDeltaPt(0x0),
166 fhnJetsEta(0x0),
167 fhnJetsPhi(0x0),
168 fhnJetsArea(0x0),
169 fhnJets3Branches(0x0),
170 fhnJetsBeforeCut1(0x0),
171 fhnJetsBeforeCut2(0x0)
3c6a60f7 172{
d5ec3aef 173 // Constructor
3c6a60f7 174
d5ec3aef 175 fJetBranchName[0] = "";
176 fJetBranchName[1] = "";
177 fJetBranchName[2] = "";
3c6a60f7 178
d5ec3aef 179 fListJets[0] = new TList;
180 fListJets[1] = new TList;
181 fListJets[2] = new TList;
3c6a60f7 182
d5ec3aef 183 DefineOutput(1, TList::Class());
3c6a60f7 184}
185
186AliAnalysisTaskJetResponseV2::~AliAnalysisTaskJetResponseV2()
187{
d5ec3aef 188 delete fListJets[0];
189 delete fListJets[1];
190 delete fListJets[2];
3c6a60f7 191}
192
d5ec3aef 193void AliAnalysisTaskJetResponseV2::SetBranchNames(const TString &branch1, const TString &branch2, const TString &branch3)
3c6a60f7 194{
d5ec3aef 195 fJetBranchName[0] = branch1;
196 fJetBranchName[1] = branch2;
197 fJetBranchName[2] = branch3;
198
199 if(strlen(fJetBranchName[2].Data()) )
200 fkNbranches = 3;
3c6a60f7 201}
202
203void AliAnalysisTaskJetResponseV2::Init()
204{
205
d5ec3aef 206 // check for jet branches
207 if(!strlen(fJetBranchName[0].Data()) || !strlen(fJetBranchName[1].Data())){
208 AliError("Jet branch name not set.");
209 }
3c6a60f7 210
211}
212
213void AliAnalysisTaskJetResponseV2::UserCreateOutputObjects()
214{
d5ec3aef 215 // Create histograms
216 // Called once
217 OpenFile(1);
218 if(!fOutputList) fOutputList = new TList;
219 fOutputList->SetOwner(kTRUE);
220
221 Bool_t oldStatus = TH1::AddDirectoryStatus();
222 TH1::AddDirectory(kFALSE);
223
224
225 fHistEvtSelection = new TH1I("fHistEvtSelection", "event selection", 6, -0.5, 5.5);
226 fHistEvtSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
227 fHistEvtSelection->GetXaxis()->SetBinLabel(2,"events IN");
228 fHistEvtSelection->GetXaxis()->SetBinLabel(3,"event selection (rejected)");
229 fHistEvtSelection->GetXaxis()->SetBinLabel(4,"vertex cut (rejected)");
230 fHistEvtSelection->GetXaxis()->SetBinLabel(5,"centrality (rejected)");
231 fHistEvtSelection->GetXaxis()->SetBinLabel(6,"multiplicity (rejected)");
232
233 fHistJetSelection = new TH1I("fHistJetSelection", "jet selection", 8, -0.5, 7.5);
234 fHistJetSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
235 fHistJetSelection->GetXaxis()->SetBinLabel(2,"probes IN");
236 fHistJetSelection->GetXaxis()->SetBinLabel(3,"no matching jet");
237 fHistJetSelection->GetXaxis()->SetBinLabel(4,"not in list");
238 fHistJetSelection->GetXaxis()->SetBinLabel(5,"fraction cut");
239 fHistJetSelection->GetXaxis()->SetBinLabel(6,"acceptance cut");
240 fHistJetSelection->GetXaxis()->SetBinLabel(7,"p_{T} cut");
241 fHistJetSelection->GetXaxis()->SetBinLabel(8,"trigger exclude mask");
242
243 fh2JetSelection = new TH2F("fh2JetSelection", "jet selection", 8, -0.5, 7.5,100,0.,200.);
244 fh2JetSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
245 fh2JetSelection->GetXaxis()->SetBinLabel(2,"probes IN");
246 fh2JetSelection->GetXaxis()->SetBinLabel(3,"no matching jet");
247 fh2JetSelection->GetXaxis()->SetBinLabel(4,"not in list");
248 fh2JetSelection->GetXaxis()->SetBinLabel(5,"fraction cut");
249 fh2JetSelection->GetXaxis()->SetBinLabel(6,"acceptance cut");
250 fh2JetSelection->GetXaxis()->SetBinLabel(7,"p_{T} cut");
251 fh2JetSelection->GetXaxis()->SetBinLabel(8,"trigger exclude mask");
252
253
254 UInt_t entries = 0; // bit coded, see GetDimParams() below
255 UInt_t opt = 0; // bit coded, default (0) or high resolution (1)
256
257 if(fbEvent){
258 entries = 1<<0 | 1<<1 | 1<<2 | 1<<26; // cent : nInpTrks : rp psi : pT hard bin
259 opt = 1<<0 | 1<<1; // centrality and nInpTrks in high resolution
260 fhnEvent = NewTHnSparseF("fhnEvent", entries, opt);
261 }
31b9d515 262
d5ec3aef 263 if(fbJetsMismatch1){ // real mismatch (no related rec jet found)
264 // cent : nInpTrks : rp bins : probe pt : probe eta : probe phi : probe area : pT hard bin
265 entries = 1<<0 | 1<<1 | 1<<3 | 1<<6 | 1<<8 | 1<<10 | 1<<12 | 1<<26;
266 opt = 1<<6 | 1<<8 | 1<<10;
267 fhnJetsMismatch1 = NewTHnSparseF("fhnJetsMismatch1", entries, opt);
268 }
269
270 if(fbJetsMismatch2){ // acceptance + fraction cut
271 // cent : nInpTrks : rp bins : jetPt(2x) : jetEta(2x) : deltaEta : deltaR : fraction : pT hard bin
272 entries = 1<<0 | 1<<1 | 1<<3 | 1<<6 | 1<<7 | 1<<8 | 1<<9 | 1<<15 | 1<<17 | 1<<19 | 1<<26;
273 opt = 1<<6 | 1<<7 | 1<<8 | 1<<9;
274 fhnJetsMismatch2 = NewTHnSparseF("fhnJetsMismatch2", entries, opt);
275
276 }
d314de7b 277
d314de7b 278
d5ec3aef 279 if(fbJetsRp){
280 // cent : nInpTrks : rp bins : rp wrt jet(2x) : probe pT : probe area: deltaPt : rp phi : rho : correction for RP : local rho
281 /*
282 entries = 1<<0 | 1<<1 | 1<<3 | 1<<4 | 1<<5 | 1<<6 | 1<<12 | 1<<14 | 1<<2 | 1<<22 | 1<<23 | 1<<24 | 1<<25;
283 opt = 1<<4 | 1<<5;*/
284 // cent : nInpTrks : rp bins : rp wrt jet(2x) : probe pT : deltaPt : rp phi : rho : correction for RP | pT hard bin
285 entries = 1<<0 | 1<<1 | 1<<3 | 1<<4 | 1<<5 | 1<<6 | 1<<14 | 1<<2 | 1<<22 | 1<<23 | 1<<26;
286 opt = 1<<4 | 1<<5;
287 fhnJetsRp = NewTHnSparseF("fhnJetsRp", entries, opt);
288 }
289
290 // cent : nInpTrks : rp bins: deltaPt : jetPt(2x) : deltaArea : pT hard bin (hr delta pt)
291 if(fbJetsDeltaPt){
292 entries = 1<<0 | 1<<1 | 1<<3 | 1<<14 | 1<<6 | 1<<7 | 1<<18 | 1<<26;
293 opt = 1<<1 | 1<<14 | 1<<6 | 1<<7 | 1<<18;
294 fhnJetsDeltaPt = NewTHnSparseF("fhnJetsDeltaPt", entries, opt);
295 }
296
297 // cent : nInpTrks : rp bins : deltaPt : jetPt(2x) : deltaR : deltaEta : jetEta(2x) : pT hard bin (hr for eta)
298 if(fbJetsEta){
299 entries = 1<<0 | 1<<1 | 1<<3 | 1<<14 | 1<<6 | 1<<7 | 1<<17 | 1<<15 | 1<<8 | 1<<9 | 1<<26;
300 opt = 1<<15 | 1<<8 | 1<<9;
301 fhnJetsEta = NewTHnSparseF("fhnJetsEta", entries, opt);
302 }
303
304 // cent : nInpTrks : rp bins : jetPt(2x) : jetPhi(2x) : deltaPt : deltaPhi : pT hard bin
305 if(fbJetsPhi){
306 entries = 1<<0 | 1<<1 | 1<<3 | 1<<6 | 1<<7 | 1<<10 | 1<<11 | 1<<14 | 1<<16 | 1<<26;
307 opt = 1<<10 | 1<<11;
308 fhnJetsPhi = NewTHnSparseF("fhnJetsPhi", entries, opt);
309 }
310
311 // cent : nInpTrks : rp bins : deltaArea : jetArea(2x) : deltaR : fraction : distance next rec jet : pT next jet : deltaPt : jetPt(2x) : pT hard bin (hr for area)
312 if(fbJetsArea){
313 entries = 1<<0 | 1<<1 | 1<<3 | 1<<18 | 1<<12 | 1<<13 | 1<<17 | 1<<19 | 1<<20 | 1<<21 | 1<<14 | 1<<6 | 1<<7 | 1<<26;
314 opt = 1<<18 | 1<<12 | 1<<13;
315 fhnJetsArea = NewTHnSparseF("fhnJetsArea", entries, opt);
316 }
317
318 // cent : nInpTrks : jetPt(3x) : deltaPt : delta : jetArea(3x) : fraction(2x) : pT hard bin
319 if(fbJets3Branches){
320 entries = 1<<0 | 1<<1 | 1<<6 | 1<<7 | 1<<7 | 1<<14 | 1<<14 | 1<<12 | 1<<13 | 1<<13 | 1<<19 | 1<<19 | 1<<26;
321 opt = 1<<6 | 1<<7 | 1<<14;
322 fhnJets3Branches = NewTHnSparseF("fhnJets3Branches", entries, opt);
323 }
324
325
326
327
328 //before cut
329
330 // cent : nInpTrks : rp bins : fraction : jetPt(2x) : jetEta(2x) : jetPhi(2x) (low resolution) (with fraction, eta, phi, pt cuts possible)
331 if(fbJetsBeforeCut1){
332 entries = 1<<0 | 1<<1 | 1<<3 | 1<<19 | 1<<6 | 1<<7 | 1<<8 | 1<<9 | 1<<10 | 1<<11;
333 opt = 0;
334 fhnJetsBeforeCut1 = NewTHnSparseF("fhnJetsBeforeCut1", entries, opt);
335 }
336
337 // cent : nInpTrks : rp bins : deltaPt : jetPt(2x) : deltaR : deltaEta : jetEta(2x) (low resolution)
338 if(fbJetsBeforeCut2){
339 entries = 1<<0 | 1<<1 | 1<<3 | 1<<14 | 1<<6 | 1<<7 | 1<<17 | 1<<15 | 1<<8 | 1<<9;
340 opt = 0;
341 fhnJetsBeforeCut2 = NewTHnSparseF("fhnJetsBeforeCut2", entries, opt);
342 }
343
344 fOutputList->Add(fHistEvtSelection);
345 fOutputList->Add(fHistJetSelection);
346 fOutputList->Add(fh2JetSelection);
347 if(fbEvent) fOutputList->Add(fhnEvent);
348 if(fbJetsMismatch1) fOutputList->Add(fhnJetsMismatch1);
349 if(fbJetsMismatch2) fOutputList->Add(fhnJetsMismatch2);
350 if(fbJetsRp) fOutputList->Add(fhnJetsRp);
351 if(fbJetsDeltaPt) fOutputList->Add(fhnJetsDeltaPt);
352 if(fbJetsEta) fOutputList->Add(fhnJetsEta);
353 if(fbJetsPhi) fOutputList->Add(fhnJetsPhi);
354 if(fbJetsArea) fOutputList->Add(fhnJetsArea);
355 if(fbJets3Branches) fOutputList->Add(fhnJets3Branches);
356 if(fbJetsBeforeCut1) fOutputList->Add(fhnJetsBeforeCut1);
357 if(fbJetsBeforeCut2) fOutputList->Add(fhnJetsBeforeCut2);
358
359 // =========== Switch on Sumw2 for all histos ===========
360 for (Int_t i=0; i<fOutputList->GetEntries(); ++i) {
361 TH1 *h1 = dynamic_cast<TH1*>(fOutputList->At(i));
362 if (h1){
363 h1->Sumw2();
364 continue;
365 }
366 THnSparse *hn = dynamic_cast<THnSparse*>(fOutputList->At(i));
367 if (hn){
368 hn->Sumw2();
369 }
370 }
371 TH1::AddDirectory(oldStatus);
372
373 PostData(1, fOutputList);
3c6a60f7 374}
375
376void AliAnalysisTaskJetResponseV2::UserExec(Option_t *)
377{
d5ec3aef 378 // load events, apply event cuts, then compare leading jets
379
380 if(!strlen(fJetBranchName[0].Data()) || !strlen(fJetBranchName[1].Data())){
381 AliError("Jet branch name not set.");
382 return;
383 }
384
385 fESD=dynamic_cast<AliESDEvent*>(InputEvent());
386 if (!fESD) {
387 AliError("ESD not available");
388 fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
389 } else {
390 fAOD = dynamic_cast<AliAODEvent*>(AODEvent());
391 }
392 if (!fAOD) {
393 AliError("AOD not available");
394 return;
395 }
396
397 // -- event selection --
398 fHistEvtSelection->Fill(1); // number of events before event selection
399
400 // physics selection
401 AliInputEventHandler* inputHandler = (AliInputEventHandler*)
402 ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
403 if(!(inputHandler->IsEventSelected() & fOfflineTrgMask)){
404 if(fDebug) Printf(" Trigger Selection: event REJECTED ... ");
405 fHistEvtSelection->Fill(2);
406 PostData(1, fOutputList);
407 return;
408 }
409
410 // vertex selection
411 AliAODVertex* primVtx = fAOD->GetPrimaryVertex();
412 Int_t nTracksPrim = primVtx->GetNContributors();
413 if ((nTracksPrim < fMinContribVtx) ||
414 (primVtx->GetZ() < fVtxZMin) ||
415 (primVtx->GetZ() > fVtxZMax) ){
416 if(fDebug) Printf("%s:%d primary vertex z = %f: event REJECTED...",(char*)__FILE__,__LINE__,primVtx->GetZ());
417 fHistEvtSelection->Fill(3);
418 PostData(1, fOutputList);
419 return;
420 }
421
422 // event class selection (from jet helper task)
423 Int_t eventClass = AliAnalysisHelperJetTasks::EventClass();
424 if(fDebug) Printf("Event class %d", eventClass);
425 if (eventClass < fEvtClassMin || eventClass > fEvtClassMax){
426 fHistEvtSelection->Fill(4);
427 PostData(1, fOutputList);
428 return;
429 }
430
431 // centrality selection
432 AliCentrality *cent = 0x0;
433 Float_t centValue = 0.;
434 if(fESD) cent = fESD->GetCentrality();
435 if(cent) centValue = cent->GetCentralityPercentile("V0M");
436 if(fDebug) printf("centrality: %f\n", centValue);
437 if (centValue < fCentMin || centValue > fCentMax){
438 fHistEvtSelection->Fill(4);
439 PostData(1, fOutputList);
440 return;
441 }
442
443
444 // multiplicity due to input tracks
445 Int_t nInputTracks = GetNInputTracks();
446
447 if (nInputTracks < fNInputTracksMin || (fNInputTracksMax > -1 && nInputTracks > fNInputTracksMax)){
448 fHistEvtSelection->Fill(5);
449 PostData(1, fOutputList);
450 return;
451 }
31b9d515 452
453
d5ec3aef 454 fHistEvtSelection->Fill(0); // accepted events
455 // -- end event selection --
456
457 // pt hard
458 Double_t pthard = AliAnalysisTaskFastEmbedding::GetPtHard();
459 Int_t pthardbin = GetPtHardBin(pthard);
460
461 // reaction plane
462 Double_t rp = AliAnalysisHelperJetTasks::ReactionPlane(kFALSE);
463 if(fbEvent){
464 Double_t eventEntries[3] = { (Double_t)centValue, (Double_t)nInputTracks, rp };
465 fhnEvent->Fill(eventEntries);
466 }
467
468
469 // get background
470 AliAODJetEventBackground* externalBackground = 0;
471 if(!externalBackground&&fBackgroundBranch.Length()){
472 externalBackground = (AliAODJetEventBackground*)(fAOD->FindListObject(fBackgroundBranch.Data()));
473 //if(!externalBackground)Printf("%s:%d Background branch not found %s",(char*)__FILE__,__LINE__,fBackgroundBranch.Data());;
474 }
475 Float_t rho = 0;
476 if(externalBackground)rho = externalBackground->GetBackground(0);
477
478
479 // fetch jets
480 TClonesArray *aodJets[3];
481 aodJets[0] = dynamic_cast<TClonesArray*>(fAOD->FindListObject(fJetBranchName[0].Data())); // in general: embedded jet
482 aodJets[1] = dynamic_cast<TClonesArray*>(fAOD->FindListObject(fJetBranchName[1].Data())); // in general: embedded jet + UE version1
483 if( strlen(fJetBranchName[2].Data()) )
484 aodJets[2] = dynamic_cast<TClonesArray*>(fAOD->FindListObject(fJetBranchName[2].Data())); // in general: embedded jet + UE version2
485
486 for (Int_t iJetType = 0; iJetType < fkNbranches; iJetType++) {
487 fListJets[iJetType]->Clear();
488 if (!aodJets[iJetType]) continue;
31b9d515 489
d5ec3aef 490 if(fDebug) Printf("%s: %d jets",fJetBranchName[iJetType].Data(),aodJets[iJetType]->GetEntriesFast());
491
492 for (Int_t iJet = 0; iJet < aodJets[iJetType]->GetEntriesFast(); iJet++) {
493 AliAODJet *jet = dynamic_cast<AliAODJet*>((*aodJets[iJetType])[iJet]);
494 if (jet) fListJets[iJetType]->Add(jet);
495 }
496 fListJets[iJetType]->Sort();
497 }
31b9d515 498
499
d5ec3aef 500 // jet matching
501 static TArrayI aMatchIndex(fListJets[0]->GetEntries());
502 static TArrayF aPtFraction(fListJets[0]->GetEntries());
503 if(aMatchIndex.GetSize()<fListJets[0]->GetEntries()) aMatchIndex.Set(fListJets[0]->GetEntries());
504 if(aPtFraction.GetSize()<fListJets[0]->GetEntries()) aPtFraction.Set(fListJets[0]->GetEntries());
505
506 // stores matched jets in 'aMatchIndex' and fraction of pT in 'aPtFraction'
507 AliAnalysisHelperJetTasks::GetJetMatching(fListJets[0], TMath::Min((Int_t)fNMatchJets,(Int_t)fListJets[0]->GetEntries()),
508 fListJets[1], TMath::Min((Int_t)fNMatchJets,(Int_t)fListJets[1]->GetEntries()),
509 aMatchIndex, aPtFraction, fDebug, fMatchMaxDist, fIsPbPb?1:2);
510 static TArrayI aMatchIndexv2(fListJets[0]->GetEntries());
511 static TArrayF aPtFractionv2(fListJets[0]->GetEntries());
512 if( strlen(fJetBranchName[2].Data()) ) {
513 if(aMatchIndexv2.GetSize()<fListJets[0]->GetEntries()) aMatchIndexv2.Set(fListJets[0]->GetEntries());
514 if(aPtFractionv2.GetSize()<fListJets[0]->GetEntries()) aPtFractionv2.Set(fListJets[0]->GetEntries());
515 AliAnalysisHelperJetTasks::GetJetMatching(fListJets[0], TMath::Min((Int_t)fNMatchJets,(Int_t)fListJets[0]->GetEntries()), fListJets[2], TMath::Min((Int_t)fNMatchJets,(Int_t)fListJets[2]->GetEntries()), aMatchIndexv2, aPtFractionv2, fDebug, fMatchMaxDist, fIsPbPb?1:2);
516 }
3c6a60f7 517
d5ec3aef 518 // loop over matched jets
519 Int_t ir = -1; // index of matched reconstruced jet
520 Int_t ir2 = -1; // index of matched reconstruced jet
521 Float_t fraction = -1.;
522 Float_t fraction2 = -1.;
523 AliAODJet *jet[3] = { 0x0, 0x0, 0x0 };
524 Float_t jetEta[3] = { -990., -990., -990. };
525 Float_t jetPhi[3] = { -990., -990., -990. };
526 Float_t jetPt[3] = { -990., -990., -990. };
527 Float_t jetArea[3] = { -990., -990., -990. };
528 Float_t rpJet[3] = { -990., -990., -990. };
529 Int_t rpBin = -990;
3c6a60f7 530
d5ec3aef 531 for(Int_t ig=0; ig<fListJets[0]->GetEntries(); ++ig){
532 ir = aMatchIndex[ig];
533 ir2 = aMatchIndexv2[ig];
534
535 //fetch jets
536 jet[0] = (AliAODJet*)(fListJets[0]->At(ig));
537 if(ir>=0) jet[1] = (AliAODJet*)(fListJets[1]->At(ir));
538 else jet[1] = 0x0;
539 if(ir2>=0) jet[2] = (AliAODJet*)(fListJets[2]->At(ir2));
540 else jet[2] = 0x0;
541
542 for(Int_t i=0; i<fkNbranches; ++i){
543 if(!jet[i]){
544 jetEta[i] = -990;
545 jetPhi[i] = -990.;
546 jetPt[i] = -990.;
547 jetArea[i] = -990.;
548 rpJet[i] = -990.;
549 if(i==1) rpBin = -990;
550 }
551 else {
552 jetEta[i] = jet[i]->Eta();
553 jetPhi[i] = jet[i]->Phi();
554 jetPt[i] = GetPt(jet[i], i);
555 jetArea[i] = jet[i]->EffectiveAreaCharged();
556 rpJet[i] = TVector2::Phi_mpi_pi(rp-jetPhi[i]);
557 if(i==1) rpBin = AliAnalysisHelperJetTasks::GetPhiBin(TVector2::Phi_mpi_pi(rp-jetPhi[i]), 3);
d314de7b 558 }
d5ec3aef 559 }
560 fraction = aPtFraction[ig];
561 fraction2 = aPtFractionv2[ig];
31b9d515 562
d5ec3aef 563 // jet statistics
564 fHistJetSelection->Fill(1); // all probe jets
565 if(jet[0]) fh2JetSelection->Fill(1.,jetPt[0]); // all probe jets
31b9d515 566
d5ec3aef 567 if(ir<0){
568 fHistJetSelection->Fill(2);
569 if(jet[0]) fh2JetSelection->Fill(2.,jetPt[0]);
31b9d515 570
d5ec3aef 571 if(fbJetsMismatch1){
572 if(!jet[0]) continue;
573 Double_t jetEntriesMismatch1[7] = {
574 (Double_t)centValue, (Double_t)nInputTracks, (Double_t)rpBin,
575 (Double_t)jetPt[0], (Double_t)jetEta[0], (Double_t)jetPhi[0], (Double_t)jetArea[0]
576 };
577 fhnJetsMismatch1->Fill(jetEntriesMismatch1);
31b9d515 578 }
d5ec3aef 579
580 continue;
581 }
31b9d515 582
d5ec3aef 583 if(!jet[0] || !jet[1]){
584 fHistJetSelection->Fill(3);
585 if(jet[0]) fh2JetSelection->Fill(3.,jetPt[0]);
586 continue;
587 }
31b9d515 588
d5ec3aef 589 // look for distance to next rec jet
590 Float_t distNextJet = -0.01; // no neighbor
591 Float_t ptNextJet = -1.;
592 for(Int_t r=0; r<fListJets[1]->GetEntries(); ++r){
593 if(r==ir) continue;
594 Float_t tmpDeltaR = jet[1]->DeltaR((AliAODJet*)fListJets[1]->At(r));
595 if(distNextJet<0. || distNextJet>tmpDeltaR){
596 distNextJet = tmpDeltaR;
597 if(fKeepJets) ptNextJet = ((AliAODJet*)fListJets[1]->At(r))->GetPtSubtracted(0);
598 else ptNextJet = ((AliAODJet*)fListJets[1]->At(r))->Pt();
d314de7b 599 }
d5ec3aef 600 }
31b9d515 601
602
603
d5ec3aef 604 //Float_t localRho = jetArea[1]>0. ? (jetPt[1]+rho*jetArea[1] - jetPt[0]) / jetArea[1] : 0.;
605 //Float_t relRho = rho>0. ? localRho / rho : 0.;
31b9d515 606
607
d5ec3aef 608 // calculate parameters of associated jets
609 /* from Leticia, kT clusterizer
610 Float_t par0[4] = { 0.00409, 0.01229, 0.05, 0.26 };
611 Float_t par1[4] = { -2.97035e-03, -2.03182e-03, -1.25702e-03, -9.95107e-04 };
612 Float_t par2[4] = { 1.02865e-01, 1.49039e-01, 1.53910e-01, 1.51109e-01 };
613 */
614 // own, from embedded tracks
615 Float_t par0[4] = { 0.02841, 0.05039, 0.09092, 0.24089 };
616 Float_t par1[4] = { -4.26725e-04, -1.15273e-03, -1.56827e-03, -3.08003e-03 };
617 Float_t par2[4] = { 4.95415e-02, 9.79538e-02, 1.32814e-01, 1.71743e-01 };
31b9d515 618
d5ec3aef 619 Float_t rpCorr = 0.;
31b9d515 620
d5ec3aef 621 if(eventClass>0&&eventClass<4){
622 rpCorr = par0[eventClass-1] + par1[eventClass-1] * 2*TMath::Cos(rpJet[1]) + par2[eventClass-1] * 2*TMath::Cos(2*rpJet[1]);
623 rpCorr *= rho * jetArea[1];
624 }
625
626 Float_t delta = jetPt[2]-jetPt[1];
627 Float_t deltaPt = jetPt[1]-jetPt[0];
628 Float_t deltaEta = jetEta[1]-jetEta[0];
629 Float_t deltaPhi = TVector2::Phi_mpi_pi(jetPhi[1]-jetPhi[0]);
630 Float_t deltaR = TMath::Sqrt(deltaEta*deltaEta+deltaPhi*deltaPhi);
631 Float_t deltaArea = jetArea[1]-jetArea[0];
31b9d515 632
633
d5ec3aef 634 // fill thnsparse before acceptance cut
635 if(fbJetsBeforeCut1){
636 Double_t jetBeforeCutEntries1[10] = {
637 (Double_t)centValue, (Double_t)nInputTracks, (Double_t)rpBin,
638 (Double_t)jetPt[0], (Double_t)jetPt[1], (Double_t)jetEta[0], (Double_t)jetEta[1],
639 (Double_t)jetPhi[0], (Double_t)jetPhi[1], (Double_t)fraction
640 };
641 fhnJetsBeforeCut1->Fill(jetBeforeCutEntries1);
642 }
31b9d515 643
d5ec3aef 644 if(fbJetsBeforeCut2){
645 Double_t jetBeforeCutEntries2[10] = {
646 (Double_t)centValue, (Double_t)nInputTracks, (Double_t)rpBin,
647 (Double_t)jetPt[0], (Double_t)jetPt[1],
648 (Double_t)jetEta[0], (Double_t)jetEta[1],
649 (Double_t)deltaPt, (Double_t)deltaEta, (Double_t)deltaR
650 };
651 fhnJetsBeforeCut2->Fill(jetBeforeCutEntries2);
652 }
31b9d515 653
d5ec3aef 654 // jet selection
655 Bool_t jetAccepted = kTRUE;
656 // minimum fraction required
657 if(fraction<fJetPtFractionMin){
658 fHistJetSelection->Fill(4);
659 fh2JetSelection->Fill(4.,jetPt[0]);
660 jetAccepted = kFALSE;
661 }
31b9d515 662
d5ec3aef 663 if(jetAccepted){
664 // jet acceptance + minimum pT check
665 if(jetEta[0]>fJetEtaMax || jetEta[0]<fJetEtaMin ||
666 jetEta[1]>fJetEtaMax || jetEta[1]<fJetEtaMin){
31b9d515 667
d5ec3aef 668 if(fDebug){
669 Printf("Jet not in eta acceptance.");
670 Printf("[0]: jet %d eta %.2f", ig, jetEta[0]);
671 Printf("[1]: jet %d eta %.2f", ir, jetEta[1]);
672 }
673 fHistJetSelection->Fill(5);
674 fh2JetSelection->Fill(5.,jetPt[0]);
675 jetAccepted = kFALSE;
31b9d515 676 }
d5ec3aef 677 }
678 if(jetAccepted){
679 if(jetPt[0] < fJetPtMin){
680 if(fDebug) Printf("Jet %d (pT %.1f GeV/c) has less than required pT.", ir, jetPt[0]);
681 fHistJetSelection->Fill(6);
682 fh2JetSelection->Fill(6.,jetPt[0]);
683 jetAccepted = kFALSE;
31b9d515 684 }
d5ec3aef 685 }
31b9d515 686
d5ec3aef 687 if(jetAccepted){
688 if(jet[1]->Trigger()&fJetTriggerExcludeMask){
689 fHistJetSelection->Fill(7);
690 fh2JetSelection->Fill(7.,jetPt[0]);
691 jetAccepted = kFALSE;
31b9d515 692 }
d5ec3aef 693
694 }
31b9d515 695
d5ec3aef 696 if(!jetAccepted){
697 if(fbJetsMismatch2){
698 Double_t jetEntriesMismatch2[10] = {
699 (Double_t)centValue, (Double_t)nInputTracks, (Double_t)rpBin,
700 (Double_t)jetPt[0], (Double_t)jetPt[1],
701 (Double_t)jetEta[0], (Double_t)jetEta[1],
702 (Double_t)deltaEta, (Double_t)deltaR,
703 (Double_t)fraction
704 };
705 fhnJetsMismatch2->Fill(jetEntriesMismatch2);
31b9d515 706 }
d5ec3aef 707 continue;
708 }
31b9d515 709
d5ec3aef 710 // all accepted jets
711 fHistJetSelection->Fill(0);
712 fh2JetSelection->Fill(0.,jetPt[0]);
31b9d515 713
d5ec3aef 714 // fill thnsparse
715 if(fbJetsRp){
716 Double_t jetEntriesRp[11] = {
717 (Double_t)centValue, (Double_t)nInputTracks, (Double_t) rp,
718 (Double_t)rpBin, (Double_t)rpJet[0], (Double_t)rpJet[1],
719 (Double_t)jetPt[0], (Double_t)deltaPt, (Double_t)rho, (Double_t)rpCorr,
720 (Double_t)pthardbin
721 };
722 fhnJetsRp->Fill(jetEntriesRp);
723 }
31b9d515 724
d5ec3aef 725 if(fbJetsDeltaPt){
726 Double_t jetEntriesDeltaPt[8] = {
727 (Double_t)centValue, (Double_t)nInputTracks, (Double_t)rpBin,
728 (Double_t)jetPt[0], (Double_t)jetPt[1], (Double_t)deltaPt, (Double_t)deltaArea,
729 (Double_t)pthardbin
730 };
731 fhnJetsDeltaPt->Fill(jetEntriesDeltaPt);
732 }
31b9d515 733
d5ec3aef 734 if(fbJetsEta){
735 Double_t jetEntriesEta[11] = {
736 (Double_t)centValue, (Double_t)nInputTracks, (Double_t)rpBin,
737 (Double_t)jetPt[0], (Double_t)jetPt[1], (Double_t)jetEta[0], (Double_t)jetEta[1],
738 (Double_t)deltaPt, (Double_t)deltaEta, (Double_t)deltaR,
739 (Double_t)pthardbin
740 };
741 fhnJetsEta->Fill(jetEntriesEta);
742 }
31b9d515 743
d5ec3aef 744 if(fbJetsPhi){
745 Double_t jetEntriesPhi[10] = {
746 (Double_t)centValue, (Double_t)nInputTracks, (Double_t)rpBin,
747 (Double_t)jetPt[0], (Double_t)jetPt[1], (Double_t)jetPhi[0], (Double_t)jetPhi[1],
748 (Double_t)deltaPt, (Double_t)deltaPhi,
749 (Double_t)pthardbin
750 };
751 fhnJetsPhi->Fill(jetEntriesPhi);
752 }
31b9d515 753
d5ec3aef 754 if(fbJetsArea){
755 Double_t jetEntriesArea[14] = {
756 (Double_t)centValue, (Double_t)nInputTracks, (Double_t)rpBin,
757 (Double_t)jetPt[0], (Double_t)jetPt[1], (Double_t)jetArea[0], (Double_t)jetArea[1],
758 (Double_t)deltaPt, (Double_t)deltaR, (Double_t)deltaArea,
759 (Double_t)fraction, (Double_t)distNextJet, (Double_t)ptNextJet,
760 (Double_t)pthardbin
761 };
762 fhnJetsArea->Fill(jetEntriesArea);
763 }
764
765 if(fbJets3Branches){
766 Double_t jetEntries3Branches[14] = {
767 (Double_t)centValue, (Double_t)nInputTracks,
768 (Double_t)jetPt[0], (Double_t)jetPt[1], (Double_t)jetPt[2],
769 (Double_t)deltaPt, (Double_t)delta,
770 (Double_t)jetArea[0], (Double_t)jetArea[1],(Double_t)jetArea[3],
771 (Double_t)fraction, (Double_t)fraction2,
772 (Double_t)pthardbin
773 };
774 fhnJets3Branches->Fill(jetEntries3Branches);
775 }
31b9d515 776
d5ec3aef 777 }
31b9d515 778
d5ec3aef 779 PostData(1, fOutputList);
3c6a60f7 780}
781
782void AliAnalysisTaskJetResponseV2::Terminate(const Option_t *)
783{
d5ec3aef 784 // Draw result to the screen
785 // Called once at the end of the query
3c6a60f7 786
d5ec3aef 787 if (!GetOutputData(1))
788 return;
3c6a60f7 789}
790
791Int_t AliAnalysisTaskJetResponseV2::GetNInputTracks()
792{
d5ec3aef 793 // number of tracks in the event, obtained via jet finder
794
795 Int_t nInputTracks = 0;
796
797 TString jbname(fJetBranchName[1]);
798 //needs complete event, use jets without background subtraction
799 for(Int_t i=1; i<=3; ++i){
800 if(jbname.Contains(Form("B%d",i))) jbname.ReplaceAll(Form("B%d",i),"B0");
801 }
802 // use only HI event
803 if(jbname.Contains("AODextraonly")) jbname.ReplaceAll("AODextraonly","AOD");
804 if(jbname.Contains("AODextra")) jbname.ReplaceAll("AODextra","AOD");
805
806 if(fDebug) Printf("Multiplicity from jet branch %s", jbname.Data());
807 TClonesArray *tmpAODjets = dynamic_cast<TClonesArray*>(fAOD->FindListObject(jbname.Data()));
808 if(!tmpAODjets){
809 Printf("Jet branch %s not found", jbname.Data());
810 Printf("AliAnalysisTaskJetResponseV2::GetNInputTracks FAILED");
811 return -1;
812 }
31b9d515 813
d5ec3aef 814 for (Int_t iJet=0; iJet<tmpAODjets->GetEntriesFast(); iJet++){
815 AliAODJet *jet = dynamic_cast<AliAODJet*>((*tmpAODjets)[iJet]);
816 if(!jet) continue;
817 TRefArray *trackList = jet->GetRefTracks();
818 Int_t nTracks = trackList->GetEntriesFast();
819 nInputTracks += nTracks;
820 if(fDebug) Printf("#jet%d: %d tracks", iJet, nTracks);
821 }
822 if(fDebug) Printf("---> input tracks: %d", nInputTracks);
823
824 return nInputTracks;
3c6a60f7 825}
d314de7b 826
827THnSparse* AliAnalysisTaskJetResponseV2::NewTHnSparseF(const char* name, UInt_t entries, UInt_t opt)
828{
d5ec3aef 829 // generate new THnSparseF, axes are defined in GetDimParams()
830
831 Int_t count = 0;
832 UInt_t tmp = entries;
833 while(tmp!=0){
834 count++;
835 tmp = tmp &~ -tmp; // clear lowest bit
836 }
837
838 TString hnTitle(name);
839 const Int_t dim = count;
840 Int_t nbins[dim];
841 Double_t xmin[dim];
842 Double_t xmax[dim];
843
844 Int_t i=0;
845 Int_t c=0;
846 while(c<dim && i<32){
847 if(entries&(1<<i)){
848 Bool_t highres = opt&(1<<i);
849 TString label("");
850 GetDimParams(i, highres, label, nbins[c], xmin[c], xmax[c]);
851 hnTitle += Form(";%s",label.Data());
852 c++;
853 }
31b9d515 854
d5ec3aef 855 i++;
856 }
857 hnTitle += ";";
31b9d515 858
d5ec3aef 859 return new THnSparseF(name, hnTitle.Data(), dim, nbins, xmin, xmax);
d314de7b 860}
861
862void AliAnalysisTaskJetResponseV2::GetDimParams(Int_t iEntry, Bool_t hr, TString &label, Int_t &nbins, Double_t &xmin, Double_t &xmax)
863{
d5ec3aef 864 // stores label and binning of axis for THnSparse
d314de7b 865
d5ec3aef 866 const Double_t pi = TMath::Pi();
d314de7b 867
d5ec3aef 868 switch(iEntry){
31b9d515 869
d5ec3aef 870 case 0:
871 label = "V0 centrality (%)";
872 if(hr){
873 nbins = 100;
874 xmin = 0.;
875 xmax = 100.;
876 } else {
877 nbins = 8;
878 xmin = 0.;
879 xmax = 80.;
880 }
881 break;
d314de7b 882
883
d5ec3aef 884 case 1:
885 label = "nb. of input tracks";
886 if(fIsPbPb){
887 if(hr){
888 nbins = 400;
889 xmin = 0.;
890 xmax = 4000.;
31b9d515 891 } else {
d5ec3aef 892 nbins = 40;
893 xmin = 0.;
894 xmax = 4000.;
31b9d515 895 }
d5ec3aef 896 } else {
897 nbins = 40;
898 xmin = 0.;
899 xmax = 400.;
900 }
901 break;
d314de7b 902
903
d5ec3aef 904 case 2:
905 label = "event plane #psi";
906 if(hr){
907 nbins = 30;
908 xmin = 0.;
909 xmax = pi;
910 } else {
911 nbins = 30;
912 xmin = 0.;
913 xmax = pi;
914 }
915 break;
d314de7b 916
917
d5ec3aef 918 case 3:
919 label = "event plane bin";
920 nbins = 3;
921 xmin = -.5;
922 xmax = 2.5;
923 break;
d314de7b 924
925
d5ec3aef 926 case 4:
927 case 5:
928 if(iEntry==4)label = "#Delta#phi(RP-jet) (probe)";
929 if(iEntry==5)label = "#Delta#phi(RP-jet) (rec)";
930 nbins = 48;
931 xmin = -pi;
932 xmax = pi;
933 break;
d314de7b 934
935
d5ec3aef 936 case 6:
937 case 7:
938 if(iEntry==6)label = "probe p_{T} (GeV/c)";
939 if(iEntry==7)label = "rec p_{T} (GeV/c)";
940 if(hr){
941 nbins = 300;
942 xmin = -50.;
943 xmax = 250.;
944 } else {
945 nbins = 50;
946 xmin = 0.;
947 xmax = 250.;
948 }
949 break;
d314de7b 950
951
d5ec3aef 952 case 8:
953 case 9:
954 if(iEntry==8)label = "probe #eta";
955 if(iEntry==9)label = "rec #eta";
956 if(hr){
957 nbins = 56;
958 xmin = -.7;
959 xmax = .7;
960 } else {
961 nbins = 28;
962 xmin = -.7;
963 xmax = .7;
964 }
965 break;
d314de7b 966
967
d5ec3aef 968 case 10:
969 case 11:
970 if(iEntry==10)label = "probe #phi";
971 if(iEntry==11)label = "rec #phi";
972 if(hr){
973 nbins = 90; // modulo 18 (sectors)
974 xmin = 0.;
975 xmax = 2*pi;
976 } else {
977 nbins = 25;
978 xmin = 0.;
979 xmax = 2*pi;
980 }
981 break;
d314de7b 982
983
d5ec3aef 984 case 12:
985 case 13:
986 if(iEntry==12)label = "probe area";
987 if(iEntry==13)label = "rec area";
988 if(hr){
989 nbins = 100;
990 xmin = 0.;
991 xmax = 1.;
992 } else {
993 nbins = 25;
994 xmin = 0.;
995 xmax = 1.;
996 }
997 break;
d314de7b 998
d5ec3aef 999 case 14:
1000 label = "#Delta p_{T}";
1001 if(hr){
1002 nbins = 241;
1003 xmin = -120.5;
1004 xmax = 120.5;
1005 } else {
1006 nbins = 101;
1007 xmin = -101.;
1008 xmax = 101.;
1009 }
1010 break;
d314de7b 1011
d5ec3aef 1012 case 15:
1013 label = "#Delta#eta";
1014 if(hr){
1015 nbins = 51;
1016 xmin = -1.02;
1017 xmax = 1.02;
1018 } else {
1019 nbins = 51;
1020 xmin = -1.02;
1021 xmax = 1.02;
1022 }
1023 break;
d314de7b 1024
1025
d5ec3aef 1026 case 16:
1027 label = "#Delta#phi";
1028 if(hr){
1029 nbins = 45;
1030 xmin = -pi;
1031 xmax = pi;
1032 } else {
1033 nbins = 45;
1034 xmin = -pi;
1035 xmax = pi;
1036 }
1037 break;
d314de7b 1038
1039
d5ec3aef 1040 case 17:
1041 label = "#DeltaR";
1042 if(hr){
1043 nbins = 50;
1044 xmin = 0.;
1045 xmax = 1.;
1046 } else {
1047 nbins = 50;
1048 xmin = 0.;
1049 xmax = 1.;
1050 }
1051 break;
d314de7b 1052
1053
d5ec3aef 1054 case 18:
1055 label = "#Deltaarea";
1056 if(hr){
1057 nbins = 81;
1058 xmin = -.81;
1059 xmax = .81;
1060 } else {
1061 nbins = 33;
1062 xmin = -.825;
1063 xmax = .825;
1064 }
1065 break;
d314de7b 1066
1067
d5ec3aef 1068 case 19:
1069 label = "fraction";
1070 if(hr){
1071 nbins = 52;
1072 xmin = 0.;
1073 xmax = 1.04;
1074 } else {
1075 nbins = 52;
1076 xmin = 0.;
1077 xmax = 1.04;
1078 }
1079 break;
d314de7b 1080
1081
d5ec3aef 1082 case 20:
1083 label = "distance to closest rec jet";
1084 if(hr){
1085 nbins = 51;
1086 xmin = -0.02;
1087 xmax = 1.;
1088 } else {
1089 nbins = 51;
1090 xmin = -0.02;
1091 xmax = 1.;
1092 }
1093 break;
d314de7b 1094
1095
d5ec3aef 1096 case 21:
1097 label = "p_{T} of closest rec jet";
1098 nbins = 100;
1099 xmin = 0.;
1100 xmax = 200.;
1101 break;
d314de7b 1102
d5ec3aef 1103 case 22:
1104 label = "#rho";
1105 nbins = 125;
1106 xmin = 0.;
1107 xmax = 250.;
1108 break;
d314de7b 1109
d5ec3aef 1110 case 23:
1111 label = "abs. correction of #rho for RP";
1112 nbins = 51;
1113 xmin = -51.;
1114 xmax = 51.;
1115 break;
d314de7b 1116
d5ec3aef 1117 case 24:
1118 label = "local #rho";
1119 nbins = 125;
1120 xmin = 0.;
1121 xmax = 250.;
1122 break;
d314de7b 1123
d5ec3aef 1124 case 25:
1125 label = "local #rho / #rho";
1126 nbins = 500;
1127 xmin = 0.;
1128 xmax = 5.;
1129 break;
31b9d515 1130
d5ec3aef 1131 case 26:
1132 label = "p_{T,hard} bin";
1133 nbins = 10;
1134 xmin = -.5;
1135 xmax = 9.5;
1136 break;
31b9d515 1137
d5ec3aef 1138 }
31b9d515 1139
1140}
1141
1142//____________________________________________________________________________
1143Int_t AliAnalysisTaskJetResponseV2::GetPtHardBin(Double_t ptHard){
1144
d5ec3aef 1145 // returns pt hard bin (for LHC10e14, LHC11a1x, LHC11a2x simulations)
7a88ef30 1146
d5ec3aef 1147 const Int_t nBins = 10;
1148 Double_t binLimits[nBins] = { 5., 11., 21., 36., 57., 84., 117., 156., 200., 249. }; // lower limits
d314de7b 1149
d5ec3aef 1150 Int_t bin = -1;
1151 while(bin<nBins-1 && binLimits[bin+1]<ptHard){
1152 bin++;
1153 }
31b9d515 1154
d5ec3aef 1155 return bin;
31b9d515 1156}
1157
1158//____________________________________________________________________________
1159Double_t AliAnalysisTaskJetResponseV2::GetPt(AliAODJet *j, Int_t mode=0){
1160
d5ec3aef 1161 // returns jet pt, also negative pt after background subtraction if available
7a88ef30 1162
d5ec3aef 1163 Double_t pt = 0.;
31b9d515 1164
d5ec3aef 1165 if(fKeepJets && mode>0){ // background subtracted pt, can be negative
1166 pt = j->GetPtSubtracted(0);
1167 }
1168 else{
1169 pt = j->Pt(); // if negative, pt=0.01
1170 }
d314de7b 1171
d5ec3aef 1172 return pt;
d314de7b 1173}