]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGJE/EMCALJetTasks/AliJetResponseMaker.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
update from Salvatore
[u/mrichter/AliRoot.git] / PWGJE / EMCALJetTasks / AliJetResponseMaker.cxx
CommitLineData
7549c451 1// $Id$
2949a2ec 2//
3// Emcal jet response matrix maker task.
4//
5// Author: S. Aiola
6
7#include "AliJetResponseMaker.h"
8
9#include <TChain.h>
10#include <TClonesArray.h>
11#include <TH1F.h>
12#include <TH2F.h>
13#include <TList.h>
14#include <TLorentzVector.h>
15
16#include "AliAnalysisManager.h"
17#include "AliCentrality.h"
18#include "AliFJWrapper.h"
19#include "AliVCluster.h"
20#include "AliVTrack.h"
21#include "AliEmcalJet.h"
1ee1b5b8 22#include "AliGenPythiaEventHeader.h"
2949a2ec 23#include "AliMCEvent.h"
24
25ClassImp(AliJetResponseMaker)
26
27//________________________________________________________________________
28AliJetResponseMaker::AliJetResponseMaker() :
e44e8726 29 AliAnalysisTaskEmcalJet("AliJetResponseMaker", kTRUE),
7549c451 30 fMCTracksName("MCParticles"),
31 fMCJetsName("MCJets"),
99c67c1b 32 fMaxDistance(0.25),
1ee1b5b8 33 fDoWeighting(kFALSE),
4643d2e8 34 fEventWeightHist(kFALSE),
91bee8bc 35 fMCFiducial(kFALSE),
36 fMCminEta(0),
37 fMCmaxEta(0),
38 fMCminPhi(0),
39 fMCmaxPhi(0),
4643d2e8 40 fSelectPtHardBin(-999),
aa4d701c 41 fDoMatching(kTRUE),
1ee1b5b8 42 fPythiaHeader(0),
43 fEventWeight(0),
44 fPtHardBin(0),
45 fNTrials(0),
7549c451 46 fMCTracks(0),
2949a2ec 47 fMCJets(0),
1ee1b5b8 48 fHistNTrials(0),
1ee1b5b8 49 fHistEvents(0),
2949a2ec 50 fHistMCJetPhiEta(0),
51 fHistMCJetsPt(0),
65bb5510 52 fHistMCJetPhiEtaFiducial(0),
53 fHistMCJetsPtFiducial(0),
2949a2ec 54 fHistMCJetsNEFvsPt(0),
55 fHistMCJetsZvsPt(0),
56 fHistJetPhiEta(0),
57 fHistJetsPt(0),
2949a2ec 58 fHistJetsNEFvsPt(0),
1b76c28f 59 fHistJetsZvsPt(0),
60 fHistClosestDistance(0),
61 fHistClosestDeltaPhi(0),
62 fHistClosestDeltaEta(0),
63 fHistClosestDeltaPt(0),
99c67c1b 64 fHistNonMatchedMCJetPt(0),
65 fHistNonMatchedJetPt(0),
2bddb6ae 66 fHistPartvsDetecPt(0),
67 fHistMissedMCJets(0)
2949a2ec 68{
69 // Default constructor.
4643d2e8 70
71 for (Int_t i = 0; i < 11; i++) {
72 fHistEventWeight[i] = 0;
73 }
a487deae 74
75 SetMakeGeneralHistograms(kTRUE);
2949a2ec 76}
77
78//________________________________________________________________________
79AliJetResponseMaker::AliJetResponseMaker(const char *name) :
e44e8726 80 AliAnalysisTaskEmcalJet(name, kTRUE),
7549c451 81 fMCTracksName("MCParticles"),
82 fMCJetsName("MCJets"),
99c67c1b 83 fMaxDistance(0.25),
1ee1b5b8 84 fDoWeighting(kFALSE),
4643d2e8 85 fEventWeightHist(kFALSE),
91bee8bc 86 fMCFiducial(kFALSE),
87 fMCminEta(0),
88 fMCmaxEta(0),
89 fMCminPhi(0),
90 fMCmaxPhi(0),
4643d2e8 91 fSelectPtHardBin(-999),
aa4d701c 92 fDoMatching(kTRUE),
1ee1b5b8 93 fPythiaHeader(0),
94 fEventWeight(0),
95 fPtHardBin(0),
96 fNTrials(0),
7549c451 97 fMCTracks(0),
2949a2ec 98 fMCJets(0),
1ee1b5b8 99 fHistNTrials(0),
1ee1b5b8 100 fHistEvents(0),
2949a2ec 101 fHistMCJetPhiEta(0),
102 fHistMCJetsPt(0),
65bb5510 103 fHistMCJetPhiEtaFiducial(0),
104 fHistMCJetsPtFiducial(0),
2949a2ec 105 fHistMCJetsNEFvsPt(0),
106 fHistMCJetsZvsPt(0),
107 fHistJetPhiEta(0),
108 fHistJetsPt(0),
2949a2ec 109 fHistJetsNEFvsPt(0),
1b76c28f 110 fHistJetsZvsPt(0),
111 fHistClosestDistance(0),
112 fHistClosestDeltaPhi(0),
113 fHistClosestDeltaEta(0),
114 fHistClosestDeltaPt(0),
99c67c1b 115 fHistNonMatchedMCJetPt(0),
116 fHistNonMatchedJetPt(0),
2bddb6ae 117 fHistPartvsDetecPt(0),
118 fHistMissedMCJets(0)
2949a2ec 119{
120 // Standard constructor.
4643d2e8 121
122 for (Int_t i = 0; i < 11; i++) {
123 fHistEventWeight[i] = 0;
124 }
a487deae 125
126 SetMakeGeneralHistograms(kTRUE);
2949a2ec 127}
128
129//________________________________________________________________________
130AliJetResponseMaker::~AliJetResponseMaker()
131{
132 // Destructor
133}
134
135//________________________________________________________________________
136void AliJetResponseMaker::UserCreateOutputObjects()
137{
138 // Create user objects.
139
a487deae 140 AliAnalysisTaskEmcalJet::UserCreateOutputObjects();
2949a2ec 141
1ee1b5b8 142 const Int_t ptHardLo[11] = { 0, 5,11,21,36,57, 84,117,152,191,234};
143 const Int_t ptHardHi[11] = { 5,11,21,36,57,84,117,152,191,234,1000000};
144
145 fHistNTrials = new TH1F("fHistNTrials", "fHistNTrials", 11, 0, 11);
146 fHistNTrials->GetXaxis()->SetTitle("p_{T} hard bin");
147 fHistNTrials->GetYaxis()->SetTitle("trials");
148 fOutput->Add(fHistNTrials);
149
1ee1b5b8 150 fHistEvents = new TH1F("fHistEvents", "fHistEvents", 11, 0, 11);
151 fHistEvents->GetXaxis()->SetTitle("p_{T} hard bin");
152 fHistEvents->GetYaxis()->SetTitle("total events");
153 fOutput->Add(fHistEvents);
154
4643d2e8 155 if (fEventWeightHist) {
156 for (Int_t i = 0; i < 11; i++) {
157 TString name(Form("fHistEventWeight_%d", i+1));
158 fHistEventWeight[i] = new TH1F(name, name, 10, 0, 10);
159 fOutput->Add(fHistEventWeight[i]);
160 }
161 }
162
1ee1b5b8 163 for (Int_t i = 1; i < 12; i++) {
164 fHistNTrials->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
1ee1b5b8 165 fHistEvents->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
166 }
167
7549c451 168 fHistJetPhiEta = new TH2F("fHistJetPhiEta", "fHistJetPhiEta", 20, -2, 2, 32, 0, 6.4);
169 fHistJetPhiEta->GetXaxis()->SetTitle("#eta");
170 fHistJetPhiEta->GetYaxis()->SetTitle("#phi");
171 fOutput->Add(fHistJetPhiEta);
172
6fd5039f 173 fHistJetsPt = new TH1F("fHistJetsPt", "fHistJetsPt", fNbins, fMinBinPt, fMaxBinPt);
2949a2ec 174 fHistJetsPt->GetXaxis()->SetTitle("p_{T} [GeV/c]");
175 fHistJetsPt->GetYaxis()->SetTitle("counts");
176 fOutput->Add(fHistJetsPt);
177
6fd5039f 178 fHistJetsZvsPt = new TH2F("fHistJetsZvsPt", "fHistJetsZvsPt", fNbins, 0, 1.2, fNbins, fMinBinPt, fMaxBinPt);
2949a2ec 179 fHistJetsZvsPt->GetXaxis()->SetTitle("Z");
180 fHistJetsZvsPt->GetYaxis()->SetTitle("p_{T} [GeV/c]");
181 fOutput->Add(fHistJetsZvsPt);
182
cd588263 183 fHistJetsNEFvsPt = new TH2F("fHistJetsNEFvsPt", "fHistJetsNEFvsPt", fNbins, 0, 1.2, fNbins, fMinBinPt, fMaxBinPt);
184 fHistJetsNEFvsPt->GetXaxis()->SetTitle("NEF");
185 fHistJetsNEFvsPt->GetYaxis()->SetTitle("p_{T} [GeV/c]");
186 fOutput->Add(fHistJetsNEFvsPt);
7549c451 187
188 fHistMCJetPhiEta = new TH2F("fHistMCJetPhiEta", "fHistMCJetPhiEta", 20, -2, 2, 32, 0, 6.4);
189 fHistMCJetPhiEta->GetXaxis()->SetTitle("#eta");
190 fHistMCJetPhiEta->GetYaxis()->SetTitle("#phi");
191 fOutput->Add(fHistMCJetPhiEta);
192
6fd5039f 193 fHistMCJetsPt = new TH1F("fHistMCJetsPt", "fHistMCJetsPt", fNbins, fMinBinPt, fMaxBinPt);
7549c451 194 fHistMCJetsPt->GetXaxis()->SetTitle("p_{T} [GeV/c]");
195 fHistMCJetsPt->GetYaxis()->SetTitle("counts");
196 fOutput->Add(fHistMCJetsPt);
87203fd3 197
198 fHistMCJetPhiEtaFiducial = new TH2F("fHistMCJetPhiEtaFiducial", "fHistMCJetPhiEtaFiducial", 20, -2, 2, 32, 0, 6.4);
199 fHistMCJetPhiEtaFiducial->GetXaxis()->SetTitle("#eta");
200 fHistMCJetPhiEtaFiducial->GetYaxis()->SetTitle("#phi");
201 fOutput->Add(fHistMCJetPhiEtaFiducial);
202
203 fHistMCJetsPtFiducial = new TH1F("fHistMCJetsPtFiducial", "fHistMCJetsPtFiducial", fNbins, fMinBinPt, fMaxBinPt);
204 fHistMCJetsPtFiducial->GetXaxis()->SetTitle("p_{T} [GeV/c]");
205 fHistMCJetsPtFiducial->GetYaxis()->SetTitle("counts");
206 fOutput->Add(fHistMCJetsPtFiducial);
7549c451 207
6fd5039f 208 fHistMCJetsZvsPt = new TH2F("fHistMCJetsZvsPt", "fHistMCJetsZvsPt", fNbins, 0, 1.2, fNbins, fMinBinPt, fMaxBinPt);
7549c451 209 fHistMCJetsZvsPt->GetXaxis()->SetTitle("Z");
210 fHistMCJetsZvsPt->GetYaxis()->SetTitle("p_{T} [GeV/c]");
211 fOutput->Add(fHistMCJetsZvsPt);
cd588263 212
213 fHistMCJetsNEFvsPt = new TH2F("fHistMCJetsNEFvsPt", "fHistMCJetsNEFvsPt", fNbins, 0, 1.2, fNbins, fMinBinPt, fMaxBinPt);
214 fHistMCJetsNEFvsPt->GetXaxis()->SetTitle("NEF");
215 fHistMCJetsNEFvsPt->GetYaxis()->SetTitle("p_{T} [GeV/c]");
216 fOutput->Add(fHistMCJetsNEFvsPt);
7549c451 217
1b76c28f 218 fHistClosestDistance = new TH1F("fHistClosestDistance", "fHistClosestDistance", 50, 0, fMaxDistance * 1.2);
219 fHistClosestDistance->GetXaxis()->SetTitle("d");
220 fHistClosestDistance->GetYaxis()->SetTitle("counts");
221 fOutput->Add(fHistClosestDistance);
222
223 fHistClosestDeltaPhi = new TH1F("fHistClosestDeltaPhi", "fHistClosestDeltaPhi", 64, -1.6, 4.8);
224 fHistClosestDeltaPhi->GetXaxis()->SetTitle("#Delta#phi");
225 fHistClosestDeltaPhi->GetYaxis()->SetTitle("counts");
226 fOutput->Add(fHistClosestDeltaPhi);
227
a487deae 228 fHistClosestDeltaEta = new TH1F("fHistClosestDeltaEta", "fHistClosestDeltaEta", TMath::CeilNint(fJetMaxEta - fJetMinEta) * 20, fJetMinEta * 2, fJetMaxEta * 2);
1b76c28f 229 fHistClosestDeltaEta->GetXaxis()->SetTitle("#Delta#eta");
230 fHistClosestDeltaEta->GetYaxis()->SetTitle("counts");
231 fOutput->Add(fHistClosestDeltaEta);
232
6fd5039f 233 fHistClosestDeltaPt = new TH1F("fHistClosestDeltaPt", "fHistClosestDeltaPt", fNbins, -fMaxBinPt / 2, fMaxBinPt / 2);
1b76c28f 234 fHistClosestDeltaPt->GetXaxis()->SetTitle("#Delta p_{T}");
235 fHistClosestDeltaPt->GetYaxis()->SetTitle("counts");
236 fOutput->Add(fHistClosestDeltaPt);
237
99c67c1b 238 fHistNonMatchedMCJetPt = new TH1F("fHistNonMatchedMCJetPt", "fHistNonMatchedMCJetPt", fNbins, fMinBinPt, fMaxBinPt);
239 fHistNonMatchedMCJetPt->GetXaxis()->SetTitle("p_{T} [GeV/c]");
240 fHistNonMatchedMCJetPt->GetYaxis()->SetTitle("counts");
241 fOutput->Add(fHistNonMatchedMCJetPt);
242
243 fHistNonMatchedJetPt = new TH1F("fHistNonMatchedJetPt", "fHistNonMatchedJetPt", fNbins, fMinBinPt, fMaxBinPt);
244 fHistNonMatchedJetPt->GetXaxis()->SetTitle("p_{T} [GeV/c]");
245 fHistNonMatchedJetPt->GetYaxis()->SetTitle("counts");
246 fOutput->Add(fHistNonMatchedJetPt);
247
6fd5039f 248 fHistPartvsDetecPt = new TH2F("fHistPartvsDetecPt", "fHistPartvsDetecPt", fNbins, fMinBinPt, fMaxBinPt, fNbins, fMinBinPt, fMaxBinPt);
99c67c1b 249 fHistPartvsDetecPt->GetXaxis()->SetTitle("p_{T}^{rec}");
250 fHistPartvsDetecPt->GetYaxis()->SetTitle("p_{T}^{gen}");
1b76c28f 251 fOutput->Add(fHistPartvsDetecPt);
252
2bddb6ae 253 fHistMissedMCJets = new TH1F("fHistMissedMCJets", "fHistMissedMCJets", fNbins, fMinBinPt, fMaxBinPt);
254 fHistMissedMCJets->GetXaxis()->SetTitle("p_{T} [GeV/c]");
255 fHistMissedMCJets->GetYaxis()->SetTitle("counts");
256 fOutput->Add(fHistMissedMCJets);
257
99c67c1b 258 PostData(1, fOutput); // Post data for ALL output slots > 0 here, to get at least an empty histogram
2949a2ec 259}
260
ceefbfbc 261//________________________________________________________________________
a487deae 262Bool_t AliJetResponseMaker::AcceptJet(AliEmcalJet *jet) const
ceefbfbc 263{
264 // Return true if jet is accepted.
265
266 if (jet->Pt() <= fJetPtCut)
267 return kFALSE;
268 if (jet->Area() <= fJetAreaCut)
269 return kFALSE;
65bb5510 270
271 return kTRUE;
ceefbfbc 272}
273
b3ccdfe2 274//________________________________________________________________________
275void AliJetResponseMaker::ExecOnce()
276{
277 // Execute once.
278
279 if (!fMCJetsName.IsNull() && !fMCJets) {
280 fMCJets = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fMCJetsName));
281 if (!fMCJets) {
282 AliError(Form("%s: Could not retrieve mc jets %s!", GetName(), fMCJetsName.Data()));
283 return;
284 }
285 else if (!fMCJets->GetClass()->GetBaseClass("AliEmcalJet")) {
286 AliError(Form("%s: Collection %s does not contain AliEmcalJet objects!", GetName(), fMCJetsName.Data()));
287 fMCJets = 0;
288 return;
289 }
290 }
291
292 if (!fMCTracksName.IsNull() && !fMCTracks) {
293 fMCTracks = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fMCTracksName));
294 if (!fMCTracks) {
295 AliError(Form("%s: Could not retrieve mc tracks %s!", GetName(), fMCTracksName.Data()));
296 return;
297 }
298 else {
299 TClass *cl = fMCTracks->GetClass();
300 if (!cl->GetBaseClass("AliVParticle") && !cl->GetBaseClass("AliEmcalParticle")) {
301 AliError(Form("%s: Collection %s does not contain AliVParticle nor AliEmcalParticle objects!", GetName(), fMCTracksName.Data()));
302 fMCTracks = 0;
303 return;
304 }
305 }
306 }
307
308 AliAnalysisTaskEmcalJet::ExecOnce();
91bee8bc 309
310 if (fMCFiducial) {
a487deae 311 fMCminEta = fJetMinEta;
312 fMCmaxEta = fJetMaxEta;
313 fMCminPhi = fJetMinPhi;
314 fMCmaxPhi = fJetMaxPhi;
91bee8bc 315 }
316 else {
a487deae 317 fMCminEta = fJetMinEta - fJetRadius;
318 fMCmaxEta = fJetMaxEta + fJetRadius;
319 fMCminPhi = fJetMinPhi - fJetRadius;
320 fMCmaxPhi = fJetMaxPhi + fJetRadius;
91bee8bc 321 }
b3ccdfe2 322}
323
4643d2e8 324//________________________________________________________________________
325Bool_t AliJetResponseMaker::IsEventSelected()
326{
327 // Check if event is selected
328
329 if (fSelectPtHardBin != -999 && fSelectPtHardBin != fPtHardBin)
330 return kFALSE;
331
332 return AliAnalysisTaskEmcalJet::IsEventSelected();
333}
334
b3ccdfe2 335//________________________________________________________________________
336Bool_t AliJetResponseMaker::RetrieveEventObjects()
337{
338 // Retrieve event objects.
339
340 if (!AliAnalysisTaskEmcalJet::RetrieveEventObjects())
341 return kFALSE;
342
343 fPythiaHeader = dynamic_cast<AliGenPythiaEventHeader*>(MCEvent()->GenEventHeader());
344
345 if (!fPythiaHeader)
346 return kFALSE;
347
348 if (fDoWeighting)
349 fEventWeight = fPythiaHeader->EventWeight();
350 else
351 fEventWeight = 1;
352
353 const Int_t ptHardLo[11] = { 0, 5,11,21,36,57, 84,117,152,191,234};
354 const Int_t ptHardHi[11] = { 5,11,21,36,57,84,117,152,191,234,1000000};
355
356 Double_t pthard = fPythiaHeader->GetPtHard();
91bee8bc 357
358 for (fPtHardBin = 0; fPtHardBin <= 11; fPtHardBin++) {
b3ccdfe2 359 if (pthard >= ptHardLo[fPtHardBin] && pthard < ptHardHi[fPtHardBin])
360 break;
361 }
b3ccdfe2 362
363 fNTrials = fPythiaHeader->Trials();
364
365 return kTRUE;
366}
367
368//________________________________________________________________________
369Bool_t AliJetResponseMaker::Run()
370{
371 // Find the closest jets
372
aa4d701c 373 if (!fDoMatching)
374 return kTRUE;
375
b3ccdfe2 376 DoJetLoop(fJets, fMCJets, kFALSE);
377 DoJetLoop(fMCJets, fJets, kTRUE);
378
aa4d701c 379 const Int_t nMCJets = fMCJets->GetEntriesFast();
380
381 for (Int_t i = 0; i < nMCJets; i++) {
382
383 AliEmcalJet* jet = static_cast<AliEmcalJet*>(fMCJets->At(i));
384
385 if (!jet) {
386 AliError(Form("Could not receive jet %d", i));
387 continue;
388 }
389
390 if (!AcceptJet(jet))
391 continue;
392
393 if (jet->Eta() < fMCminEta || jet->Eta() > fMCmaxEta || jet->Phi() < fMCminPhi || jet->Phi() > fMCmaxPhi)
394 continue;
395
396 if (jet->Pt() > fMaxBinPt)
397 continue;
398
399 if (jet->ClosestJet() && jet->ClosestJet()->ClosestJet() == jet &&
400 jet->ClosestJetDistance() < fMaxDistance) { // Matched jet found
401 jet->SetMatchedToClosest();
402 jet->ClosestJet()->SetMatchedToClosest();
403 }
404 }
405
b3ccdfe2 406 return kTRUE;
407}
408
2949a2ec 409//________________________________________________________________________
44476be7 410void AliJetResponseMaker::DoJetLoop(TClonesArray *jets1, TClonesArray *jets2, Bool_t mc)
2949a2ec 411{
44476be7 412 // Do the jet loop.
1ee1b5b8 413
44476be7 414 Int_t nJets1 = jets1->GetEntriesFast();
415 Int_t nJets2 = jets2->GetEntriesFast();
1ee1b5b8 416
44476be7 417 for (Int_t i = 0; i < nJets1; i++) {
1ee1b5b8 418
44476be7 419 AliEmcalJet* jet1 = static_cast<AliEmcalJet*>(jets1->At(i));
2949a2ec 420
44476be7 421 if (!jet1) {
422 AliError(Form("Could not receive jet %d", i));
423 continue;
424 }
425
ceefbfbc 426 if (!AcceptJet(jet1))
44476be7 427 continue;
428
ceefbfbc 429 if (!mc) {
a487deae 430 if (jet1->Eta() < fJetMinEta || jet1->Eta() > fJetMaxEta || jet1->Phi() < fJetMinPhi || jet1->Phi() > fJetMaxPhi)
ceefbfbc 431 continue;
432 }
433 else {
91bee8bc 434 if (jet1->Eta() < fMCminEta || jet1->Eta() > fMCmaxEta || jet1->Phi() < fMCminPhi || jet1->Phi() > fMCmaxPhi)
ceefbfbc 435 continue;
436 }
437
44476be7 438 for (Int_t j = 0; j < nJets2; j++) {
439
440 AliEmcalJet* jet2 = static_cast<AliEmcalJet*>(jets2->At(j));
441
442 if (!jet2) {
443 AliError(Form("Could not receive jet %d", j));
444 continue;
445 }
446
ceefbfbc 447 if (!AcceptJet(jet2))
44476be7 448 continue;
ceefbfbc 449
450 if (mc) {
a487deae 451 if (jet2->Eta() < fJetMinEta || jet2->Eta() > fJetMaxEta || jet2->Phi() < fJetMinPhi || jet2->Phi() > fJetMaxPhi)
ceefbfbc 452 continue;
453 }
454 else {
91bee8bc 455 if (jet1->Eta() < fMCminEta || jet1->Eta() > fMCmaxEta || jet1->Phi() < fMCminPhi || jet1->Phi() > fMCmaxPhi)
ceefbfbc 456 continue;
457 }
44476be7 458
459 Double_t deta = jet2->Eta() - jet1->Eta();
460 Double_t dphi = jet2->Phi() - jet1->Phi();
461 Double_t d = TMath::Sqrt(deta * deta + dphi * dphi);
462
463 if (d < jet1->ClosestJetDistance()) {
464 jet1->SetSecondClosestJet(jet1->ClosestJet(), jet1->ClosestJetDistance());
465 jet1->SetClosestJet(jet2, d);
466 }
467 else if (d < jet1->SecondClosestJetDistance()) {
468 jet1->SetSecondClosestJet(jet2, d);
469 }
470 }
471 }
1ee1b5b8 472}
473
1ee1b5b8 474//________________________________________________________________________
475Bool_t AliJetResponseMaker::FillHistograms()
476{
477 // Fill histograms.
478
4643d2e8 479 fHistEvents->SetBinContent(fPtHardBin + 1, fHistEvents->GetBinContent(fPtHardBin + 1) + 1);
480 fHistNTrials->SetBinContent(fPtHardBin + 1, fHistNTrials->GetBinContent(fPtHardBin + 1) + fNTrials);
481 if (fEventWeightHist)
482 fHistEventWeight[fPtHardBin]->Fill(fPythiaHeader->EventWeight());
4643d2e8 483
2bddb6ae 484 const Int_t nMCJets = fMCJets->GetEntriesFast();
2949a2ec 485
2bddb6ae 486 for (Int_t i = 0; i < nMCJets; i++) {
2949a2ec 487
e44e8726 488 AliEmcalJet* jet = static_cast<AliEmcalJet*>(fMCJets->At(i));
2949a2ec 489
490 if (!jet) {
1b76c28f 491 AliError(Form("Could not receive jet %d", i));
2949a2ec 492 continue;
493 }
494
ceefbfbc 495 if (!AcceptJet(jet))
496 continue;
497
91bee8bc 498 if (jet->Eta() < fMCminEta || jet->Eta() > fMCmaxEta || jet->Phi() < fMCminPhi || jet->Phi() > fMCmaxPhi)
2bddb6ae 499 continue;
500
501 if (jet->Pt() > fMaxBinPt)
2949a2ec 502 continue;
503
aa4d701c 504 if (jet->MatchedJet()) {
505
506 if (!AcceptBiasJet(jet->MatchedJet()) ||
507 jet->MatchedJet()->MaxTrackPt() > fMaxTrackPt || jet->MatchedJet()->MaxClusterPt() > fMaxClusterPt ||
508 jet->MatchedJet()->Pt() > fMaxBinPt) {
1ee1b5b8 509 fHistMissedMCJets->Fill(jet->Pt(), fEventWeight);
2bddb6ae 510 }
511 else {
1ee1b5b8 512 fHistClosestDistance->Fill(jet->ClosestJetDistance(), fEventWeight);
513
2bddb6ae 514 Double_t deta = jet->MatchedJet()->Eta() - jet->Eta();
1ee1b5b8 515 fHistClosestDeltaEta->Fill(deta, fEventWeight);
516
2bddb6ae 517 Double_t dphi = jet->MatchedJet()->Phi() - jet->Phi();
1ee1b5b8 518 fHistClosestDeltaPhi->Fill(dphi, fEventWeight);
519
2bddb6ae 520 Double_t dpt = jet->MatchedJet()->Pt() - jet->Pt();
1ee1b5b8 521 fHistClosestDeltaPt->Fill(dpt, fEventWeight);
522
523 fHistPartvsDetecPt->Fill(jet->MatchedJet()->Pt(), jet->Pt(), fEventWeight);
2bddb6ae 524 }
2949a2ec 525 }
99c67c1b 526 else {
1ee1b5b8 527 fHistNonMatchedMCJetPt->Fill(jet->Pt(), fEventWeight);
528 fHistMissedMCJets->Fill(jet->Pt(), fEventWeight);
99c67c1b 529 }
2949a2ec 530
1ee1b5b8 531 fHistMCJetsPt->Fill(jet->Pt(), fEventWeight);
1ee1b5b8 532 fHistMCJetPhiEta->Fill(jet->Eta(), jet->Phi(), fEventWeight);
2949a2ec 533
cd588263 534 fHistMCJetsNEFvsPt->Fill(jet->NEF(), jet->Pt(), fEventWeight);
535
2949a2ec 536 for (Int_t it = 0; it < jet->GetNumberOfTracks(); it++) {
2bddb6ae 537 AliVParticle *track = jet->TrackAt(it, fMCTracks);
2949a2ec 538 if (track)
1ee1b5b8 539 fHistMCJetsZvsPt->Fill(track->Pt() / jet->Pt(), jet->Pt(), fEventWeight);
1b76c28f 540 }
ceefbfbc 541
542 if (!AcceptBiasJet(jet))
543 continue;
a487deae 544 if (jet->Eta() < fJetMinEta || jet->Eta() > fJetMaxEta || jet->Phi() < fJetMinPhi || jet->Phi() > fJetMaxPhi)
ceefbfbc 545 continue;
546
547 fHistMCJetsPtFiducial->Fill(jet->Pt(), fEventWeight);
548 fHistMCJetPhiEtaFiducial->Fill(jet->Eta(), jet->Phi(), fEventWeight);
1b76c28f 549 }
550
2bddb6ae 551 const Int_t nJets = fJets->GetEntriesFast();
1b76c28f 552
2bddb6ae 553 for (Int_t i = 0; i < nJets; i++) {
1b76c28f 554
e44e8726 555 AliEmcalJet* jet = static_cast<AliEmcalJet*>(fJets->At(i));
1b76c28f 556
557 if (!jet) {
558 AliError(Form("Could not receive mc jet %d", i));
2949a2ec 559 continue;
1b76c28f 560 }
99c67c1b 561
2bddb6ae 562 if (!AcceptJet(jet))
1b76c28f 563 continue;
91bee8bc 564 if (!AcceptBiasJet(jet))
565 continue;
566 if (jet->MaxTrackPt() > fMaxTrackPt || jet->MaxClusterPt() > fMaxClusterPt)
567 continue;
a487deae 568 if (jet->Eta() < fJetMinEta || jet->Eta() > fJetMaxEta || jet->Phi() < fJetMinPhi || jet->Phi() > fJetMaxPhi)
91bee8bc 569 continue;
2949a2ec 570
99c67c1b 571 if (!jet->MatchedJet()) {
1ee1b5b8 572 fHistNonMatchedJetPt->Fill(jet->Pt(), fEventWeight);
99c67c1b 573 }
1b76c28f 574
1ee1b5b8 575 fHistJetsPt->Fill(jet->Pt(), fEventWeight);
1b76c28f 576
1ee1b5b8 577 fHistJetPhiEta->Fill(jet->Eta(), jet->Phi(), fEventWeight);
1b76c28f 578
cd588263 579 fHistJetsNEFvsPt->Fill(jet->NEF(), jet->Pt(), fEventWeight);
1b76c28f 580
581 for (Int_t it = 0; it < jet->GetNumberOfTracks(); it++) {
2bddb6ae 582 AliVParticle *track = jet->TrackAt(it, fTracks);
1b76c28f 583 if (track)
1ee1b5b8 584 fHistJetsZvsPt->Fill(track->Pt() / jet->Pt(), jet->Pt(), fEventWeight);
2bddb6ae 585 }
586
587 for (Int_t ic = 0; ic < jet->GetNumberOfClusters(); ic++) {
588 AliVCluster *cluster = jet->ClusterAt(ic, fCaloClusters);
589 if (cluster) {
590 TLorentzVector nP;
591 cluster->GetMomentum(nP, fVertex);
1ee1b5b8 592 fHistJetsZvsPt->Fill(nP.Pt() / jet->Pt(), jet->Pt(), fEventWeight);
2bddb6ae 593 }
1b76c28f 594 }
595 }
6fd5039f 596
597 return kTRUE;
1b76c28f 598}