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