]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGJE/AliAnalysisTaskJetResponse.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Transition PWG0 -> PWGUD
[u/mrichter/AliRoot.git] / PWGJE / AliAnalysisTaskJetResponse.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//
17// task compares jets in two branches,
18// written for analysis of jet embedding in HI events
19//
20// newer class: AliAnalysisTaskJetResponseV2
21//
22
2a436fa1 23#include "TChain.h"
24#include "TTree.h"
25#include "TMath.h"
26#include "TH1F.h"
27#include "TH2F.h"
46465e39 28#include "TH3F.h"
2a436fa1 29#include "TCanvas.h"
30
31#include "AliLog.h"
32
33#include "AliAnalysisTask.h"
34#include "AliAnalysisManager.h"
35
36#include "AliVEvent.h"
37#include "AliESDEvent.h"
38#include "AliESDInputHandler.h"
39#include "AliCentrality.h"
40#include "AliAnalysisHelperJetTasks.h"
41#include "AliInputEventHandler.h"
42
43#include "AliAODEvent.h"
44#include "AliAODJet.h"
45
46#include "AliAnalysisTaskJetResponse.h"
47
48ClassImp(AliAnalysisTaskJetResponse)
49
50AliAnalysisTaskJetResponse::AliAnalysisTaskJetResponse() :
51 AliAnalysisTaskSE(),
52 fESD(0x0),
53 fAOD(0x0),
54 fOfflineTrgMask(AliVEvent::kAny),
55 fMinContribVtx(1),
56 fVtxZMin(-8.),
57 fVtxZMax(8.),
46465e39 58 fEvtClassMin(1),
59 fEvtClassMax(4),
2a436fa1 60 fCentMin(0.),
61 fCentMax(100.),
3c6a60f7 62 fNInputTracksMin(0),
63 fNInputTracksMax(-1),
64 fReactionPlaneBin(-1),
2a436fa1 65 fJetEtaMin(-.5),
66 fJetEtaMax(.5),
46465e39 67 fJetPtMin(20.),
68 fJetPtFractionMin(0.5),
69 fNMatchJets(4),
70 //fJetDeltaEta(.4),
71 //fJetDeltaPhi(.4),
3c6a60f7 72 fEvtClassMode(0),
2a436fa1 73 fkNbranches(2),
3c6a60f7 74 fkEvtClasses(12),
2a436fa1 75 fOutputList(0x0),
76 fHistEvtSelection(0x0),
46465e39 77 fHistEvtClass(0x0),
78 fHistCentrality(0x0),
3c6a60f7 79 fHistNInputTracks(0x0),
80 //fHistNInputTracks2(0x0),
81 fHistCentVsTracks(0x0),
82 fHistReactionPlane(0x0),
83 fHistReactionPlaneWrtJets(0x0),
46465e39 84 fHistPtJet(new TH1F*[fkNbranches]),
85 fHistEtaPhiJet(new TH2F*[fkNbranches]),
86 fHistEtaPhiJetCut(new TH2F*[fkNbranches]),
87 fHistDeltaEtaDeltaPhiJet(new TH2F*[fkEvtClasses]),
88 fHistDeltaEtaDeltaPhiJetCut(new TH2F*[fkEvtClasses]),
3c6a60f7 89 //fHistDeltaEtaDeltaPhiJetNOMatching(new TH2F*[fkEvtClasses]),
46465e39 90 fHistDeltaEtaEtaJet(new TH2F*[fkEvtClasses]),
91 fHistDeltaPtEtaJet(new TH2F*[fkEvtClasses]),
92 fHistPtFraction(new TH2F*[fkEvtClasses]),
2a436fa1 93 fHistPtResponse(new TH2F*[fkEvtClasses]),
a11972e9 94 fHistPtSmearing(new TH2F*[fkEvtClasses]),
46465e39 95 fHistDeltaR(new TH2F*[fkEvtClasses]),
3c6a60f7 96 fHistArea(new TH3F*[fkEvtClasses]),
97 //fHistJetPtArea(new TH2F*[fkEvtClasses]),
98 fHistDeltaArea(new TH2F*[fkEvtClasses]),
99 fHistJetPtDeltaArea(new TH2F*[fkEvtClasses])
2a436fa1 100{
101 // default Constructor
102
103 fJetBranchName[0] = "";
104 fJetBranchName[1] = "";
105
106 fListJets[0] = new TList;
107 fListJets[1] = new TList;
108}
109
110AliAnalysisTaskJetResponse::AliAnalysisTaskJetResponse(const char *name) :
111 AliAnalysisTaskSE(name),
112 fESD(0x0),
113 fAOD(0x0),
114 fOfflineTrgMask(AliVEvent::kAny),
115 fMinContribVtx(1),
116 fVtxZMin(-8.),
117 fVtxZMax(8.),
46465e39 118 fEvtClassMin(1),
119 fEvtClassMax(4),
2a436fa1 120 fCentMin(0.),
121 fCentMax(100.),
3c6a60f7 122 fNInputTracksMin(0),
123 fNInputTracksMax(-1),
124 fReactionPlaneBin(-1),
2a436fa1 125 fJetEtaMin(-.5),
126 fJetEtaMax(.5),
46465e39 127 fJetPtMin(20.),
128 fJetPtFractionMin(0.5),
129 fNMatchJets(4),
3c6a60f7 130 fEvtClassMode(0),
2a436fa1 131 fkNbranches(2),
3c6a60f7 132 fkEvtClasses(12),
2a436fa1 133 fOutputList(0x0),
134 fHistEvtSelection(0x0),
46465e39 135 fHistEvtClass(0x0),
136 fHistCentrality(0x0),
3c6a60f7 137 fHistNInputTracks(0x0),
138 //fHistNInputTracks2(0x0),
139 fHistCentVsTracks(0x0),
140 fHistReactionPlane(0x0),
141 fHistReactionPlaneWrtJets(0x0),
46465e39 142 fHistPtJet(new TH1F*[fkNbranches]),
143 fHistEtaPhiJet(new TH2F*[fkNbranches]),
144 fHistEtaPhiJetCut(new TH2F*[fkNbranches]),
145 fHistDeltaEtaDeltaPhiJet(new TH2F*[fkEvtClasses]),
146 fHistDeltaEtaDeltaPhiJetCut(new TH2F*[fkEvtClasses]),
3c6a60f7 147 //fHistDeltaEtaDeltaPhiJetNOMatching(new TH2F*[fkEvtClasses]),
46465e39 148 fHistDeltaEtaEtaJet(new TH2F*[fkEvtClasses]),
149 fHistDeltaPtEtaJet(new TH2F*[fkEvtClasses]),
150 fHistPtFraction(new TH2F*[fkEvtClasses]),
2a436fa1 151 fHistPtResponse(new TH2F*[fkEvtClasses]),
a11972e9 152 fHistPtSmearing(new TH2F*[fkEvtClasses]),
46465e39 153 fHistDeltaR(new TH2F*[fkEvtClasses]),
3c6a60f7 154 fHistArea(new TH3F*[fkEvtClasses]),
155 //fHistJetPtArea(new TH2F*[fkEvtClasses]),
156 fHistDeltaArea(new TH2F*[fkEvtClasses]),
157 fHistJetPtDeltaArea(new TH2F*[fkEvtClasses])
2a436fa1 158{
159 // Constructor
160
161 fJetBranchName[0] = "";
162 fJetBranchName[1] = "";
163
164 fListJets[0] = new TList;
165 fListJets[1] = new TList;
166
167 DefineOutput(1, TList::Class());
168}
169
170AliAnalysisTaskJetResponse::~AliAnalysisTaskJetResponse()
171{
172 delete fListJets[0];
173 delete fListJets[1];
174}
175
176void AliAnalysisTaskJetResponse::SetBranchNames(const TString &branch1, const TString &branch2)
177{
178 fJetBranchName[0] = branch1;
179 fJetBranchName[1] = branch2;
180}
181
182void AliAnalysisTaskJetResponse::Init()
183{
184
185 // check for jet branches
186 if(!strlen(fJetBranchName[0].Data()) || !strlen(fJetBranchName[1].Data())){
187 AliError("Jet branch name not set.");
188 }
189
190}
191
192void AliAnalysisTaskJetResponse::UserCreateOutputObjects()
193{
194 // Create histograms
195 // Called once
46465e39 196 OpenFile(1);
197 if(!fOutputList) fOutputList = new TList;
198 fOutputList->SetOwner(kTRUE);
199
200 Bool_t oldStatus = TH1::AddDirectoryStatus();
201 TH1::AddDirectory(kFALSE);
202
2a436fa1 203
3c6a60f7 204 fHistEvtSelection = new TH1I("fHistEvtSelection", "event selection", 6, -0.5, 5.5);
2a436fa1 205 fHistEvtSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
206 fHistEvtSelection->GetXaxis()->SetBinLabel(2,"events IN");
207 fHistEvtSelection->GetXaxis()->SetBinLabel(3,"event selection (rejected)");
208 fHistEvtSelection->GetXaxis()->SetBinLabel(4,"vertex cut (rejected)");
209 fHistEvtSelection->GetXaxis()->SetBinLabel(5,"centrality (rejected)");
3c6a60f7 210 fHistEvtSelection->GetXaxis()->SetBinLabel(6,"multiplicity (rejected)");
46465e39 211
3c6a60f7 212 fHistEvtClass = new TH1I("fHistEvtClass", "event classes (centrality) ", fkEvtClasses, -0.5, fkEvtClasses-0.5);
46465e39 213 fHistCentrality = new TH1F("fHistCentrality", "event centrality", 100, 0., 100.);
3c6a60f7 214 fHistNInputTracks = new TH1F("fHistNInputTracks", "nb. of input tracks", 400, 0., 4000.);
215 //fHistNInputTracks2 = new TH2F("fHistNInputTracks2", "nb. of input tracks;from B0;from Bx", 400, 0., 4000., 400, 0., 4000.);
216 fHistCentVsTracks = new TH2F("fHistCentVsTracks", "centrality vs. nb. of input tracks;centrality;input tracks", 100, 0., 100., 400, 0., 4000.);
217 fHistReactionPlane = new TH1F("fHistReactionPlane", "reaction plane", 30, 0., TMath::Pi());
218 fHistReactionPlaneWrtJets = new TH1F("fHistReactionPlaneWrtJets", "reaction plane wrt jets", 48, -TMath::Pi(), TMath::Pi());
2a436fa1 219
220 for (Int_t iBranch = 0; iBranch < fkNbranches; iBranch++) {
46465e39 221 fHistPtJet[iBranch] = new TH1F(Form("fHistPtJet%i", iBranch),
222 Form("p_{T} of jets from branch %i;p_{T} (GeV/c);counts", iBranch),
223 250, 0., 250.);
224 fHistPtJet[iBranch]->SetMarkerStyle(kFullCircle);
225
226 fHistEtaPhiJet[iBranch] = new TH2F(Form("fHistEtaPhiJet%i", iBranch),
227 Form("#eta - #phi of jets from branch %i (before cut);#eta;#phi", iBranch),
3c6a60f7 228 28, -0.7, 0.7, 100, 0., 2*TMath::Pi());
46465e39 229 fHistEtaPhiJetCut[iBranch] = new TH2F(Form("fHistEtaPhiJetCut%i", iBranch),
230 Form("#eta - #phi of jets from branch %i (before cut);#eta;#phi", iBranch),
3c6a60f7 231 28, -0.7, 0.7, 100, 0., 2*TMath::Pi());
2a436fa1 232
233 }
46465e39 234
2a436fa1 235 for (Int_t iEvtClass =0; iEvtClass < fkEvtClasses; iEvtClass++) {
46465e39 236 fHistDeltaEtaDeltaPhiJet[iEvtClass] = new TH2F(Form("fHistDeltaEtaDeltaPhiJet%i",iEvtClass),
237 "#Delta#eta - #Delta#phi of jet;#Delta#eta;#Delta#phi",
238 101, -1.01, 1.01, 101, -1.01, 1.01);
239 fHistDeltaEtaDeltaPhiJetCut[iEvtClass] = new TH2F(Form("fHistDeltaEtaDeltaPhiJetCut%i",iEvtClass),
240 "#Delta#eta - #Delta#phi of jet;#Delta#eta;#Delta#phi",
241 101, -1.01, 1.01, 101, -1.01, 1.01);
3c6a60f7 242 //fHistDeltaEtaDeltaPhiJetNOMatching[iEvtClass] = new TH2F("fHistDeltaEtaDeltaPhiJetNOMatching",
243 // "#Delta#eta - #Delta#phi of jets which do not match;#Delta#eta;#Delta#phi",
244 // 100, -2., 2., 100, TMath::Pi(), TMath::Pi());
2a436fa1 245 fHistPtResponse[iEvtClass] = new TH2F(Form("pt_response%i",iEvtClass), Form("pt_response%i;p_{T} (GeV/c);p_{T} (GeV/c)",iEvtClass),
46465e39 246 250, 0., 250., 250, 0., 250.);
247 fHistPtSmearing[iEvtClass] = new TH2F(Form("pt_smearing%i",iEvtClass), Form("pt_smearing%i;#Deltap_{T} (GeV/c);p_{T} (GeV/c)",iEvtClass),
3c6a60f7 248 161, -80.5, 80.5, 250, 0., 250.);
249 fHistDeltaR[iEvtClass] = new TH2F(Form("hist_DeltaR%i",iEvtClass), "#DeltaR of matched jets;#Deltap_{T};#DeltaR", 161, -80.5,80.5, 85, 0.,3.4);
250 fHistArea[iEvtClass] = new TH3F(Form("hist_Area%i",iEvtClass), "jet area;probe p_{T};#Deltap_{T};jet area", 250, 0., 250., 161, -80.5,80.5, 100, 0.,1.0);
251 //fHistJetPtArea[iEvtClass] = new TH2F(Form("hist_JetPtArea%i",iEvtClass), "jet area vs. probe p_{T};jet p_{T};jet area", 250, 0.,250., 100, 0.,1.0);
252 fHistDeltaArea[iEvtClass] = new TH2F(Form("hist_DeltaArea%i",iEvtClass), "#DeltaArea of matched jets;#Deltap_{T};#Deltaarea", 161, -80.5, 80.5, 32, -.8, .8);
253 fHistJetPtDeltaArea[iEvtClass] = new TH2F(Form("hist_JetPtDeltaArea%i",iEvtClass), "#DeltaArea of matched jets vs. probe jet p_{T};#Deltap_{T};#Deltaarea", 250, 0., 250., 32, -.8, .8);
46465e39 254
255 fHistDeltaEtaEtaJet[iEvtClass] = new TH2F(Form("fHistDeltaEtaEtaJet%i", iEvtClass),
256 Form("#eta - #Delta#eta of matched jets from event class %i;#eta;#Delta#eta", iEvtClass),
3c6a60f7 257 28, -.7, .7, 101, -1.01, 1.01);
46465e39 258 fHistDeltaPtEtaJet[iEvtClass] = new TH2F(Form("fHistDeltaPtEtaJet%i", iEvtClass),
259 Form("#eta - #Deltap_{T} of matched jets from event class %i;#eta;#Deltap_{T}", iEvtClass),
3c6a60f7 260 28, -.7, .7, 161, -80.5, 80.5);
46465e39 261 fHistPtFraction[iEvtClass] = new TH2F(Form("fHistPtFraction%i", iEvtClass),
262 Form("fraction from embedded jet in reconstructed jet;fraction (event class %i);p_{T}^{emb}", iEvtClass),
3c6a60f7 263 52, 0., 1.04, 250, 0., 250.);
2a436fa1 264 }
265
46465e39 266
2a436fa1 267 fOutputList->Add(fHistEvtSelection);
46465e39 268 fOutputList->Add(fHistEvtClass);
269 fOutputList->Add(fHistCentrality);
3c6a60f7 270 fOutputList->Add(fHistNInputTracks);
271 fOutputList->Add(fHistCentVsTracks);
272 fOutputList->Add(fHistReactionPlane);
273 fOutputList->Add(fHistReactionPlaneWrtJets);
274
46465e39 275
2a436fa1 276 for (Int_t iBranch = 0; iBranch < fkNbranches; iBranch++) {
46465e39 277 fOutputList->Add(fHistPtJet[iBranch]);
278 fOutputList->Add(fHistEtaPhiJet[iBranch]);
279 fOutputList->Add(fHistEtaPhiJetCut[iBranch]);
2a436fa1 280 }
3c6a60f7 281
2a436fa1 282 for (Int_t iEvtClass =0; iEvtClass < fkEvtClasses; iEvtClass++) {
46465e39 283 fOutputList->Add(fHistDeltaEtaDeltaPhiJet[iEvtClass]);
284 fOutputList->Add(fHistDeltaEtaDeltaPhiJetCut[iEvtClass]);
3c6a60f7 285 //fOutputList->Add(fHistDeltaEtaDeltaPhiJetNOMatching[iEvtClass]);
46465e39 286 fOutputList->Add(fHistDeltaEtaEtaJet[iEvtClass]);
287 fOutputList->Add(fHistDeltaPtEtaJet[iEvtClass]);
288 fOutputList->Add(fHistPtFraction[iEvtClass]);
289
2a436fa1 290 fOutputList->Add(fHistPtResponse[iEvtClass]);
291 fOutputList->Add(fHistPtSmearing[iEvtClass]);
46465e39 292 fOutputList->Add(fHistDeltaR[iEvtClass]);
a11972e9 293 fOutputList->Add(fHistArea[iEvtClass]);
3c6a60f7 294 //fOutputList->Add(fHistJetPtArea[iEvtClass]);
46465e39 295 fOutputList->Add(fHistDeltaArea[iEvtClass]);
3c6a60f7 296 fOutputList->Add(fHistJetPtDeltaArea[iEvtClass]);
2a436fa1 297 }
46465e39 298
299 // =========== Switch on Sumw2 for all histos ===========
300 for (Int_t i=0; i<fOutputList->GetEntries(); ++i) {
301 TH1 *h1 = dynamic_cast<TH1*>(fOutputList->At(i));
302 if (h1){
303 h1->Sumw2();
304 continue;
305 }
306 }
307 TH1::AddDirectory(oldStatus);
308
309 PostData(1, fOutputList);
2a436fa1 310}
311
312void AliAnalysisTaskJetResponse::UserExec(Option_t *)
313{
314 // load events, apply event cuts, then compare leading jets
2a436fa1 315 if(!strlen(fJetBranchName[0].Data()) || !strlen(fJetBranchName[1].Data())){
316 AliError("Jet branch name not set.");
317 return;
318 }
319
320 fESD=dynamic_cast<AliESDEvent*>(InputEvent());
321 if (!fESD) {
322 AliError("ESD not available");
323 return;
324 }
325 fAOD = dynamic_cast<AliAODEvent*>(AODEvent());
326 if (!fAOD) {
327 AliError("AOD not available");
328 return;
329 }
330
46465e39 331 // -- event selection --
2a436fa1 332 fHistEvtSelection->Fill(1); // number of events before event selection
333
334 // physics selection
335 AliInputEventHandler* inputHandler = (AliInputEventHandler*)
336 ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
337 if(!(inputHandler->IsEventSelected() & fOfflineTrgMask)){
46465e39 338 if(fDebug) Printf(" Trigger Selection: event REJECTED ... ");
2a436fa1 339 fHistEvtSelection->Fill(2);
340 PostData(1, fOutputList);
341 return;
342 }
343
344 // vertex selection
345 AliAODVertex* primVtx = fAOD->GetPrimaryVertex();
346 Int_t nTracksPrim = primVtx->GetNContributors();
347 if ((nTracksPrim < fMinContribVtx) ||
348 (primVtx->GetZ() < fVtxZMin) ||
349 (primVtx->GetZ() > fVtxZMax) ){
46465e39 350 if(fDebug) Printf("%s:%d primary vertex z = %f: event REJECTED...",(char*)__FILE__,__LINE__,primVtx->GetZ());
2a436fa1 351 fHistEvtSelection->Fill(3);
352 PostData(1, fOutputList);
353 return;
354 }
355
356 // event class selection (from jet helper task)
357 Int_t eventClass = AliAnalysisHelperJetTasks::EventClass();
46465e39 358 if(fDebug) Printf("Event class %d", eventClass);
2a436fa1 359 if (eventClass < fEvtClassMin || eventClass > fEvtClassMax){
360 fHistEvtSelection->Fill(4);
361 PostData(1, fOutputList);
362 return;
363 }
364
365 // centrality selection
366 AliCentrality *cent = fESD->GetCentrality();
3c6a60f7 367 Float_t centValue = cent->GetCentralityPercentile("V0M");
46465e39 368 if(fDebug) printf("centrality: %f\n", centValue);
2a436fa1 369 if (centValue < fCentMin || centValue > fCentMax){
46465e39 370 fHistEvtSelection->Fill(4);
2a436fa1 371 PostData(1, fOutputList);
372 return;
373 }
374
3c6a60f7 375
376 // multiplicity due to input tracks
377 Int_t nInputTracks = 0;
378
379 TString jbname(fJetBranchName[1]);
380 for(Int_t i=1; i<=3; ++i){
381 if(jbname.Contains(Form("B%d",i))) jbname.ReplaceAll(Form("B%d",i),"B0");
382 }
383 // use only HI event
384 if(jbname.Contains("AODextraonly")) jbname.ReplaceAll("AODextraonly","AOD");
385 if(jbname.Contains("AODextra")) jbname.ReplaceAll("AODextra","AOD");
386
387 if(fDebug) Printf("Multiplicity from jet branch %s", jbname.Data());
388 TClonesArray *tmpAODjets = dynamic_cast<TClonesArray*>(fAOD->FindListObject(jbname.Data()));
389
74348808 390 if (tmpAODjets) {
391 for (Int_t iJet=0; iJet<tmpAODjets->GetEntriesFast(); iJet++){
3c6a60f7 392 AliAODJet *jet = dynamic_cast<AliAODJet*>((*tmpAODjets)[iJet]);
393 if(!jet) continue;
394 TRefArray *trackList = jet->GetRefTracks();
395 Int_t nTracks = trackList->GetEntriesFast();
396 nInputTracks += nTracks;
397 if(fDebug) Printf("#jet%d: %d tracks", iJet, nTracks);
74348808 398 }
3c6a60f7 399 }
400 if(fDebug) Printf("---> input tracks: %d", nInputTracks);
401
402 if (nInputTracks < fNInputTracksMin || (fNInputTracksMax > -1 && nInputTracks > fNInputTracksMax)){
403 fHistEvtSelection->Fill(5);
404 PostData(1, fOutputList);
405 return;
406 }
407
408 if(fEvtClassMode==1){ //multiplicity
409 fHistEvtClass->SetTitle("event classes (multiplicity)");
410 eventClass = fkEvtClasses-1 - (Int_t)(nInputTracks/250.);
411 if(eventClass<0) eventClass = 0;
412 }
413
414
415
2a436fa1 416 fHistEvtSelection->Fill(0); // accepted events
46465e39 417 fHistEvtClass->Fill(eventClass);
418 fHistCentrality->Fill(centValue);
3c6a60f7 419 fHistNInputTracks->Fill(nInputTracks);
420 fHistCentVsTracks->Fill(centValue,nInputTracks);
421 Double_t rp = AliAnalysisHelperJetTasks::ReactionPlane(kFALSE);
46465e39 422 // -- end event selection --
2a436fa1 423
3c6a60f7 424
425
46465e39 426 // fetch jets
2a436fa1 427 TClonesArray *aodJets[2];
a11972e9 428 aodJets[0] = dynamic_cast<TClonesArray*>(fAOD->FindListObject(fJetBranchName[0].Data())); // in general: embedded jet
429 aodJets[1] = dynamic_cast<TClonesArray*>(fAOD->FindListObject(fJetBranchName[1].Data())); // in general: embedded jet + UE
2a436fa1 430
431 for (Int_t iJetType = 0; iJetType < 2; iJetType++) {
432 fListJets[iJetType]->Clear();
46465e39 433 if (!aodJets[iJetType]) continue;
434
2a436fa1 435 for (Int_t iJet = 0; iJet < aodJets[iJetType]->GetEntriesFast(); iJet++) {
436 AliAODJet *jet = dynamic_cast<AliAODJet*>((*aodJets[iJetType])[iJet]);
46465e39 437 if (jet) fListJets[iJetType]->Add(jet);
2a436fa1 438 }
439 fListJets[iJetType]->Sort();
2a436fa1 440 }
3c6a60f7 441
46465e39 442 // jet matching
443 static TArrayI aMatchIndex(fListJets[0]->GetEntries());
444 static TArrayF aPtFraction(fListJets[0]->GetEntries());
445 if(aMatchIndex.GetSize()<fListJets[0]->GetEntries()) aMatchIndex.Set(fListJets[0]->GetEntries());
446 if(aPtFraction.GetSize()<fListJets[0]->GetEntries()) aPtFraction.Set(fListJets[0]->GetEntries());
447
448 // stores matched jets in 'aMatchIndex' and fraction of pT in 'aPtFraction'
449 AliAnalysisHelperJetTasks::GetJetMatching(fListJets[0], TMath::Min((Int_t)fNMatchJets,(Int_t)fListJets[0]->GetEntries()),
450 fListJets[1], TMath::Min((Int_t)fNMatchJets,(Int_t)fListJets[1]->GetEntries()),
451 aMatchIndex, aPtFraction, fDebug);
452
453 // loop over matched jets
454 Int_t ir = -1; // index of matched reconstruced jet
455 Float_t fraction = 0.;
456 AliAODJet *jet[2] = { 0x0, 0x0 };
457 Float_t jetEta[2] = { 0., 0. };
458 Float_t jetPhi[2] = { 0., 0. };
459 Float_t jetPt[2] = { 0., 0. };
460 Float_t jetArea[2] = { 0., 0. };
461
462 for(Int_t ig=0; ig<fListJets[0]->GetEntries(); ++ig){
463 ir = aMatchIndex[ig];
464 if(ir<0) continue;
465 fraction = aPtFraction[ig];
466
467 // fetch jets
468 jet[0] = (AliAODJet*)(fListJets[0]->At(ig));
469 jet[1] = (AliAODJet*)(fListJets[1]->At(ir));
470 if(!jet[0] || !jet[1]) continue;
471
472 for(Int_t i=0; i<fkNbranches; ++i){
473 jetEta[i] = jet[i]->Eta();
474 jetPhi[i] = jet[i]->Phi();
475 jetPt[i] = jet[i]->Pt();
476 jetArea[i] = jet[i]->EffectiveAreaCharged();
477 }
3c6a60f7 478
479 // reaction plane;
480 if(fReactionPlaneBin>=0){
481 Int_t rpBin = AliAnalysisHelperJetTasks::GetPhiBin(TVector2::Phi_mpi_pi(rp-jetPhi[1]), 3);
482 if(rpBin!=fReactionPlaneBin) continue;
483 }
484 fHistReactionPlane->Fill(rp);
485 fHistReactionPlaneWrtJets->Fill(TVector2::Phi_mpi_pi(rp-jetPhi[1]));
486
487
488
46465e39 489 if(eventClass > -1 && eventClass < fkEvtClasses){
490 fHistPtFraction[eventClass] -> Fill(fraction, jetPt[0]);
491 }
492
493 if(fraction<fJetPtFractionMin) continue;
494
495 // calculate parameters of associated jets
496 Float_t deltaPt = jetPt[1]-jetPt[0];
497 Float_t deltaEta = jetEta[1]-jetEta[0];
498 Float_t deltaPhi = TVector2::Phi_mpi_pi(jetPhi[1]-jetPhi[0]);
499 Float_t deltaR = TMath::Sqrt(deltaEta*deltaEta+deltaPhi*deltaPhi);
500 Float_t deltaArea = jetArea[1]-jetArea[0];
501
502 // fill histograms before acceptance cut
503 for(Int_t i=0; i<fkNbranches; ++i){
504 fHistEtaPhiJet[i] -> Fill(jetEta[i], jetPhi[i]);
505 }
506 if(eventClass > -1 && eventClass < fkEvtClasses){
507 fHistDeltaEtaDeltaPhiJet[eventClass] -> Fill(deltaEta, deltaPhi);
508 }
509
510 // jet acceptance + minimum pT check
511 if(jetEta[0]>fJetEtaMax || jetEta[0]<fJetEtaMin ||
512 jetEta[1]>fJetEtaMax || jetEta[1]<fJetEtaMin){
513
514 if(fDebug){
515 Printf("Jet not in eta acceptance.");
516 Printf("[0]: jet %d eta %.2f", ig, jetEta[0]);
517 Printf("[1]: jet %d eta %.2f", ir, jetEta[1]);
518 }
519 continue;
520 }
521 if(jetPt[1] < fJetPtMin){
522 if(fDebug) Printf("Jet %d (pT %.1f GeV/c) has less than required pT.", ir, jetPt[1]);
523 continue;
524 }
525
526
527
528 // fill histograms
529 for(Int_t i=0; i<fkNbranches; ++i){
530 fHistPtJet[i] -> Fill(jetPt[i]);
531 fHistEtaPhiJetCut[i] -> Fill(jetEta[i], jetPhi[i]);
532 }
533
46465e39 534 if(eventClass > -1 && eventClass < fkEvtClasses){
535 fHistDeltaEtaDeltaPhiJetCut[eventClass] -> Fill(deltaEta, deltaPhi);
536
537 fHistPtResponse[eventClass] -> Fill(jetPt[1], jetPt[0]);
538 fHistPtSmearing[eventClass] -> Fill(deltaPt, jetPt[0]);
539
540 fHistDeltaPtEtaJet[eventClass] -> Fill(jetEta[0], deltaPt);
541 fHistDeltaEtaEtaJet[eventClass] -> Fill(jetEta[0], deltaEta);
542
543 fHistDeltaR[eventClass] -> Fill(deltaPt, deltaR);
3c6a60f7 544 fHistArea[eventClass] -> Fill(jetPt[0], deltaPt, jetArea[1]);
545 //fHistJetPtArea[eventClass] -> Fill(jetPt[0], jetArea[1]);
46465e39 546 fHistDeltaArea[eventClass] -> Fill(deltaPt, deltaArea);
3c6a60f7 547 fHistJetPtDeltaArea[eventClass] -> Fill(jetPt[0], deltaArea);
46465e39 548
549 }
2a436fa1 550 }
551
552 PostData(1, fOutputList);
553}
554
555void AliAnalysisTaskJetResponse::Terminate(const Option_t *)
556{
557 // Draw result to the screen
558 // Called once at the end of the query
559
560 if (!GetOutputData(1))
561 return;
562}
46465e39 563