]> 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
Changes 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
2949a2ec 9#include <TClonesArray.h>
10#include <TH1F.h>
11#include <TH2F.h>
b16bb001 12#include <TH3F.h>
2949a2ec 13#include <TLorentzVector.h>
ca5c29fa 14#include <TSystem.h>
15#include <TFile.h>
16#include <TKey.h>
17#include <TProfile.h>
2949a2ec 18
ca5c29fa 19#include "AliAnalysisManager.h"
2949a2ec 20#include "AliVCluster.h"
21#include "AliVTrack.h"
22#include "AliEmcalJet.h"
1ee1b5b8 23#include "AliGenPythiaEventHeader.h"
2949a2ec 24#include "AliMCEvent.h"
787a3c4f 25#include "AliLog.h"
cd6431de 26#include "AliRhoParameter.h"
5be3857d 27#include "AliNamedArrayI.h"
2949a2ec 28
29ClassImp(AliJetResponseMaker)
30
31//________________________________________________________________________
32AliJetResponseMaker::AliJetResponseMaker() :
e44e8726 33 AliAnalysisTaskEmcalJet("AliJetResponseMaker", kTRUE),
cd6431de 34 fTracks2Name(""),
35 fCalo2Name(""),
36 fJets2Name(""),
37 fRho2Name(""),
38 fPtBiasJet2Track(0),
39 fPtBiasJet2Clus(0),
40 fAreCollections1MC(kFALSE),
41 fAreCollections2MC(kTRUE),
42 fMatching(kNoMatching),
7f76e479 43 fMatchingPar1(0),
44 fMatchingPar2(0),
cd6431de 45 fJet2MinEta(-999),
46 fJet2MaxEta(-999),
47 fJet2MinPhi(-999),
48 fJet2MaxPhi(-999),
4643d2e8 49 fSelectPtHardBin(-999),
1ee1b5b8 50 fPythiaHeader(0),
1ee1b5b8 51 fPtHardBin(0),
52 fNTrials(0),
cd6431de 53 fTracks2(0),
54 fCaloClusters2(0),
55 fJets2(0),
56 fRho2(0),
57 fRho2Val(0),
cfc2ac24 58 fTracks2Map(0),
ca5c29fa 59 fHistTrialsAfterSel(0),
60 fHistEventsAfterSel(0),
61 fHistTrials(0),
62 fHistXsection(0),
1ee1b5b8 63 fHistEvents(0),
cd6431de 64 fHistJets1PhiEta(0),
65 fHistJets1PtArea(0),
66 fHistJets1CorrPtArea(0),
ca5c29fa 67 fHistLeadingJets1PtArea(0),
68 fHistLeadingJets1CorrPtArea(0),
63fac07f 69 fHistJets1NEFvsPt(0),
70 fHistJets1CEFvsCEFPt(0),
71 fHistJets1ZvsPt(0),
cd6431de 72 fHistJets2PhiEta(0),
73 fHistJets2PtArea(0),
74 fHistJets2CorrPtArea(0),
ca5c29fa 75 fHistLeadingJets2PtArea(0),
76 fHistLeadingJets2CorrPtArea(0),
cfc2ac24 77 fHistJets2PhiEtaAcceptance(0),
78 fHistJets2PtAreaAcceptance(0),
79 fHistJets2CorrPtAreaAcceptance(0),
ca5c29fa 80 fHistLeadingJets2PtAreaAcceptance(0),
81 fHistLeadingJets2CorrPtAreaAcceptance(0),
63fac07f 82 fHistJets2NEFvsPt(0),
83 fHistJets2CEFvsCEFPt(0),
84 fHistJets2ZvsPt(0),
ca5c29fa 85 fHistCommonEnergy1vsJet1Pt(0),
86 fHistCommonEnergy2vsJet2Pt(0),
87 fHistDistancevsJet1Pt(0),
88 fHistDistancevsJet2Pt(0),
6a20534a 89 fHistDistancevsCommonEnergy1(0),
90 fHistDistancevsCommonEnergy2(0),
c560b734 91 fHistJet2PtOverJet1PtvsJet2Pt(0),
7030f36f 92 fHistJet1PtOverJet2PtvsJet1Pt(0),
ca5c29fa 93 fHistDeltaEtaPhi(0),
7030f36f 94 fHistDeltaPtvsJet1Pt(0),
cfc2ac24 95 fHistDeltaPtvsJet2Pt(0),
96 fHistDeltaPtvsMatchingLevel(0),
7030f36f 97 fHistDeltaCorrPtvsJet1Pt(0),
cfc2ac24 98 fHistDeltaCorrPtvsJet2Pt(0),
99 fHistDeltaCorrPtvsMatchingLevel(0),
cd6431de 100 fHistNonMatchedJets1PtArea(0),
101 fHistNonMatchedJets2PtArea(0),
102 fHistNonMatchedJets1CorrPtArea(0),
103 fHistNonMatchedJets2CorrPtArea(0),
104 fHistJet1PtvsJet2Pt(0),
105 fHistJet1CorrPtvsJet2CorrPt(0),
106 fHistMissedJets2PtArea(0)
2949a2ec 107{
108 // Default constructor.
4643d2e8 109
a487deae 110 SetMakeGeneralHistograms(kTRUE);
2949a2ec 111}
112
113//________________________________________________________________________
114AliJetResponseMaker::AliJetResponseMaker(const char *name) :
e44e8726 115 AliAnalysisTaskEmcalJet(name, kTRUE),
cd6431de 116 fTracks2Name("MCParticles"),
117 fCalo2Name(""),
118 fJets2Name("MCJets"),
119 fRho2Name(""),
120 fPtBiasJet2Track(0),
121 fPtBiasJet2Clus(0),
122 fAreCollections1MC(kFALSE),
123 fAreCollections2MC(kTRUE),
124 fMatching(kNoMatching),
7f76e479 125 fMatchingPar1(0),
126 fMatchingPar2(0),
cd6431de 127 fJet2MinEta(-999),
128 fJet2MaxEta(-999),
129 fJet2MinPhi(-999),
130 fJet2MaxPhi(-999),
4643d2e8 131 fSelectPtHardBin(-999),
1ee1b5b8 132 fPythiaHeader(0),
1ee1b5b8 133 fPtHardBin(0),
134 fNTrials(0),
cd6431de 135 fTracks2(0),
136 fCaloClusters2(0),
137 fJets2(0),
138 fRho2(0),
139 fRho2Val(0),
cfc2ac24 140 fTracks2Map(0),
ca5c29fa 141 fHistTrialsAfterSel(0),
142 fHistEventsAfterSel(0),
143 fHistTrials(0),
144 fHistXsection(0),
1ee1b5b8 145 fHistEvents(0),
cd6431de 146 fHistJets1PhiEta(0),
147 fHistJets1PtArea(0),
148 fHistJets1CorrPtArea(0),
ca5c29fa 149 fHistLeadingJets1PtArea(0),
150 fHistLeadingJets1CorrPtArea(0),
63fac07f 151 fHistJets1NEFvsPt(0),
152 fHistJets1CEFvsCEFPt(0),
153 fHistJets1ZvsPt(0),
cd6431de 154 fHistJets2PhiEta(0),
155 fHistJets2PtArea(0),
156 fHistJets2CorrPtArea(0),
ca5c29fa 157 fHistLeadingJets2PtArea(0),
158 fHistLeadingJets2CorrPtArea(0),
cfc2ac24 159 fHistJets2PhiEtaAcceptance(0),
160 fHistJets2PtAreaAcceptance(0),
161 fHistJets2CorrPtAreaAcceptance(0),
ca5c29fa 162 fHistLeadingJets2PtAreaAcceptance(0),
163 fHistLeadingJets2CorrPtAreaAcceptance(0),
63fac07f 164 fHistJets2NEFvsPt(0),
165 fHistJets2CEFvsCEFPt(0),
166 fHistJets2ZvsPt(0),
ca5c29fa 167 fHistCommonEnergy1vsJet1Pt(0),
168 fHistCommonEnergy2vsJet2Pt(0),
169 fHistDistancevsJet1Pt(0),
170 fHistDistancevsJet2Pt(0),
6a20534a 171 fHistDistancevsCommonEnergy1(0),
172 fHistDistancevsCommonEnergy2(0),
c560b734 173 fHistJet2PtOverJet1PtvsJet2Pt(0),
7030f36f 174 fHistJet1PtOverJet2PtvsJet1Pt(0),
ca5c29fa 175 fHistDeltaEtaPhi(0),
7030f36f 176 fHistDeltaPtvsJet1Pt(0),
cfc2ac24 177 fHistDeltaPtvsJet2Pt(0),
178 fHistDeltaPtvsMatchingLevel(0),
7030f36f 179 fHistDeltaCorrPtvsJet1Pt(0),
cfc2ac24 180 fHistDeltaCorrPtvsJet2Pt(0),
181 fHistDeltaCorrPtvsMatchingLevel(0),
cd6431de 182 fHistNonMatchedJets1PtArea(0),
183 fHistNonMatchedJets2PtArea(0),
184 fHistNonMatchedJets1CorrPtArea(0),
185 fHistNonMatchedJets2CorrPtArea(0),
186 fHistJet1PtvsJet2Pt(0),
187 fHistJet1CorrPtvsJet2CorrPt(0),
188 fHistMissedJets2PtArea(0)
2949a2ec 189{
190 // Standard constructor.
4643d2e8 191
a487deae 192 SetMakeGeneralHistograms(kTRUE);
2949a2ec 193}
194
195//________________________________________________________________________
196AliJetResponseMaker::~AliJetResponseMaker()
197{
198 // Destructor
199}
200
ca5c29fa 201//________________________________________________________________________
202Bool_t AliJetResponseMaker::PythiaInfoFromFile(const char* currFile, Float_t &fXsec, Float_t &fTrials, Int_t &pthard)
203{
204 //
205 // Get the cross section and the trails either from pyxsec.root or from pysec_hists.root
206 // Get the pt hard bin from the file path
207 // This is to called in Notify and should provide the path to the AOD/ESD file
208 // (Partially copied from AliAnalysisHelperJetTasks)
209
210 TString file(currFile);
211 fXsec = 0;
212 fTrials = 1;
213
214 if(file.Contains("root_archive.zip#")){
215 Ssiz_t pos1 = file.Index("root_archive",12,0,TString::kExact);
216 Ssiz_t pos = file.Index("#",1,pos1,TString::kExact);
217 Ssiz_t pos2 = file.Index(".root",5,TString::kExact);
218 file.Replace(pos+1,pos2-pos1,"");
219 }
220 else {
221 // not an archive take the basename....
222 file.ReplaceAll(gSystem->BaseName(file.Data()),"");
223 }
224 Printf("%s",file.Data());
225
226 // Get the pt hard bin
227 TString strPthard(file);
228 strPthard.Remove(strPthard.Last('/'));
229 strPthard.Remove(strPthard.Last('/'));
230 strPthard.Remove(0,strPthard.Last('/')+1);
231 if (strPthard.IsDec())
232 pthard = strPthard.Atoi();
233 else
234 AliWarning(Form("Could not extract file number from path %s", strPthard.Data()));
235
236 TFile *fxsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec.root")); // problem that we cannot really test the existance of a file in a archive so we have to lvie with open error message from root
237 if(!fxsec){
238 // next trial fetch the histgram file
239 fxsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec_hists.root"));
240 if(!fxsec){
241 // not a severe condition but inciate that we have no information
242 return kFALSE;
243 }
244 else{
245 // find the tlist we want to be independtent of the name so use the Tkey
246 TKey* key = (TKey*)fxsec->GetListOfKeys()->At(0);
247 if(!key){
248 fxsec->Close();
249 return kFALSE;
250 }
251 TList *list = dynamic_cast<TList*>(key->ReadObj());
252 if(!list){
253 fxsec->Close();
254 return kFALSE;
255 }
256 fXsec = ((TProfile*)list->FindObject("h1Xsec"))->GetBinContent(1);
257 fTrials = ((TH1F*)list->FindObject("h1Trials"))->GetBinContent(1);
258 fxsec->Close();
259 }
260 } // no tree pyxsec.root
261 else {
262 TTree *xtree = (TTree*)fxsec->Get("Xsection");
263 if(!xtree){
264 fxsec->Close();
265 return kFALSE;
266 }
267 UInt_t ntrials = 0;
268 Double_t xsection = 0;
269 xtree->SetBranchAddress("xsection",&xsection);
270 xtree->SetBranchAddress("ntrials",&ntrials);
271 xtree->GetEntry(0);
272 fTrials = ntrials;
273 fXsec = xsection;
274 fxsec->Close();
275 }
276 return kTRUE;
277}
278
279//________________________________________________________________________
280Bool_t AliJetResponseMaker::UserNotify()
281{
282 TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
283 if (!tree) {
284 AliError(Form("%s - UserNotify: No current tree!",GetName()));
285 return kFALSE;
286 }
287
288 Float_t xsection = 0;
289 Float_t trials = 0;
290 Int_t pthard = 0;
291
292 TFile *curfile = tree->GetCurrentFile();
293 if (!curfile) {
294 AliError(Form("%s - UserNotify: No current file!",GetName()));
295 return kFALSE;
296 }
297
298 PythiaInfoFromFile(curfile->GetName(), xsection, trials, pthard);
299
300 fHistTrials->SetBinContent(pthard + 1, fHistTrials->GetBinContent(pthard + 1) + trials);
301 fHistXsection->SetBinContent(pthard + 1, fHistXsection->GetBinContent(pthard + 1) + xsection);
302 fHistEvents->SetBinContent(pthard + 1, fHistEvents->GetBinContent(pthard + 1) + 1);
303
304 return kTRUE;
305}
306
2949a2ec 307//________________________________________________________________________
308void AliJetResponseMaker::UserCreateOutputObjects()
309{
310 // Create user objects.
311
a487deae 312 AliAnalysisTaskEmcalJet::UserCreateOutputObjects();
2949a2ec 313
ca5c29fa 314 fHistTrialsAfterSel = new TH1F("fHistTrialsAfterSel", "fHistTrialsAfterSel", 11, 0, 11);
315 fHistTrialsAfterSel->GetXaxis()->SetTitle("p_{T} hard bin");
316 fHistTrialsAfterSel->GetYaxis()->SetTitle("trials");
317 fOutput->Add(fHistTrialsAfterSel);
318
319 fHistEventsAfterSel = new TH1F("fHistEventsAfterSel", "fHistEventsAfterSel", 11, 0, 11);
320 fHistEventsAfterSel->GetXaxis()->SetTitle("p_{T} hard bin");
321 fHistEventsAfterSel->GetYaxis()->SetTitle("total events");
322 fOutput->Add(fHistEventsAfterSel);
1ee1b5b8 323
ca5c29fa 324 fHistTrials = new TH1F("fHistTrials", "fHistTrials", 11, 0, 11);
325 fHistTrials->GetXaxis()->SetTitle("p_{T} hard bin");
326 fHistTrials->GetYaxis()->SetTitle("trials");
327 fOutput->Add(fHistTrials);
328
329 fHistXsection = new TH1F("fHistXsection", "fHistXsection", 11, 0, 11);
330 fHistXsection->GetXaxis()->SetTitle("p_{T} hard bin");
331 fHistXsection->GetYaxis()->SetTitle("xsection");
332 fOutput->Add(fHistXsection);
1ee1b5b8 333
1ee1b5b8 334 fHistEvents = new TH1F("fHistEvents", "fHistEvents", 11, 0, 11);
335 fHistEvents->GetXaxis()->SetTitle("p_{T} hard bin");
336 fHistEvents->GetYaxis()->SetTitle("total events");
337 fOutput->Add(fHistEvents);
338
ca5c29fa 339 const Int_t ptHardLo[11] = { 0, 5,11,21,36,57, 84,117,152,191,234};
340 const Int_t ptHardHi[11] = { 5,11,21,36,57,84,117,152,191,234,1000000};
341
1ee1b5b8 342 for (Int_t i = 1; i < 12; i++) {
ca5c29fa 343 fHistTrialsAfterSel->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
344 fHistEventsAfterSel->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
345
346 fHistTrials->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
347 fHistXsection->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
1ee1b5b8 348 fHistEvents->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
349 }
350
cfc2ac24 351 fHistJets1PhiEta = new TH2F("fHistJets1PhiEta", "fHistJets1PhiEta", 40, -1, 1, 40, 0, TMath::Pi()*2);
cd6431de 352 fHistJets1PhiEta->GetXaxis()->SetTitle("#eta");
353 fHistJets1PhiEta->GetYaxis()->SetTitle("#phi");
354 fOutput->Add(fHistJets1PhiEta);
2949a2ec 355
cd6431de 356 fHistJets1PtArea = new TH2F("fHistJets1PtArea", "fHistJets1PtArea", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, fNbins, fMinBinPt, fMaxBinPt);
357 fHistJets1PtArea->GetXaxis()->SetTitle("area");
358 fHistJets1PtArea->GetYaxis()->SetTitle("p_{T,1} (GeV/c)");
359 fHistJets1PtArea->GetZaxis()->SetTitle("counts");
360 fOutput->Add(fHistJets1PtArea);
361
ca5c29fa 362 fHistLeadingJets1PtArea = new TH2F("fHistLeadingJets1PtArea", "fHistLeadingJets1PtArea", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, fNbins, fMinBinPt, fMaxBinPt);
363 fHistLeadingJets1PtArea->GetXaxis()->SetTitle("area");
364 fHistLeadingJets1PtArea->GetYaxis()->SetTitle("p_{T,1} (GeV/c)");
365 fHistLeadingJets1PtArea->GetZaxis()->SetTitle("counts");
366 fOutput->Add(fHistLeadingJets1PtArea);
367
cd6431de 368 if (!fRhoName.IsNull()) {
369 fHistJets1CorrPtArea = new TH2F("fHistJets1CorrPtArea", "fHistJets1CorrPtArea", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, 2*fNbins, -fMaxBinPt, fMaxBinPt);
370 fHistJets1CorrPtArea->GetXaxis()->SetTitle("area");
371 fHistJets1CorrPtArea->GetYaxis()->SetTitle("p_{T,1}^{corr} (GeV/c)");
372 fHistJets1CorrPtArea->GetZaxis()->SetTitle("counts");
373 fOutput->Add(fHistJets1CorrPtArea);
ca5c29fa 374
375 fHistLeadingJets1CorrPtArea = new TH2F("fHistLeadingJets1CorrPtArea", "fHistLeadingJets1CorrPtArea", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, 2*fNbins, -fMaxBinPt, fMaxBinPt);
376 fHistLeadingJets1CorrPtArea->GetXaxis()->SetTitle("area");
377 fHistLeadingJets1CorrPtArea->GetYaxis()->SetTitle("p_{T,1}^{corr} (GeV/c)");
378 fHistLeadingJets1CorrPtArea->GetZaxis()->SetTitle("counts");
379 fOutput->Add(fHistLeadingJets1CorrPtArea);
cd6431de 380 }
381
63fac07f 382 fHistJets1ZvsPt = new TH2F("fHistJets1ZvsPt", "fHistJets1ZvsPt", 120, 0, 1.2, fNbins, fMinBinPt, fMaxBinPt);
383 fHistJets1ZvsPt->GetXaxis()->SetTitle("Z");
384 fHistJets1ZvsPt->GetYaxis()->SetTitle("p_{T,1} (GeV/c)");
385 fHistJets1ZvsPt->GetZaxis()->SetTitle("counts");
386 fOutput->Add(fHistJets1ZvsPt);
387
388 fHistJets1NEFvsPt = new TH2F("fHistJets1NEFvsPt", "fHistJets1NEFvsPt", 120, 0, 1.2, fNbins, fMinBinPt, fMaxBinPt);
389 fHistJets1NEFvsPt->GetXaxis()->SetTitle("NEF");
390 fHistJets1NEFvsPt->GetYaxis()->SetTitle("p_{T,1} (GeV/c)");
391 fHistJets1NEFvsPt->GetZaxis()->SetTitle("counts");
392 fOutput->Add(fHistJets1NEFvsPt);
393
394 fHistJets1CEFvsCEFPt = new TH2F("fHistJets1CEFvsCEFPt", "fHistJets1CEFvsCEFPt", 120, 0, 1.2, fNbins, fMinBinPt, fMaxBinPt);
395 fHistJets1CEFvsCEFPt->GetXaxis()->SetTitle("1-NEF");
396 fHistJets1CEFvsCEFPt->GetYaxis()->SetTitle("(1-NEF)*p_{T,1} (GeV/c)");
397 fHistJets1CEFvsCEFPt->GetZaxis()->SetTitle("counts");
398 fOutput->Add(fHistJets1CEFvsCEFPt);
399
cfc2ac24 400 fHistJets2PhiEta = new TH2F("fHistJets2PhiEta", "fHistJets2PhiEta", 40, -1, 1, 40, 0, TMath::Pi()*2);
cd6431de 401 fHistJets2PhiEta->GetXaxis()->SetTitle("#eta");
402 fHistJets2PhiEta->GetYaxis()->SetTitle("#phi");
403 fOutput->Add(fHistJets2PhiEta);
cfc2ac24 404
405 fHistJets2PhiEtaAcceptance = new TH2F("fHistJets2PhiEtaAcceptance", "fHistJets2PhiEtaAcceptance", 40, -1, 1, 40, 0, TMath::Pi()*2);
406 fHistJets2PhiEtaAcceptance->GetXaxis()->SetTitle("#eta");
407 fHistJets2PhiEtaAcceptance->GetYaxis()->SetTitle("#phi");
408 fOutput->Add(fHistJets2PhiEtaAcceptance);
7549c451 409
cd6431de 410 fHistJets2PtArea = new TH2F("fHistJets2PtArea", "fHistJets2PtArea", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, fNbins, fMinBinPt, fMaxBinPt);
411 fHistJets2PtArea->GetXaxis()->SetTitle("area");
412 fHistJets2PtArea->GetYaxis()->SetTitle("p_{T,2} (GeV/c)");
413 fHistJets2PtArea->GetZaxis()->SetTitle("counts");
414 fOutput->Add(fHistJets2PtArea);
415
ca5c29fa 416 fHistLeadingJets2PtArea = new TH2F("fHistLeadingJets2PtArea", "fHistLeadingJets2PtArea", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, fNbins, fMinBinPt, fMaxBinPt);
417 fHistLeadingJets2PtArea->GetXaxis()->SetTitle("area");
418 fHistLeadingJets2PtArea->GetYaxis()->SetTitle("p_{T,2} (GeV/c)");
419 fHistLeadingJets2PtArea->GetZaxis()->SetTitle("counts");
420 fOutput->Add(fHistLeadingJets2PtArea);
421
cfc2ac24 422 fHistJets2PtAreaAcceptance = new TH2F("fHistJets2PtAreaAcceptance", "fHistJets2PtAreaAcceptance", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, fNbins, fMinBinPt, fMaxBinPt);
423 fHistJets2PtAreaAcceptance->GetXaxis()->SetTitle("area");
424 fHistJets2PtAreaAcceptance->GetYaxis()->SetTitle("p_{T,2} (GeV/c)");
425 fHistJets2PtAreaAcceptance->GetZaxis()->SetTitle("counts");
426 fOutput->Add(fHistJets2PtAreaAcceptance);
427
ca5c29fa 428 fHistLeadingJets2PtAreaAcceptance = new TH2F("fHistLeadingJets2PtAreaAcceptance", "fHistLeadingJets2PtAreaAcceptance", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, fNbins, fMinBinPt, fMaxBinPt);
429 fHistLeadingJets2PtAreaAcceptance->GetXaxis()->SetTitle("area");
430 fHistLeadingJets2PtAreaAcceptance->GetYaxis()->SetTitle("p_{T,2} (GeV/c)");
431 fHistLeadingJets2PtAreaAcceptance->GetZaxis()->SetTitle("counts");
432 fOutput->Add(fHistLeadingJets2PtAreaAcceptance);
433
cd6431de 434 if (!fRho2Name.IsNull()) {
435 fHistJets2CorrPtArea = new TH2F("fHistJets2CorrPtArea", "fHistJets2CorrPtArea", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, 2*fNbins, -fMaxBinPt, fMaxBinPt);
436 fHistJets2CorrPtArea->GetXaxis()->SetTitle("area");
437 fHistJets2CorrPtArea->GetYaxis()->SetTitle("p_{T,2}^{corr} (GeV/c)");
438 fHistJets2CorrPtArea->GetZaxis()->SetTitle("counts");
439 fOutput->Add(fHistJets2CorrPtArea);
cfc2ac24 440
ca5c29fa 441 fHistLeadingJets2CorrPtArea = new TH2F("fHistLeadingJets2CorrPtArea", "fHistLeadingJets2CorrPtArea", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, 2*fNbins, -fMaxBinPt, fMaxBinPt);
442 fHistLeadingJets2CorrPtArea->GetXaxis()->SetTitle("area");
443 fHistLeadingJets2CorrPtArea->GetYaxis()->SetTitle("p_{T,2}^{corr} (GeV/c)");
444 fHistLeadingJets2CorrPtArea->GetZaxis()->SetTitle("counts");
445 fOutput->Add(fHistLeadingJets2CorrPtArea);
446
cfc2ac24 447 fHistJets2CorrPtAreaAcceptance = new TH2F("fHistJets2CorrPtAreaAcceptance", "fHistJets2CorrPtAreaAcceptance", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, 2*fNbins, -fMaxBinPt, fMaxBinPt);
448 fHistJets2CorrPtAreaAcceptance->GetXaxis()->SetTitle("area");
449 fHistJets2CorrPtAreaAcceptance->GetYaxis()->SetTitle("p_{T,2}^{corr} (GeV/c)");
450 fHistJets2CorrPtAreaAcceptance->GetZaxis()->SetTitle("counts");
451 fOutput->Add(fHistJets2CorrPtAreaAcceptance);
ca5c29fa 452
453 fHistLeadingJets2CorrPtAreaAcceptance = new TH2F("fHistLeadingJets2CorrPtAreaAcceptance", "fHistLeadingJets2CorrPtAreaAcceptance", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, 2*fNbins, -fMaxBinPt, fMaxBinPt);
454 fHistLeadingJets2CorrPtAreaAcceptance->GetXaxis()->SetTitle("area");
455 fHistLeadingJets2CorrPtAreaAcceptance->GetYaxis()->SetTitle("p_{T,2}^{corr} (GeV/c)");
456 fHistLeadingJets2CorrPtAreaAcceptance->GetZaxis()->SetTitle("counts");
457 fOutput->Add(fHistLeadingJets2CorrPtAreaAcceptance);
cd6431de 458 }
459
63fac07f 460 fHistJets2ZvsPt = new TH2F("fHistJets2ZvsPt", "fHistJets2ZvsPt", 120, 0, 1.2, fNbins, fMinBinPt, fMaxBinPt);
461 fHistJets2ZvsPt->GetXaxis()->SetTitle("Z");
462 fHistJets2ZvsPt->GetYaxis()->SetTitle("p_{T,2} (GeV/c)");
463 fHistJets2ZvsPt->GetZaxis()->SetTitle("counts");
464 fOutput->Add(fHistJets2ZvsPt);
465
466 fHistJets2NEFvsPt = new TH2F("fHistJets2NEFvsPt", "fHistJets2NEFvsPt", 120, 0, 1.2, fNbins, fMinBinPt, fMaxBinPt);
467 fHistJets2NEFvsPt->GetXaxis()->SetTitle("NEF");
468 fHistJets2NEFvsPt->GetYaxis()->SetTitle("p_{T,2} (GeV/c)");
469 fHistJets2NEFvsPt->GetZaxis()->SetTitle("counts");
470 fOutput->Add(fHistJets2NEFvsPt);
471
472 fHistJets2CEFvsCEFPt = new TH2F("fHistJets2CEFvsCEFPt", "fHistJets2CEFvsCEFPt", 120, 0, 1.2, fNbins, fMinBinPt, fMaxBinPt);
473 fHistJets2CEFvsCEFPt->GetXaxis()->SetTitle("1-NEF");
474 fHistJets2CEFvsCEFPt->GetYaxis()->SetTitle("(1-NEF)*p_{T,2} (GeV/c)");
475 fHistJets2CEFvsCEFPt->GetZaxis()->SetTitle("counts");
476 fOutput->Add(fHistJets2CEFvsCEFPt);
477
ca5c29fa 478 fHistCommonEnergy1vsJet1Pt = new TH2F("fHistCommonEnergy1vsJet1Pt", "fHistCommonEnergy1vsJet1Pt", fNbins/2, 0, 1.2, fNbins/2, fMinBinPt, fMaxBinPt);
479 fHistCommonEnergy1vsJet1Pt->GetXaxis()->SetTitle("Common energy 1 (%)");
480 fHistCommonEnergy1vsJet1Pt->GetYaxis()->SetTitle("p_{T,1}");
481 fHistCommonEnergy1vsJet1Pt->GetZaxis()->SetTitle("counts");
482 fOutput->Add(fHistCommonEnergy1vsJet1Pt);
483
484 fHistCommonEnergy2vsJet2Pt = new TH2F("fHistCommonEnergy2vsJet2Pt", "fHistCommonEnergy2vsJet2Pt", fNbins/2, 0, 1.2, fNbins/2, fMinBinPt, fMaxBinPt);
485 fHistCommonEnergy2vsJet2Pt->GetXaxis()->SetTitle("Common energy 2 (%)");
486 fHistCommonEnergy2vsJet2Pt->GetYaxis()->SetTitle("p_{T,2}");
487 fHistCommonEnergy2vsJet2Pt->GetZaxis()->SetTitle("counts");
488 fOutput->Add(fHistCommonEnergy2vsJet2Pt);
489
490 fHistDistancevsJet1Pt = new TH2F("fHistDistancevsJet1Pt", "fHistDistancevsJet1Pt", fNbins/2, 0, 1.2, fNbins/2, fMinBinPt, fMaxBinPt);
491 fHistDistancevsJet1Pt->GetXaxis()->SetTitle("Distance");
492 fHistDistancevsJet1Pt->GetYaxis()->SetTitle("p_{T,1}");
493 fHistDistancevsJet1Pt->GetZaxis()->SetTitle("counts");
494 fOutput->Add(fHistDistancevsJet1Pt);
495
496 fHistDistancevsJet2Pt = new TH2F("fHistDistancevsJet2Pt", "fHistDistancevsJet2Pt", fNbins/2, 0, 1.2, fNbins/2, fMinBinPt, fMaxBinPt);
497 fHistDistancevsJet2Pt->GetXaxis()->SetTitle("Distance");
498 fHistDistancevsJet2Pt->GetYaxis()->SetTitle("p_{T,2}");
499 fHistDistancevsJet2Pt->GetZaxis()->SetTitle("counts");
500 fOutput->Add(fHistDistancevsJet2Pt);
cd6431de 501
6a20534a 502 fHistDistancevsCommonEnergy1 = new TH2F("fHistDistancevsCommonEnergy1", "fHistDistancevsCommonEnergy1", fNbins/2, 0, 1.2, fNbins/2, 0, 1.2);
503 fHistDistancevsCommonEnergy1->GetXaxis()->SetTitle("Distance");
504 fHistDistancevsCommonEnergy1->GetYaxis()->SetTitle("Common energy 1 (%)");
505 fHistDistancevsCommonEnergy1->GetZaxis()->SetTitle("counts");
506 fOutput->Add(fHistDistancevsCommonEnergy1);
507
508 fHistDistancevsCommonEnergy2 = new TH2F("fHistDistancevsCommonEnergy2", "fHistDistancevsCommonEnergy2", fNbins/2, 0, 1.2, fNbins/2, 0, 1.2);
509 fHistDistancevsCommonEnergy2->GetXaxis()->SetTitle("Distance");
510 fHistDistancevsCommonEnergy2->GetYaxis()->SetTitle("Common energy 2 (%)");
511 fHistDistancevsCommonEnergy2->GetZaxis()->SetTitle("counts");
512 fOutput->Add(fHistDistancevsCommonEnergy2);
cfc2ac24 513
7030f36f 514 fHistJet2PtOverJet1PtvsJet2Pt = new TH2F("fHistJet2PtOverJet1PtvsJet2Pt", "fHistJet2PtOverJet1PtvsJet2Pt", fNbins, fMinBinPt, fMaxBinPt, 300, 0, 1.5);
c560b734 515 fHistJet2PtOverJet1PtvsJet2Pt->GetXaxis()->SetTitle("p_{T,2}");
516 fHistJet2PtOverJet1PtvsJet2Pt->GetYaxis()->SetTitle("p_{T,2} / p_{T,1}");
517 fHistJet2PtOverJet1PtvsJet2Pt->GetZaxis()->SetTitle("counts");
518 fOutput->Add(fHistJet2PtOverJet1PtvsJet2Pt);
519
7030f36f 520 fHistJet1PtOverJet2PtvsJet1Pt = new TH2F("fHistJet1PtOverJet2PtvsJet1Pt", "fHistJet1PtOverJet2PtvsJet1Pt", fNbins, fMinBinPt, fMaxBinPt, 300, 0, 1.5);
521 fHistJet1PtOverJet2PtvsJet1Pt->GetXaxis()->SetTitle("p_{T,1}");
522 fHistJet1PtOverJet2PtvsJet1Pt->GetYaxis()->SetTitle("p_{T,1} / p_{T,2}");
523 fHistJet1PtOverJet2PtvsJet1Pt->GetZaxis()->SetTitle("counts");
524 fOutput->Add(fHistJet1PtOverJet2PtvsJet1Pt);
525
ca5c29fa 526 fHistDeltaEtaPhi = new TH2F("fHistDeltaEtaPhi", "fHistDeltaEtaPhi", 200, -1, 1, 250, -1.6, 4.8);
527 fHistDeltaEtaPhi->GetXaxis()->SetTitle("#Delta#eta");
528 fHistDeltaEtaPhi->GetYaxis()->SetTitle("#Delta#phi");
529 fHistDeltaEtaPhi->GetZaxis()->SetTitle("counts");
530 fOutput->Add(fHistDeltaEtaPhi);
cfc2ac24 531
7030f36f 532 fHistDeltaPtvsJet1Pt = new TH2F("fHistDeltaPtvsJet1Pt", "fHistDeltaPtvsJet1Pt", fNbins, fMinBinPt, fMaxBinPt, 2*fNbins, -fMaxBinPt, fMaxBinPt);
533 fHistDeltaPtvsJet1Pt->GetXaxis()->SetTitle("p_{T,1}");
534 fHistDeltaPtvsJet1Pt->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
535 fHistDeltaPtvsJet1Pt->GetZaxis()->SetTitle("counts");
536 fOutput->Add(fHistDeltaPtvsJet1Pt);
537
c560b734 538 fHistDeltaPtvsJet2Pt = new TH2F("fHistDeltaPtvsJet2Pt", "fHistDeltaPtvsJet2Pt", fNbins, fMinBinPt, fMaxBinPt, 2*fNbins, -fMaxBinPt, fMaxBinPt);
cfc2ac24 539 fHistDeltaPtvsJet2Pt->GetXaxis()->SetTitle("p_{T,2}");
540 fHistDeltaPtvsJet2Pt->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
541 fHistDeltaPtvsJet2Pt->GetZaxis()->SetTitle("counts");
542 fOutput->Add(fHistDeltaPtvsJet2Pt);
543
544 fHistDeltaPtvsMatchingLevel = new TH2F("fHistDeltaPtvsMatchingLevel", "fHistDeltaPtvsMatchingLevel", fNbins/2, 0, 1.2, 2*fNbins, -fMaxBinPt, fMaxBinPt);
545 fHistDeltaPtvsMatchingLevel->GetXaxis()->SetTitle("Matching level");
546 fHistDeltaPtvsMatchingLevel->GetYaxis()->SetTitle("#Deltap_{T} (GeV/c)");
547 fHistDeltaPtvsMatchingLevel->GetZaxis()->SetTitle("counts");
548 fOutput->Add(fHistDeltaPtvsMatchingLevel);
cd6431de 549
550 if (!fRhoName.IsNull() || !fRho2Name.IsNull()) {
7030f36f 551 fHistDeltaCorrPtvsJet1Pt = new TH2F("fHistDeltaCorrPtvsJet1Pt", "fHistDeltaCorrPtvsJet1Pt", fNbins/2, fMinBinPt, fMaxBinPt, 2*fNbins, -fMaxBinPt, fMaxBinPt);
552 fHistDeltaCorrPtvsJet1Pt->GetXaxis()->SetTitle("p_{T,1}");
553 fHistDeltaCorrPtvsJet1Pt->GetYaxis()->SetTitle("#Deltap_{T}^{corr} (GeV/c)");
554 fHistDeltaCorrPtvsJet1Pt->GetZaxis()->SetTitle("counts");
555 fOutput->Add(fHistDeltaCorrPtvsJet1Pt);
556
cfc2ac24 557 fHistDeltaCorrPtvsJet2Pt = new TH2F("fHistDeltaCorrPtvsJet2Pt", "fHistDeltaCorrPtvsJet2Pt", fNbins/2, fMinBinPt, fMaxBinPt, 2*fNbins, -fMaxBinPt, fMaxBinPt);
558 fHistDeltaCorrPtvsJet2Pt->GetXaxis()->SetTitle("p_{T,2}");
559 fHistDeltaCorrPtvsJet2Pt->GetYaxis()->SetTitle("#Deltap_{T}^{corr} (GeV/c)");
560 fHistDeltaCorrPtvsJet2Pt->GetZaxis()->SetTitle("counts");
561 fOutput->Add(fHistDeltaCorrPtvsJet2Pt);
562
563 fHistDeltaCorrPtvsMatchingLevel = new TH2F("fHistDeltaCorrPtvsMatchingLevel", "fHistDeltaCorrPtvsMatchingLevel", fNbins/2, 0, 1.2, 2*fNbins, -fMaxBinPt, fMaxBinPt);
564 fHistDeltaCorrPtvsMatchingLevel->GetXaxis()->SetTitle("Matching level");
565 fHistDeltaCorrPtvsMatchingLevel->GetYaxis()->SetTitle("#Deltap_{T}^{corr} (GeV/c)");
566 fHistDeltaCorrPtvsMatchingLevel->GetZaxis()->SetTitle("counts");
787a3c4f 567 fOutput->Add(fHistDeltaCorrPtvsMatchingLevel);
cd6431de 568 }
569
570 fHistNonMatchedJets1PtArea = new TH2F("fHistNonMatchedJets1PtArea", "fHistNonMatchedJets1PtArea", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, fNbins, fMinBinPt, fMaxBinPt);
571 fHistNonMatchedJets1PtArea->GetXaxis()->SetTitle("area");
572 fHistNonMatchedJets1PtArea->GetYaxis()->SetTitle("p_{T,1} (GeV/c)");
573 fHistNonMatchedJets1PtArea->GetZaxis()->SetTitle("counts");
574 fOutput->Add(fHistNonMatchedJets1PtArea);
575
576 fHistNonMatchedJets2PtArea = new TH2F("fHistNonMatchedJets2PtArea", "fHistNonMatchedJets2PtArea", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, fNbins, fMinBinPt, fMaxBinPt);
577 fHistNonMatchedJets2PtArea->GetXaxis()->SetTitle("area");
578 fHistNonMatchedJets2PtArea->GetYaxis()->SetTitle("p_{T,2} (GeV/c)");
579 fHistNonMatchedJets2PtArea->GetZaxis()->SetTitle("counts");
580 fOutput->Add(fHistNonMatchedJets2PtArea);
581
582 if (!fRhoName.IsNull()) {
583 fHistNonMatchedJets1CorrPtArea = new TH2F("fHistNonMatchedJets1CorrPtArea", "fHistNonMatchedJets1CorrPtArea", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, 2*fNbins, -fMaxBinPt, fMaxBinPt);
584 fHistNonMatchedJets1CorrPtArea->GetXaxis()->SetTitle("area");
585 fHistNonMatchedJets1CorrPtArea->GetYaxis()->SetTitle("p_{T,1} (GeV/c)");
586 fHistNonMatchedJets1CorrPtArea->GetZaxis()->SetTitle("counts");
587 fOutput->Add(fHistNonMatchedJets1CorrPtArea);
588 }
589
590 if (!fRho2Name.IsNull()) {
591 fHistNonMatchedJets2CorrPtArea = new TH2F("fHistNonMatchedJets2CorrPtArea", "fHistNonMatchedJets2CorrPtArea", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, 2*fNbins, -fMaxBinPt, fMaxBinPt);
592 fHistNonMatchedJets2CorrPtArea->GetXaxis()->SetTitle("area");
593 fHistNonMatchedJets2CorrPtArea->GetYaxis()->SetTitle("p_{T,2} (GeV/c)");
594 fHistNonMatchedJets2CorrPtArea->GetZaxis()->SetTitle("counts");
595 fOutput->Add(fHistNonMatchedJets2CorrPtArea);
596 }
597
598 fHistJet1PtvsJet2Pt = new TH2F("fHistJet1PtvsJet2Pt", "fHistJet1PtvsJet2Pt", fNbins, fMinBinPt, fMaxBinPt, fNbins, fMinBinPt, fMaxBinPt);
599 fHistJet1PtvsJet2Pt->GetXaxis()->SetTitle("p_{T,1}");
600 fHistJet1PtvsJet2Pt->GetYaxis()->SetTitle("p_{T,2}");
601 fHistJet1PtvsJet2Pt->GetZaxis()->SetTitle("counts");
602 fOutput->Add(fHistJet1PtvsJet2Pt);
603
604 if (!fRhoName.IsNull() || !fRho2Name.IsNull()) {
605 if (fRhoName.IsNull())
606 fHistJet1CorrPtvsJet2CorrPt = new TH2F("fHistJet1CorrPtvsJet2CorrPt", "fHistJet1CorrPtvsJet2CorrPt", fNbins, fMinBinPt, fMaxBinPt, 2*fNbins, -fMaxBinPt, fMaxBinPt);
607 else if (fRho2Name.IsNull())
608 fHistJet1CorrPtvsJet2CorrPt = new TH2F("fHistJet1CorrPtvsJet2CorrPt", "fHistJet1CorrPtvsJet2CorrPt", 2*fNbins, -fMaxBinPt, fMaxBinPt, fNbins, fMinBinPt, fMaxBinPt);
609 else
610 fHistJet1CorrPtvsJet2CorrPt = new TH2F("fHistJet1CorrPtvsJet2CorrPt", "fHistJet1CorrPtvsJet2CorrPt", 2*fNbins, -fMaxBinPt, fMaxBinPt, 2*fNbins, -fMaxBinPt, fMaxBinPt);
611 fHistJet1CorrPtvsJet2CorrPt->GetXaxis()->SetTitle("p_{T,1}^{corr}");
612 fHistJet1CorrPtvsJet2CorrPt->GetYaxis()->SetTitle("p_{T,2}^{corr}");
613 fHistJet1CorrPtvsJet2CorrPt->GetZaxis()->SetTitle("counts");
614 fOutput->Add(fHistJet1CorrPtvsJet2CorrPt);
615 }
616
617 fHistMissedJets2PtArea = new TH2F("fHistMissedJets2PtArea", "fHistMissedJets2PtArea", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, fNbins, fMinBinPt, fMaxBinPt);
618 fHistMissedJets2PtArea->GetXaxis()->SetTitle("area");
619 fHistMissedJets2PtArea->GetYaxis()->SetTitle("p_{T,2} (GeV/c)");
620 fHistMissedJets2PtArea->GetZaxis()->SetTitle("counts");
621 fOutput->Add(fHistMissedJets2PtArea);
2bddb6ae 622
99c67c1b 623 PostData(1, fOutput); // Post data for ALL output slots > 0 here, to get at least an empty histogram
2949a2ec 624}
625
ceefbfbc 626//________________________________________________________________________
a487deae 627Bool_t AliJetResponseMaker::AcceptJet(AliEmcalJet *jet) const
ceefbfbc 628{
629 // Return true if jet is accepted.
630
631 if (jet->Pt() <= fJetPtCut)
632 return kFALSE;
633 if (jet->Area() <= fJetAreaCut)
634 return kFALSE;
65bb5510 635
636 return kTRUE;
ceefbfbc 637}
638
cd6431de 639//________________________________________________________________________
640Bool_t AliJetResponseMaker::AcceptBiasJet2(AliEmcalJet *jet) const
641{
642 // Accept jet with a bias.
643
644 if (fLeadingHadronType == 0) {
645 if (jet->MaxTrackPt() < fPtBiasJet2Track) return kFALSE;
646 }
647 else if (fLeadingHadronType == 1) {
648 if (jet->MaxClusterPt() < fPtBiasJet2Clus) return kFALSE;
649 }
650 else {
651 if (jet->MaxTrackPt() < fPtBiasJet2Track && jet->MaxClusterPt() < fPtBiasJet2Clus) return kFALSE;
652 }
653
654 return kTRUE;
655}
656
b3ccdfe2 657//________________________________________________________________________
658void AliJetResponseMaker::ExecOnce()
659{
660 // Execute once.
661
cd6431de 662 if (!fJets2Name.IsNull() && !fJets2) {
663 fJets2 = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fJets2Name));
664 if (!fJets2) {
665 AliError(Form("%s: Could not retrieve jets2 %s!", GetName(), fJets2Name.Data()));
b3ccdfe2 666 return;
667 }
cd6431de 668 else if (!fJets2->GetClass()->GetBaseClass("AliEmcalJet")) {
669 AliError(Form("%s: Collection %s does not contain AliEmcalJet objects!", GetName(), fJets2Name.Data()));
670 fJets2 = 0;
b3ccdfe2 671 return;
672 }
673 }
674
cd6431de 675 if (!fTracks2Name.IsNull() && !fTracks2) {
676 fTracks2 = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fTracks2Name));
677 if (!fTracks2) {
678 AliError(Form("%s: Could not retrieve tracks2 %s!", GetName(), fTracks2Name.Data()));
b3ccdfe2 679 return;
680 }
681 else {
cd6431de 682 TClass *cl = fTracks2->GetClass();
b3ccdfe2 683 if (!cl->GetBaseClass("AliVParticle") && !cl->GetBaseClass("AliEmcalParticle")) {
cd6431de 684 AliError(Form("%s: Collection %s does not contain AliVParticle nor AliEmcalParticle objects!", GetName(), fTracks2Name.Data()));
685 fTracks2 = 0;
b3ccdfe2 686 return;
687 }
688 }
cfc2ac24 689
690 if (fAreCollections2MC) {
5be3857d 691 fTracks2Map = dynamic_cast<AliNamedArrayI*>(InputEvent()->FindListObject(fTracks2Name + "_Map"));
7f76e479 692 // this is needed to map the MC labels with the indexes of the MC particle collection
693 // if teh map is not given, the MC labels are assumed to be consistent with the indexes (which is not the case if AliEmcalMCTrackSelector is used)
cfc2ac24 694 if (!fTracks2Map) {
7f76e479 695 AliWarning(Form("%s: Could not retrieve map for tracks2 %s! Will assume MC labels consistent with indexes...", GetName(), fTracks2Name.Data()));
5be3857d 696 fTracks2Map = new AliNamedArrayI("tracksMap",9999);
7f76e479 697 for (Int_t i = 0; i < 9999; i++) {
5be3857d 698 fTracks2Map->AddAt(i,i);
7f76e479 699 }
cfc2ac24 700 }
701 }
b3ccdfe2 702 }
703
cd6431de 704 if (!fCalo2Name.IsNull() && !fCaloClusters2) {
705 fCaloClusters2 = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fCalo2Name));
706 if (!fCaloClusters2) {
707 AliError(Form("%s: Could not retrieve clusters %s!", GetName(), fCalo2Name.Data()));
708 return;
709 } else {
710 TClass *cl = fCaloClusters2->GetClass();
711 if (!cl->GetBaseClass("AliVCluster") && !cl->GetBaseClass("AliEmcalParticle")) {
712 AliError(Form("%s: Collection %s does not contain AliVCluster nor AliEmcalParticle objects!", GetName(), fCalo2Name.Data()));
713 fCaloClusters2 = 0;
714 return;
715 }
716 }
91bee8bc 717 }
cd6431de 718
719 if (!fRho2Name.IsNull() && !fRho2) {
720 fRho2 = dynamic_cast<AliRhoParameter*>(InputEvent()->FindListObject(fRho2Name));
721 if (!fRho2) {
722 AliError(Form("%s: Could not retrieve rho %s!", GetName(), fRho2Name.Data()));
723 fInitialized = kFALSE;
724 return;
725 }
91bee8bc 726 }
cd6431de 727
728 if (fJet2MinEta == -999)
729 fJet2MinEta = fJetMinEta - fJetRadius;
730 if (fJet2MaxEta == -999)
731 fJet2MaxEta = fJetMaxEta + fJetRadius;
732 if (fJet2MinPhi == -999)
733 fJet2MinPhi = fJetMinPhi - fJetRadius;
734 if (fJet2MaxPhi == -999)
735 fJet2MaxPhi = fJetMaxPhi + fJetRadius;
736
737 AliAnalysisTaskEmcalJet::ExecOnce();
b3ccdfe2 738}
739
4643d2e8 740//________________________________________________________________________
741Bool_t AliJetResponseMaker::IsEventSelected()
742{
743 // Check if event is selected
744
745 if (fSelectPtHardBin != -999 && fSelectPtHardBin != fPtHardBin)
746 return kFALSE;
747
748 return AliAnalysisTaskEmcalJet::IsEventSelected();
749}
750
b3ccdfe2 751//________________________________________________________________________
752Bool_t AliJetResponseMaker::RetrieveEventObjects()
753{
754 // Retrieve event objects.
755
756 if (!AliAnalysisTaskEmcalJet::RetrieveEventObjects())
757 return kFALSE;
cd6431de 758
759 if (fRho2)
760 fRho2Val = fRho2->GetVal();
b3ccdfe2 761
b3ccdfe2 762 const Int_t ptHardLo[11] = { 0, 5,11,21,36,57, 84,117,152,191,234};
763 const Int_t ptHardHi[11] = { 5,11,21,36,57,84,117,152,191,234,1000000};
787a3c4f 764
765 if (MCEvent())
766 fPythiaHeader = dynamic_cast<AliGenPythiaEventHeader*>(MCEvent()->GenEventHeader());
767
768 if (fPythiaHeader) {
769 Double_t pthard = fPythiaHeader->GetPtHard();
770
771 for (fPtHardBin = 0; fPtHardBin < 11; fPtHardBin++) {
772 if (pthard >= ptHardLo[fPtHardBin] && pthard < ptHardHi[fPtHardBin])
773 break;
774 }
775
776 fNTrials = fPythiaHeader->Trials();
b3ccdfe2 777 }
b3ccdfe2 778
b3ccdfe2 779 return kTRUE;
780}
781
782//________________________________________________________________________
783Bool_t AliJetResponseMaker::Run()
784{
785 // Find the closest jets
cfc2ac24 786
f660c2d6 787 if (fMatching == kNoMatching)
aa4d701c 788 return kTRUE;
f660c2d6 789 else
790 return DoJetMatching();
cfc2ac24 791}
aa4d701c 792
cfc2ac24 793//________________________________________________________________________
f660c2d6 794Bool_t AliJetResponseMaker::DoJetMatching()
cfc2ac24 795{
796 DoJetLoop(kFALSE);
cfc2ac24 797
f660c2d6 798 const Int_t nJets = fJets->GetEntriesFast();
cfc2ac24 799
f660c2d6 800 for (Int_t i = 0; i < nJets; i++) {
cfc2ac24 801
802 AliEmcalJet* jet1 = static_cast<AliEmcalJet*>(fJets->At(i));
803
804 if (!jet1) {
805 AliError(Form("Could not receive jet %d", i));
806 continue;
807 }
808
809 if (!AcceptJet(jet1))
810 continue;
811
812 if (jet1->Eta() < fJetMinEta || jet1->Eta() > fJetMaxEta || jet1->Phi() < fJetMinPhi || jet1->Phi() > fJetMaxPhi)
813 continue;
814
f660c2d6 815 if (jet1->ClosestJet() && jet1->ClosestJet()->ClosestJet() == jet1 &&
7f76e479 816 jet1->ClosestJetDistance() < fMatchingPar1 && jet1->ClosestJet()->ClosestJetDistance() < fMatchingPar2) { // Matched jet found
cfc2ac24 817 jet1->SetMatchedToClosest(fMatching);
cfc2ac24 818 jet1->ClosestJet()->SetMatchedToClosest(fMatching);
aa4d701c 819 }
820 }
821
b3ccdfe2 822 return kTRUE;
823}
824
2949a2ec 825//________________________________________________________________________
cfc2ac24 826void AliJetResponseMaker::DoJetLoop(Bool_t order)
2949a2ec 827{
44476be7 828 // Do the jet loop.
1ee1b5b8 829
cfc2ac24 830 TClonesArray *jets1 = 0;
831 TClonesArray *jets2 = 0;
832
833 if (order) {
834 jets1 = fJets2;
835 jets2 = fJets;
836 }
837 else {
838 jets1 = fJets;
839 jets2 = fJets2;
840 }
841
44476be7 842 Int_t nJets1 = jets1->GetEntriesFast();
843 Int_t nJets2 = jets2->GetEntriesFast();
1ee1b5b8 844
44476be7 845 for (Int_t i = 0; i < nJets1; i++) {
1ee1b5b8 846
44476be7 847 AliEmcalJet* jet1 = static_cast<AliEmcalJet*>(jets1->At(i));
2949a2ec 848
44476be7 849 if (!jet1) {
850 AliError(Form("Could not receive jet %d", i));
851 continue;
852 }
853
ceefbfbc 854 if (!AcceptJet(jet1))
44476be7 855 continue;
856
cfc2ac24 857 if (order) {
cd6431de 858 if (jet1->Eta() < fJet2MinEta || jet1->Eta() > fJet2MaxEta || jet1->Phi() < fJet2MinPhi || jet1->Phi() > fJet2MaxPhi)
ceefbfbc 859 continue;
860 }
861 else {
cd6431de 862 if (jet1->Eta() < fJetMinEta || jet1->Eta() > fJetMaxEta || jet1->Phi() < fJetMinPhi || jet1->Phi() > fJetMaxPhi)
ceefbfbc 863 continue;
864 }
865
44476be7 866 for (Int_t j = 0; j < nJets2; j++) {
867
868 AliEmcalJet* jet2 = static_cast<AliEmcalJet*>(jets2->At(j));
869
870 if (!jet2) {
871 AliError(Form("Could not receive jet %d", j));
872 continue;
873 }
874
ceefbfbc 875 if (!AcceptJet(jet2))
44476be7 876 continue;
ceefbfbc 877
cfc2ac24 878 if (order) {
a487deae 879 if (jet2->Eta() < fJetMinEta || jet2->Eta() > fJetMaxEta || jet2->Phi() < fJetMinPhi || jet2->Phi() > fJetMaxPhi)
ceefbfbc 880 continue;
881 }
882 else {
cd6431de 883 if (jet1->Eta() < fJet2MinEta || jet1->Eta() > fJet2MaxEta || jet1->Phi() < fJet2MinPhi || jet1->Phi() > fJet2MaxPhi)
ceefbfbc 884 continue;
885 }
cd6431de 886
7f76e479 887 SetMatchingLevel(jet1, jet2, fMatching);
44476be7 888 }
889 }
1ee1b5b8 890}
891
cd6431de 892//________________________________________________________________________
7f76e479 893void AliJetResponseMaker::GetGeometricalMatchingLevel(AliEmcalJet *jet1, AliEmcalJet *jet2, Double_t &d) const
cd6431de 894{
7f76e479 895 Double_t deta = jet2->Eta() - jet1->Eta();
896 Double_t dphi = jet2->Phi() - jet1->Phi();
897 d = TMath::Sqrt(deta * deta + dphi * dphi);
898}
cd6431de 899
7f76e479 900//________________________________________________________________________
901void AliJetResponseMaker::GetMCLabelMatchingLevel(AliEmcalJet *jet1, AliEmcalJet *jet2, Double_t &d1, Double_t &d2) const
902{
903 // d1 and d2 represent the matching level: 0 = maximum level of matching, 1 = the two jets are completely unrelated
904 d1 = jet1->Pt();
905 d2 = jet2->Pt();
906 Double_t totalPt1 = d1; // the total pt of the reconstructed jet will be cleaned from the background
907
908 for (Int_t iTrack2 = 0; iTrack2 < jet2->GetNumberOfTracks(); iTrack2++) {
909 Bool_t track2Found = kFALSE;
910 Int_t index2 = jet2->TrackAt(iTrack2);
911 for (Int_t iTrack = 0; iTrack < jet1->GetNumberOfTracks(); iTrack++) {
912 AliVParticle *track = jet1->TrackAt(iTrack,fTracks);
913 if (!track) {
914 AliWarning(Form("Could not find track %d!", iTrack));
915 continue;
916 }
787a3c4f 917 Int_t MClabel = TMath::Abs(track->GetLabel());
7f76e479 918 Int_t index = -1;
919
787a3c4f 920 if (MClabel == 0) {// this is not a MC particle; remove it completely
7f76e479 921 AliDebug(3,Form("Track %d (pT = %f) is not a MC particle (MClabel = %d)!",iTrack,track->Pt(),MClabel));
922 totalPt1 -= track->Pt();
923 d1 -= track->Pt();
924 continue;
925 }
5be3857d 926 else if (MClabel < fTracks2Map->GetSize()) {
927 index = fTracks2Map->At(MClabel);
7f76e479 928 }
929
930 if (index < 0) {
931 AliDebug(2,Form("Track %d (pT = %f) does not have an associated MC particle (MClabel = %d)!",iTrack,track->Pt(),MClabel));
932 continue;
933 }
934
935 if (index2 == index) { // found common particle
936 track2Found = kTRUE;
937 d1 -= track->Pt();
938 AliVParticle *MCpart = static_cast<AliVParticle*>(fTracks2->At(index2));
939 AliDebug(3,Form("Track %d (pT = %f, eta = %f, phi = %f) is associated with the MC particle %d (pT = %f, eta = %f, phi = %f)!",
940 iTrack,track->Pt(),track->Eta(),track->Phi(),MClabel,MCpart->Pt(),MCpart->Eta(),MCpart->Phi()));
941 d2 -= MCpart->Pt();
942 break;
943 }
cd6431de 944 }
7f76e479 945 for (Int_t iClus = 0; iClus < jet1->GetNumberOfClusters(); iClus++) {
946 AliVCluster *clus = jet1->ClusterAt(iClus,fCaloClusters);
947 if (!clus) {
948 AliWarning(Form("Could not find cluster %d!", iClus));
949 continue;
950 }
951 TLorentzVector part;
952 clus->GetMomentum(part, const_cast<Double_t*>(fVertex));
953
954 if (fCaloCells) { // if the cell colection is available, look for cells with a matched MC particle
955 for (Int_t iCell = 0; iCell < clus->GetNCells(); iCell++) {
956 Int_t cellId = clus->GetCellAbsId(iCell);
957 Double_t cellFrac = clus->GetCellAmplitudeFraction(iCell);
958
787a3c4f 959 Int_t MClabel = TMath::Abs(fCaloCells->GetCellMCLabel(cellId));
7f76e479 960 Int_t index = -1;
961
787a3c4f 962 if (MClabel == 0) {// this is not a MC particle; remove it completely
7f76e479 963 AliDebug(3,Form("Cell %d (frac = %f) is not a MC particle (MClabel = %d)!",iCell,cellFrac,MClabel));
964 totalPt1 -= part.Pt() * cellFrac;
965 d1 -= part.Pt() * cellFrac;
966 continue;
967 }
5be3857d 968 else if (MClabel < fTracks2Map->GetSize()) {
969 index = fTracks2Map->At(MClabel);
7f76e479 970 }
971
972 if (index < 0) {
973 AliDebug(3,Form("Cell %d (frac = %f) does not have an associated MC particle (MClabel = %d)!",iCell,cellFrac,MClabel));
974 continue;
975 }
976 if (index2 == index) { // found common particle
977 d1 -= part.Pt() * cellFrac;
978
979 if (!track2Found) {// only if it is not already found among charged tracks (charged particles are most likely already found)
980 AliVParticle *MCpart = static_cast<AliVParticle*>(fTracks2->At(index2));
981 AliDebug(3,Form("Cell %d belonging to cluster %d (pT = %f, eta = %f, phi = %f) is associated with the MC particle %d (pT = %f, eta = %f, phi = %f)!",
982 iCell,iClus,part.Pt(),part.Eta(),part.Phi(),MClabel,MCpart->Pt(),MCpart->Eta(),MCpart->Phi()));
983 d2 -= MCpart->Pt() * cellFrac;
984 }
985 break;
986 }
cfc2ac24 987 }
988 }
7f76e479 989 else { //otherwise look for the first contributor to the cluster, and if matched to a MC label remove it
787a3c4f 990 Int_t MClabel = TMath::Abs(clus->GetLabel());
7f76e479 991 Int_t index = -1;
992
787a3c4f 993 if (MClabel == 0) {// this is not a MC particle; remove it completely
7f76e479 994 AliDebug(3,Form("Cluster %d (pT = %f) is not a MC particle (MClabel = %d)!",iClus,part.Pt(),MClabel));
995 totalPt1 -= part.Pt();
996 d1 -= part.Pt();
cfc2ac24 997 continue;
998 }
5be3857d 999 else if (MClabel < fTracks2Map->GetSize()) {
1000 index = fTracks2Map->At(MClabel);
f660c2d6 1001 }
7f76e479 1002
cfc2ac24 1003 if (index < 0) {
7f76e479 1004 AliDebug(3,Form("Cluster %d (pT = %f) does not have an associated MC particle (MClabel = %d)!",iClus,part.Pt(),MClabel));
cfc2ac24 1005 continue;
1006 }
7f76e479 1007 if (index2 == index) { // found common particle
1008 d1 -= part.Pt();
1009
1010 if (!track2Found) {// only if it is not already found among charged tracks (charged particles are most likely already found)
f660c2d6 1011 AliVParticle *MCpart = static_cast<AliVParticle*>(fTracks2->At(index2));
7f76e479 1012 AliDebug(3,Form("Cluster %d (pT = %f, eta = %f, phi = %f) is associated with the MC particle %d (pT = %f, eta = %f, phi = %f)!",
1013 iClus,part.Pt(),part.Eta(),part.Phi(),MClabel,MCpart->Pt(),MCpart->Eta(),MCpart->Phi()));
1014
f660c2d6 1015 d2 -= MCpart->Pt();
cfc2ac24 1016 }
7f76e479 1017 break;
cfc2ac24 1018 }
1019 }
7f76e479 1020 }
1021 }
1022 if (d1 <= 0 || totalPt1 < 1)
1023 d1 = 0;
1024 else
1025 d1 /= totalPt1;
f660c2d6 1026
7f76e479 1027 if (jet2->Pt() > 0 && d2 > 0)
1028 d2 /= jet2->Pt();
1029 else
1030 d2 = 0;
1031}
1032
1033//________________________________________________________________________
1034void AliJetResponseMaker::GetSameCollectionsMatchingLevel(AliEmcalJet *jet1, AliEmcalJet *jet2, Double_t &d1, Double_t &d2) const
1035{
1036 // d1 and d2 represent the matching level: 0 = maximum level of matching, 1 = the two jets are completely unrelated
1037 d1 = jet1->Pt();
1038 d2 = jet2->Pt();
1039
1040 if (fTracks && fTracks2) {
1041
1042 for (Int_t iTrack2 = 0; iTrack2 < jet2->GetNumberOfTracks(); iTrack2++) {
1043 Int_t index2 = jet2->TrackAt(iTrack2);
1044 for (Int_t iTrack = 0; iTrack < jet1->GetNumberOfTracks(); iTrack++) {
1045 Int_t index = jet1->TrackAt(iTrack);
1046 if (index2 == index) { // found common particle
1047 AliVParticle *part = static_cast<AliVParticle*>(fTracks->At(index));
1048 if (!part) {
1049 AliWarning(Form("Could not find track %d!", index));
f660c2d6 1050 continue;
1051 }
7f76e479 1052 AliVParticle *part2 = static_cast<AliVParticle*>(fTracks2->At(index2));
1053 if (!part2) {
1054 AliWarning(Form("Could not find track %d!", index2));
cfc2ac24 1055 continue;
1056 }
7f76e479 1057
1058 d1 -= part->Pt();
1059 d2 -= part2->Pt();
1060 break;
cfc2ac24 1061 }
1062 }
cfc2ac24 1063 }
7f76e479 1064
1065 }
1066
1067 if (fCaloClusters && fCaloClusters2) {
1068
1069 for (Int_t iClus2 = 0; iClus2 < jet2->GetNumberOfClusters(); iClus2++) {
1070 Int_t index2 = jet2->ClusterAt(iClus2);
1071 for (Int_t iClus = 0; iClus < jet1->GetNumberOfClusters(); iClus++) {
1072 Int_t index = jet1->ClusterAt(iClus);
1073 AliVCluster *clus = static_cast<AliVCluster*>(fCaloClusters->At(index));
1074 if (!clus) {
1075 AliWarning(Form("Could not find cluster %d!", index));
1076 continue;
1077 }
1078 AliVCluster *clus2 = static_cast<AliVCluster*>(fCaloClusters2->At(index2));
1079 if (!clus2) {
1080 AliWarning(Form("Could not find cluster %d!", index2));
1081 continue;
1082 }
1083 TLorentzVector part, part2;
1084 clus->GetMomentum(part, const_cast<Double_t*>(fVertex));
1085 clus2->GetMomentum(part2, const_cast<Double_t*>(fVertex));
1086
1087 d1 -= part.Pt();
1088 d2 -= part2.Pt();
1089 break;
1090 }
1091 }
1092
1093 }
1094
1095 if (jet1->Pt() > 0 && d1 > 0)
1096 d1 /= jet1->Pt();
1097 else
1098 d1 = 0;
1099
1100 if (jet2->Pt() > 0 && d2 > 0)
1101 d2 /= jet2->Pt();
1102 else
1103 d2 = 0;
1104}
1105
1106//________________________________________________________________________
1107void AliJetResponseMaker::SetMatchingLevel(AliEmcalJet *jet1, AliEmcalJet *jet2, MatchingType matching)
1108{
1109 Double_t d1 = -1;
1110 Double_t d2 = -1;
1111
1112 switch (matching) {
1113 case kGeometrical:
1114 GetGeometricalMatchingLevel(jet1,jet2,d1);
1115 d2 = d1;
1116 break;
1117 case kMCLabel: // jet1 = detector level and jet2 = particle level!
1118 GetMCLabelMatchingLevel(jet1,jet2,d1,d2);
1119 break;
1120 case kSameCollections:
1121 GetSameCollectionsMatchingLevel(jet1,jet2,d1,d2);
cd6431de 1122 break;
1123 default:
1124 ;
1125 }
f660c2d6 1126
1127 if (d1 > 0) {
1128
1129 if (d1 < jet1->ClosestJetDistance()) {
1130 jet1->SetSecondClosestJet(jet1->ClosestJet(), jet1->ClosestJetDistance());
1131 jet1->SetClosestJet(jet2, d1);
1132 }
1133 else if (d1 < jet1->SecondClosestJetDistance()) {
1134 jet1->SetSecondClosestJet(jet2, d1);
1135 }
1136 }
1137
1138 if (d2 > 0) {
1139
1140 if (d2 < jet2->ClosestJetDistance()) {
1141 jet2->SetSecondClosestJet(jet2->ClosestJet(), jet2->ClosestJetDistance());
1142 jet2->SetClosestJet(jet1, d2);
1143 }
1144 else if (d2 < jet2->SecondClosestJetDistance()) {
1145 jet2->SetSecondClosestJet(jet1, d2);
1146 }
1147 }
cd6431de 1148}
1149
1ee1b5b8 1150//________________________________________________________________________
1151Bool_t AliJetResponseMaker::FillHistograms()
1152{
1153 // Fill histograms.
1154
ca5c29fa 1155 static Int_t indexes[9999] = {-1};
1156
1157 fHistEventsAfterSel->SetBinContent(fPtHardBin + 1, fHistEventsAfterSel->GetBinContent(fPtHardBin + 1) + 1);
1158 fHistTrialsAfterSel->SetBinContent(fPtHardBin + 1, fHistTrialsAfterSel->GetBinContent(fPtHardBin + 1) + fNTrials);
1159
1160 GetSortedArray(indexes, fJets2, fRho2Val);
4643d2e8 1161
cd6431de 1162 const Int_t nJets2 = fJets2->GetEntriesFast();
2949a2ec 1163
ca5c29fa 1164 Int_t naccJets2 = 0;
1165 Int_t naccJets2Acceptance = 0;
1166
cd6431de 1167 for (Int_t i = 0; i < nJets2; i++) {
2949a2ec 1168
ca5c29fa 1169 AliEmcalJet* jet2 = static_cast<AliEmcalJet*>(fJets2->At(indexes[i]));
2949a2ec 1170
cd6431de 1171 if (!jet2) {
1172 AliError(Form("Could not receive jet2 %d", i));
2949a2ec 1173 continue;
cfc2ac24 1174 }
1175
cd6431de 1176 if (!AcceptJet(jet2))
1177 continue;
1178
cfc2ac24 1179 if (AcceptBiasJet(jet2) &&
1180 (jet2->Eta() > fJetMinEta && jet2->Eta() < fJetMaxEta && jet2->Phi() > fJetMinPhi && jet2->Phi() < fJetMaxPhi)) {
1181
1182 fHistJets2PtAreaAcceptance->Fill(jet2->Area(), jet2->Pt());
1183 fHistJets2PhiEtaAcceptance->Fill(jet2->Eta(), jet2->Phi());
1184
ca5c29fa 1185 if (naccJets2Acceptance < fNLeadingJets)
1186 fHistLeadingJets2PtAreaAcceptance->Fill(jet2->Area(), jet2->Pt());
1187
1188 if (!fRho2Name.IsNull()) {
cfc2ac24 1189 fHistJets2CorrPtAreaAcceptance->Fill(jet2->Area(), jet2->Pt() - fRho2Val * jet2->Area());
ca5c29fa 1190 if (naccJets2Acceptance < fNLeadingJets)
1191 fHistLeadingJets2CorrPtAreaAcceptance->Fill(jet2->Area(), jet2->Pt() - fRho2Val * jet2->Area());
1192 }
63fac07f 1193
1194 if (fTracks2) {
1195 for (Int_t it = 0; it < jet2->GetNumberOfTracks(); it++) {
1196 AliVParticle *track2 = jet2->TrackAt(it, fTracks2);
1197 if (track2)
1198 fHistJets2ZvsPt->Fill(track2->Pt() / jet2->Pt(), jet2->Pt());
1199 }
1200 }
1201
1202 if (fCaloClusters2) {
1203 for (Int_t ic = 0; ic < jet2->GetNumberOfClusters(); ic++) {
1204 AliVCluster *cluster2 = jet2->ClusterAt(ic, fCaloClusters2);
1205
1206 if (cluster2) {
1207 TLorentzVector nPart2;
1208 cluster2->GetMomentum(nPart2, fVertex);
1209 fHistJets2ZvsPt->Fill(nPart2.Et() / jet2->Pt(), jet2->Pt());
1210 }
1211 }
1212 }
1213
1214 fHistJets2NEFvsPt->Fill(jet2->NEF(), jet2->Pt());
1215 fHistJets2CEFvsCEFPt->Fill(1-jet2->NEF(), (1-jet2->NEF())*jet2->Pt());
ca5c29fa 1216
1217 naccJets2Acceptance++;
cfc2ac24 1218 }
1219
cd6431de 1220 if (!AcceptBiasJet2(jet2))
ceefbfbc 1221 continue;
1222
cd6431de 1223 if (jet2->Eta() < fJet2MinEta || jet2->Eta() > fJet2MaxEta || jet2->Phi() < fJet2MinPhi || jet2->Phi() > fJet2MaxPhi)
2bddb6ae 1224 continue;
ca5c29fa 1225
cfc2ac24 1226 fHistJets2PtArea->Fill(jet2->Area(), jet2->Pt());
1227 fHistJets2PhiEta->Fill(jet2->Eta(), jet2->Phi());
ca5c29fa 1228
1229 if (naccJets2 < fNLeadingJets)
1230 fHistLeadingJets2PtArea->Fill(jet2->Area(), jet2->Pt());
cfc2ac24 1231
ca5c29fa 1232 if (!fRho2Name.IsNull()) {
cfc2ac24 1233 fHistJets2CorrPtArea->Fill(jet2->Area(), jet2->Pt() - fRho2Val * jet2->Area());
ca5c29fa 1234 if (naccJets2 < fNLeadingJets)
1235 fHistLeadingJets2CorrPtArea->Fill(jet2->Area(), jet2->Pt() - fRho2Val * jet2->Area());
1236 }
1237
1238 naccJets2++;
2949a2ec 1239
cd6431de 1240 if (jet2->MatchedJet()) {
aa4d701c 1241
cd6431de 1242 if (!AcceptBiasJet(jet2->MatchedJet()) ||
ca5c29fa 1243 jet2->MatchedJet()->MaxTrackPt() > fMaxTrackPt || jet2->MatchedJet()->MaxClusterPt() > fMaxClusterPt) {
cd6431de 1244 fHistMissedJets2PtArea->Fill(jet2->Area(), jet2->Pt());
2bddb6ae 1245 }
1246 else {
ca5c29fa 1247 if (jet2->MatchedJet()->Pt() > fMaxBinPt)
1248 fHistMissedJets2PtArea->Fill(jet2->Area(), jet2->Pt());
1249
7f76e479 1250 Double_t d1=-1, d2=-1;
1251 if (jet2->GetMatchingType() == kGeometrical) {
1252
1253 if (fAreCollections2MC && !fAreCollections1MC)
1254 GetMCLabelMatchingLevel(jet2->MatchedJet(), jet2, d1, d2);
1255 else if ((fAreCollections1MC && fAreCollections2MC) || (!fAreCollections1MC && !fAreCollections2MC))
1256 GetSameCollectionsMatchingLevel(jet2->MatchedJet(), jet2, d1, d2);
1257
6a20534a 1258 fHistDistancevsCommonEnergy1->Fill(jet2->ClosestJetDistance(), d1);
1259 fHistDistancevsCommonEnergy2->Fill(jet2->ClosestJetDistance(), d2);
ca5c29fa 1260
1261 fHistDistancevsJet1Pt->Fill(jet2->ClosestJetDistance(), jet2->MatchedJet()->Pt());
1262 fHistDistancevsJet2Pt->Fill(jet2->ClosestJetDistance(), jet2->Pt());
1263
1264 fHistCommonEnergy1vsJet1Pt->Fill(d1, jet2->MatchedJet()->Pt());
1265 fHistCommonEnergy2vsJet2Pt->Fill(d2, jet2->Pt());
7f76e479 1266 }
1267 else if (jet2->GetMatchingType() == kMCLabel || jet2->GetMatchingType() == kSameCollections) {
1268 GetGeometricalMatchingLevel(jet2->MatchedJet(), jet2, d1);
ca5c29fa 1269
6a20534a 1270 fHistDistancevsCommonEnergy1->Fill(d1, jet2->MatchedJet()->ClosestJetDistance());
1271 fHistDistancevsCommonEnergy2->Fill(d1, jet2->ClosestJetDistance());
ca5c29fa 1272
1273 fHistDistancevsJet1Pt->Fill(d1, jet2->MatchedJet()->Pt());
1274 fHistDistancevsJet2Pt->Fill(d1, jet2->Pt());
1275
1276 fHistCommonEnergy1vsJet1Pt->Fill(jet2->MatchedJet()->ClosestJetDistance(), jet2->MatchedJet()->Pt());
1277 fHistCommonEnergy2vsJet2Pt->Fill(jet2->ClosestJetDistance(), jet2->Pt());
7f76e479 1278 }
cd6431de 1279
1280 Double_t deta = jet2->MatchedJet()->Eta() - jet2->Eta();
1281 Double_t dphi = jet2->MatchedJet()->Phi() - jet2->Phi();
ca5c29fa 1282 fHistDeltaEtaPhi->Fill(deta, dphi, jet2->Pt());
cd6431de 1283
1284 Double_t dpt = jet2->MatchedJet()->Pt() - jet2->Pt();
7030f36f 1285 fHistDeltaPtvsJet1Pt->Fill(jet2->MatchedJet()->Pt(), dpt);
cfc2ac24 1286 fHistDeltaPtvsJet2Pt->Fill(jet2->Pt(), dpt);
1287 fHistDeltaPtvsMatchingLevel->Fill(jet2->ClosestJetDistance(), dpt);
cd6431de 1288
7030f36f 1289 fHistJet2PtOverJet1PtvsJet2Pt->Fill(jet2->Pt(), jet2->Pt() / jet2->MatchedJet()->Pt());
1290 fHistJet1PtOverJet2PtvsJet1Pt->Fill(jet2->MatchedJet()->Pt(), jet2->MatchedJet()->Pt() / jet2->Pt());
c560b734 1291
cd6431de 1292 fHistJet1PtvsJet2Pt->Fill(jet2->MatchedJet()->Pt(), jet2->Pt());
1293
1294 if (!fRhoName.IsNull() || !fRho2Name.IsNull()) {
1295 dpt -= fRhoVal * jet2->MatchedJet()->Area() - fRho2Val * jet2->Area();
7030f36f 1296 fHistDeltaCorrPtvsJet1Pt->Fill(jet2->MatchedJet()->Pt(), dpt);
cfc2ac24 1297 fHistDeltaCorrPtvsJet2Pt->Fill(jet2->Pt(), dpt);
1298 fHistDeltaCorrPtvsMatchingLevel->Fill(jet2->ClosestJetDistance(), dpt);
cd6431de 1299 fHistJet1CorrPtvsJet2CorrPt->Fill(jet2->MatchedJet()->Pt() - fRhoVal * jet2->MatchedJet()->Area(), jet2->Pt() - fRho2Val * jet2->Area());
1300 }
2bddb6ae 1301 }
2949a2ec 1302 }
99c67c1b 1303 else {
cd6431de 1304 fHistNonMatchedJets2PtArea->Fill(jet2->Area(), jet2->Pt());
1305 fHistMissedJets2PtArea->Fill(jet2->Area(), jet2->Pt());
2949a2ec 1306
cd6431de 1307 if (!fRho2Name.IsNull())
1308 fHistNonMatchedJets2CorrPtArea->Fill(jet2->Area(), jet2->Pt() - fRhoVal * jet2->Area());
1b76c28f 1309 }
1310 }
1311
ca5c29fa 1312 GetSortedArray(indexes, fJets, fRhoVal);
1313
cd6431de 1314 const Int_t nJets1 = fJets->GetEntriesFast();
ca5c29fa 1315 Int_t naccJets1 = 0;
1b76c28f 1316
cd6431de 1317 for (Int_t i = 0; i < nJets1; i++) {
1b76c28f 1318
ca5c29fa 1319 AliDebug(2,Form("Processing jet %d", i));
1320 AliEmcalJet* jet1 = static_cast<AliEmcalJet*>(fJets->At(indexes[i]));
1b76c28f 1321
cd6431de 1322 if (!jet1) {
1323 AliError(Form("Could not receive jet %d", i));
2949a2ec 1324 continue;
1b76c28f 1325 }
99c67c1b 1326
cd6431de 1327 if (!AcceptJet(jet1))
1b76c28f 1328 continue;
cd6431de 1329
1330 if (!AcceptBiasJet(jet1))
91bee8bc 1331 continue;
cd6431de 1332
1333 if (jet1->MaxTrackPt() > fMaxTrackPt || jet1->MaxClusterPt() > fMaxClusterPt)
91bee8bc 1334 continue;
cd6431de 1335
1336 if (jet1->Eta() < fJetMinEta || jet1->Eta() > fJetMaxEta || jet1->Phi() < fJetMinPhi || jet1->Phi() > fJetMaxPhi)
91bee8bc 1337 continue;
2949a2ec 1338
cd6431de 1339 if (!jet1->MatchedJet()) {
1340 fHistNonMatchedJets1PtArea->Fill(jet1->Area(), jet1->Pt());
1341 if (!fRhoName.IsNull())
1342 fHistNonMatchedJets1CorrPtArea->Fill(jet1->Area(), jet1->Pt() - fRhoVal * jet1->Area());
99c67c1b 1343 }
1b76c28f 1344
cd6431de 1345 fHistJets1PtArea->Fill(jet1->Area(), jet1->Pt());
1346 fHistJets1PhiEta->Fill(jet1->Eta(), jet1->Phi());
1b76c28f 1347
ca5c29fa 1348 if (naccJets1 < fNLeadingJets)
1349 fHistLeadingJets1PtArea->Fill(jet1->Area(), jet1->Pt());
1350
1351 if (!fRhoName.IsNull()) {
cd6431de 1352 fHistJets1CorrPtArea->Fill(jet1->Area(), jet1->Pt() - fRhoVal * jet1->Area());
ca5c29fa 1353
1354 if (naccJets1 < fNLeadingJets)
1355 fHistLeadingJets1CorrPtArea->Fill(jet1->Area(), jet1->Pt() - fRhoVal * jet1->Area());
1356 }
1357
63fac07f 1358 if (fTracks) {
1359 for (Int_t it = 0; it < jet1->GetNumberOfTracks(); it++) {
1360 AliVParticle *track1 = jet1->TrackAt(it, fTracks2);
1361 if (track1)
1362 fHistJets1ZvsPt->Fill(track1->Pt() / jet1->Pt(), jet1->Pt());
1363 }
1364 }
1365
1366 if (fCaloClusters) {
1367 for (Int_t ic = 0; ic < jet1->GetNumberOfClusters(); ic++) {
1368 AliVCluster *cluster1 = jet1->ClusterAt(ic, fCaloClusters);
1369
1370 if (cluster1) {
1371 TLorentzVector nPart1;
1372 cluster1->GetMomentum(nPart1, fVertex);
1373 fHistJets2ZvsPt->Fill(nPart1.Et() / jet1->Pt(), jet1->Pt());
1374 }
1375 }
1376 }
1377
1378 fHistJets1NEFvsPt->Fill(jet1->NEF(), jet1->Pt());
1379 fHistJets1CEFvsCEFPt->Fill(1-jet1->NEF(), (1-jet1->NEF())*jet1->Pt());
1380
ca5c29fa 1381 naccJets1++;
1b76c28f 1382 }
6fd5039f 1383
1384 return kTRUE;
1b76c28f 1385}