Charged jets (pPb): Improved trackcut analysis
[u/mrichter/AliRoot.git] / PWGJE / EMCALJetTasks / UserTasks / AliAnalysisTaskEmcalDiJetResponse.cxx
1 //
2 // Dijet response analysis task.
3 //
4 // Author: M.Verweij
5
6 #include <TClonesArray.h>
7 #include <TH1F.h>
8 #include <TH2F.h>
9 #include <TH3F.h>
10 #include <THnSparse.h>
11 #include <TList.h>
12 #include <TLorentzVector.h>
13 #include <TProfile.h>
14 #include <TChain.h>
15 #include <TSystem.h>
16 #include <TFile.h>
17 #include <TKey.h>
18
19 #include "AliVCluster.h"
20 #include "AliVTrack.h"
21 #include "AliEmcalJet.h"
22 #include "AliRhoParameter.h"
23 #include "AliLog.h"
24 #include "AliEmcalParticle.h"
25 #include "AliMCEvent.h"
26 #include "AliGenPythiaEventHeader.h"
27 #include "AliAODMCHeader.h"
28 #include "AliMCEvent.h"
29 #include "AliAnalysisManager.h"
30 #include "AliJetContainer.h"
31
32 #include "AliAnalysisTaskEmcalDiJetResponse.h"
33
34 ClassImp(AliAnalysisTaskEmcalDiJetResponse)
35
36 //________________________________________________________________________
37 AliAnalysisTaskEmcalDiJetResponse::AliAnalysisTaskEmcalDiJetResponse() : 
38   AliAnalysisTaskEmcalDiJetBase("AliAnalysisTaskEmcalDiJetResponse"),
39   fDoMatchFullCharged(kTRUE),
40   fhnDiJetResponseCharged(0),
41   fhnDiJetResponseFullCharged(0),
42   fh1TriggersLostCharged(0),
43   fh1TriggersLostFull(0),
44   fh3AssocLostPtDeltaPhiCharged(0),
45   fh3AssocLostPtDeltaPhiFull(0),
46   fhnMatchingCharged(0),
47   fhnMatchingFull(0),
48   fnUsedResponseVar(0)
49 {
50   // Default constructor.
51   
52   SetMakeGeneralHistograms(kTRUE);
53
54   for(Int_t i = 0; i<2; i++) {
55     fh1TriggersCharged[i] = 0;
56     fh1TriggersFull[i] = 0;
57   }
58
59 }
60
61 //________________________________________________________________________
62 AliAnalysisTaskEmcalDiJetResponse::AliAnalysisTaskEmcalDiJetResponse(const char *name) : 
63   AliAnalysisTaskEmcalDiJetBase(name),
64   fDoMatchFullCharged(kTRUE),
65   fhnDiJetResponseCharged(0),
66   fhnDiJetResponseFullCharged(0),
67   fh1TriggersLostCharged(0),
68   fh1TriggersLostFull(0),
69   fh3AssocLostPtDeltaPhiCharged(0),
70   fh3AssocLostPtDeltaPhiFull(0),
71   fhnMatchingCharged(0),
72   fhnMatchingFull(0),
73   fnUsedResponseVar(0)
74 {
75   // Standard constructor.
76
77   SetMakeGeneralHistograms(kTRUE);
78
79   for(Int_t i = 0; i<2; i++) {
80     fh1TriggersCharged[i] = 0;
81     fh1TriggersFull[i] = 0;
82   }
83
84 }
85
86 //________________________________________________________________________
87 AliAnalysisTaskEmcalDiJetResponse::~AliAnalysisTaskEmcalDiJetResponse()
88 {
89   // Destructor.
90 }
91
92
93 //________________________________________________________________________
94 void AliAnalysisTaskEmcalDiJetResponse::UserCreateOutputObjects()
95 {
96   // Create user output.
97
98   AliAnalysisTaskEmcalDiJetBase::UserCreateOutputObjects();
99
100   Bool_t oldStatus = TH1::AddDirectoryStatus();
101   TH1::AddDirectory(kFALSE);
102
103   //Store dijet vars: pt,trig MC, pt,trig DET, pt,ass MC, pt,ass DET, dPhi MC, dPhi Det, kT MC, kT Det 
104   const Int_t nBinsSparse0 = 10;
105   const Int_t nBinsPt = 250;
106   const Int_t nBinsDPhi     = 36;
107   const Int_t nBinsKt       = 25;
108   const Int_t nBinsDiJetEta = 40;
109   const Int_t nBinsAj       = 50;
110   const Int_t nBinsVar[2] = {nBinsKt,nBinsDiJetEta};
111
112   const Int_t nBins0[nBinsSparse0] = {nBinsPt,nBinsPt,nBinsPt,nBinsPt,nBinsDPhi,nBinsDPhi,nBinsVar[fnUsedResponseVar],nBinsVar[fnUsedResponseVar],nBinsAj,nBinsAj};
113
114   const Double_t minPt = 0.;
115   const Double_t maxPt = 250.;
116   const Double_t minVar[2] = {   0.,-1.};
117   const Double_t maxVar[2] = { 100., 1.};
118
119   const Double_t xmin0[nBinsSparse0]  = {  minPt, minPt, minPt, minPt, 0.5*TMath::Pi(), 0.5*TMath::Pi(), minVar[fnUsedResponseVar], minVar[fnUsedResponseVar],0.,0.};
120   const Double_t xmax0[nBinsSparse0]  = {  maxPt, maxPt, maxPt, maxPt, 1.5*TMath::Pi(), 1.5*TMath::Pi(), maxVar[fnUsedResponseVar], maxVar[fnUsedResponseVar],1.,1.};
121
122   fhnDiJetResponseCharged = new THnSparseF("fhnDiJetResponseCharged","fhnDiJetResponseCharged;p_{T,trig}^{part};p_{T,trig}^{det};p_{T,ass}^{part};p_{T,ass}^{det};#Delta#varphi_{part};#Delta#varphi_{det};k_{T}^{part},k_{T}^{det};A_{j}^{part}A_{j}^{det}",nBinsSparse0,nBins0,xmin0,xmax0);
123
124   fhnDiJetResponseFullCharged = new THnSparseF("fhnDiJetResponseFullCharged","fhnDiJetResponseFullCharged;p_{T,trig}^{part};p_{T,trig}^{det};p_{T,ass}^{part};p_{T,ass}^{det};#Delta#varphi_{part};#Delta#varphi_{det};k_{T}^{part},k_{T}^{det};A_{j}^{part}A_{j}^{det}",nBinsSparse0,nBins0,xmin0,xmax0);
125
126   if(fnUsedResponseVar==1) {
127     fhnDiJetResponseCharged->SetTitle("fhnDiJetResponseCharged DiJetEta"); 
128     fhnDiJetResponseCharged->GetAxis(6)->SetTitle("#eta_{dijet}^{part}");
129     fhnDiJetResponseCharged->GetAxis(7)->SetTitle("#eta_{dijet}^{det}");
130
131     fhnDiJetResponseFullCharged->SetTitle("fhnDiJetResponseFullCharged DiJetEta"); 
132     fhnDiJetResponseFullCharged->GetAxis(6)->SetTitle("#eta_{dijet}^{part}");
133     fhnDiJetResponseFullCharged->GetAxis(7)->SetTitle("#eta_{dijet}^{det}");
134   }
135
136   fOutput->Add(fhnDiJetResponseCharged);
137   fOutput->Add(fhnDiJetResponseFullCharged);
138
139   TString strType = "";
140   for(Int_t i = 0; i<2; i++) {
141     if(i==0)      strType="Part";
142     else if(i==1) strType="Det";
143     fh1TriggersCharged[i] = new TH1F(Form("fh1TriggersCharged%s",strType.Data()),Form("fh1TriggersCharged%s;p_{T,trig}^{ch}",strType.Data()),nBinsPt,minPt,maxPt);
144     fOutput->Add(fh1TriggersCharged[i]);
145   
146     fh1TriggersFull[i] = new TH1F(Form("fh1TriggersFull%s",strType.Data()),Form("fh1TriggersFull%s;p_{T,trig}^{ch}",strType.Data()),nBinsPt,minPt,maxPt);
147     fOutput->Add(fh1TriggersFull[i]);
148   }
149
150   fh1TriggersLostCharged = new TH1F("fh1TriggersLostCharged","fh1TriggersLostCharged;p_{T,trig}^{ch}",nBinsPt,minPt,maxPt);
151   fOutput->Add(fh1TriggersLostCharged);
152
153   fh1TriggersLostFull = new TH1F("fh1TriggersLostFull","fh1TriggersLostFull;p_{T,trig}^{ch}",nBinsPt,minPt,maxPt);
154   fOutput->Add(fh1TriggersLostFull);
155
156   fh3AssocLostPtDeltaPhiCharged = new TH3F("fh3AssocLostPtDeltaPhiCharged","fh3AssocLostPtDeltaPhiCharged;p_{T,trig}^{ch};p_{T,assoc}^{ch};#Delta#varphi",nBinsPt,minPt,maxPt,nBinsPt,minPt,maxPt,nBinsDPhi,-0.5*TMath::Pi(),1.5*TMath::Pi());
157   fOutput->Add(fh3AssocLostPtDeltaPhiCharged);
158
159   fh3AssocLostPtDeltaPhiFull = new TH3F("fh3AssocLostPtDeltaPhiFull","fh3AssocLostPtDeltaPhiFull;p_{T,trig}^{ch};p_{T,assoc}^{ch};#Delta#varphi",nBinsPt,minPt,maxPt,nBinsPt,minPt,maxPt,nBinsDPhi,-0.5*TMath::Pi(),1.5*TMath::Pi());
160   fOutput->Add(fh3AssocLostPtDeltaPhiFull);
161
162   const Int_t nBinsSparseMatch = 6;
163   const Int_t nBinsDPhiMatch = 80;
164   const Int_t nBinsDEtaMatch = 80;
165   const Int_t nBinsDR = 20;
166   const Int_t nBinsType = 3;
167   const Int_t nBinsMatch[nBinsSparseMatch] = {nBinsPt,nBinsPt,nBinsDPhiMatch,nBinsDEtaMatch,nBinsDR,nBinsType};
168   //pTpart, pTdet, deltaPhi, deltaEta, deltaR, jet type (leading,subleading,other)
169   const Double_t xminMatch[nBinsSparseMatch]  = { minPt, minPt, -0.5,-0.5, 0., 0};
170   const Double_t xmaxMatch[nBinsSparseMatch]  = { maxPt, maxPt,  0.5, 0.5, 0.5,3};
171   fhnMatchingCharged = new THnSparseF("fhnMatchingCharged","fhnMatchingCharged;#it{p}_{T,part} (GeV/#it{c});#it{p}_{T,det} (GeV/#it{c});#Delta#varphi;#Delta#eta;#Delta R;type",
172                                           nBinsSparseMatch,nBinsMatch,xminMatch,xmaxMatch);
173   fOutput->Add(fhnMatchingCharged);
174
175   fhnMatchingFull = new THnSparseF("fhnMatchingFull","fhnMatchingFull;#it{p}_{T,part} (GeV/#it{c});#it{p}_{T,det} (GeV/#it{c});#Delta#varphi;#Delta#eta;#Delta R;type",
176                                           nBinsSparseMatch,nBinsMatch,xminMatch,xmaxMatch);
177   fOutput->Add(fhnMatchingFull);
178
179
180   // =========== Switch on Sumw2 for all histos ===========
181   for (Int_t i=0; i<fOutput->GetEntries(); ++i) {
182     TH1 *h1 = dynamic_cast<TH1*>(fOutput->At(i));
183     if (h1){
184       h1->Sumw2();
185       continue;
186     }
187     THnSparse *hn = dynamic_cast<THnSparse*>(fOutput->At(i));
188     if(hn)hn->Sumw2();
189   }
190
191   TH1::AddDirectory(oldStatus);
192
193   PostData(1, fOutput); // Post data for ALL output slots > 0 here.
194 }
195
196 //________________________________________________________________________
197 Bool_t AliAnalysisTaskEmcalDiJetResponse::Run()
198 {
199   // Run analysis code here, if needed. It will be executed before FillHistograms().
200       
201   //Check if event is selected (vertex & pile-up)
202   if(!SelectEvent())
203     return kFALSE;
204
205   if(fRhoType==0) {
206     fRhoFullVal = 0.;
207     fRhoChVal = 0.;
208   }
209   if(fRhoType==1) {
210     fRhoFullVal = GetRhoVal(fContainerFull);
211     fRhoChVal = GetRhoVal(fContainerCharged);
212   }
213   
214
215   //Do matching
216   MatchJetsGeo(fContainerCharged,fContainerChargedMC,0,0.3,1);
217   MatchJetsGeo(fContainerFull,fContainerFullMC,0,0.3,2);
218
219   //Fill particle-detector level matching histos
220   
221   if(fDoChargedCharged)   CorrelateJets(1);
222
223   if(fDoFullCharged) {
224     SetChargedFractionIndexMC();
225     CorrelateJets(2);
226   }
227   
228   return kTRUE;  // If return kFALSE FillHistogram() will NOT be executed.
229 }
230
231 //________________________________________________________________________
232 void AliAnalysisTaskEmcalDiJetResponse::CorrelateJets(const Int_t type) {
233   //
234   // Correlate jets and fill histos
235   //
236
237   if( fJetCorrelationType==kCorrelateAll )
238     CorrelateAllJets(type);
239   else if( fJetCorrelationType==kCorrelateTwo )
240     CorrelateTwoJets(type);
241   else if( fJetCorrelationType==kCorrelateLS )
242     AliWarning(Form("%s: leading-subleading correlation not implemented for response!",GetName()));
243
244   return;
245
246 }
247
248
249 //________________________________________________________________________
250 void AliAnalysisTaskEmcalDiJetResponse::CorrelateAllJets(const Int_t type) {
251   //
252   // Correlate jets and fill histos
253   //
254
255   Int_t typet = 0;
256   Int_t typetMC = 0;
257   Int_t typeaMC = 0;
258   if(type==0) { //full-full
259     typetMC = fContainerFullMC;
260     typeaMC = fContainerFullMC;
261     typet = fContainerFull;
262   }
263   else if(type==1) { //charged-charged
264     typetMC = fContainerChargedMC;
265     typeaMC = fContainerChargedMC;
266     typet = fContainerCharged;
267   }
268   else if(type==2) { //full-charged
269     typetMC = fContainerFullMC;
270     typeaMC = fContainerChargedMC;
271     typet = fContainerFull;
272   }
273   else {
274     AliWarning(Form("%s: type %d of dijet correlation not defined!",GetName(),type));
275     return;
276   }
277
278   Int_t nJetsTrig  = 0;
279   Int_t nJetsAssoc = 0;
280   if(type==0) {
281     nJetsTrig  = GetNJets(fContainerFullMC);
282     nJetsAssoc = nJetsTrig;
283   }
284   else if(type==1) {
285     nJetsTrig  = GetNJets(fContainerChargedMC);
286     nJetsAssoc = nJetsTrig;
287   }
288   else if(type==2) {
289     nJetsTrig  = GetNJets(fContainerFullMC);
290     nJetsAssoc = GetNJets(fContainerChargedMC);
291   }
292
293
294   for(Int_t ijt=0; ijt<nJetsTrig; ijt++) {
295
296     AliEmcalJet *jetTrigMC = static_cast<AliEmcalJet*>(GetAcceptJetFromArray(ijt, typetMC));
297     if(!jetTrigMC) continue; //jet not selected
298
299     Double_t jetTrigPtMC = GetJetPt(jetTrigMC,typetMC);
300
301     if(jetTrigPtMC<fPtMinTriggerJet)
302       continue;
303
304     if(type==1)
305       fh1TriggersCharged[0]->Fill(jetTrigPtMC);
306     if(type==2)
307       fh1TriggersFull[0]->Fill(jetTrigPtMC);
308
309     AliEmcalJet *jetTrigDet = jetTrigMC->ClosestJet();
310     if(!jetTrigDet) {
311       //trigger is lost
312       if(type==1)
313         fh1TriggersLostCharged->Fill(jetTrigPtMC);
314       if(type==2)
315         fh1TriggersLostFull->Fill(jetTrigPtMC);
316       
317       continue;
318     }
319
320     if(type==1)
321       fh1TriggersCharged[1]->Fill(GetJetPt(jetTrigDet,typet));
322     if(type==2)
323       fh1TriggersFull[1]->Fill(GetJetPt(jetTrigDet,typet));
324
325     for(Int_t ija=0; ija<nJetsAssoc; ija++) {
326       if(IsSameJet(ijt,ija,type,kTRUE)) continue;
327
328       AliEmcalJet *jetAssocMC = static_cast<AliEmcalJet*>(GetAcceptJetFromArray(ija, typeaMC));
329       if(!jetAssocMC) continue;
330
331       Double_t jetAssocPtMC = GetJetPt(jetAssocMC,typeaMC);
332
333       //Now check if jets are also there on detector level
334       AliEmcalJet *jetAssocDet = jetAssocMC->ClosestJet();
335       if(!jetAssocDet) {
336         //dijet is lost
337       if(type==1)
338         fh3AssocLostPtDeltaPhiCharged->Fill(jetTrigPtMC,jetAssocPtMC,GetDeltaPhi(jetTrigMC,jetAssocMC));
339       if(type==2)
340         fh3AssocLostPtDeltaPhiFull->Fill(jetTrigPtMC,jetAssocPtMC,GetDeltaPhi(jetTrigMC,jetAssocMC));
341         continue;
342       }
343
344       FillDiJetResponse(jetTrigMC,jetAssocMC,jetTrigDet,jetAssocDet,type);
345
346     } // associate jet loop
347   }//trigger jet loop
348
349 }
350
351 //________________________________________________________________________
352 void AliAnalysisTaskEmcalDiJetResponse::CorrelateTwoJets(const Int_t type) {
353   //
354   // Correlate jets and fill histos
355   //
356
357   Int_t typet = 0;
358   Int_t typea = 0;
359   Int_t typetMC = 0;
360   Int_t typeaMC = 0;
361   if(type==0) { //full-full
362     typetMC = fContainerFullMC;
363     typeaMC = fContainerFullMC;
364     typet = fContainerFull;
365     typea = fContainerFull;
366   }
367   else if(type==1) { //charged-charged
368     typetMC = fContainerChargedMC;
369     typeaMC = fContainerChargedMC;
370     typet = fContainerCharged;
371     typea = fContainerCharged;
372   }
373   else if(type==2) { //full-charged
374     typetMC = fContainerFullMC;
375     typeaMC = fContainerChargedMC;
376     typet = fContainerFull;
377     typea = fContainerCharged;
378   }
379   else {
380     AliWarning(Form("%s: type %d of dijet correlation not defined!",GetName(),type));
381     return;
382   }
383
384   Int_t nJetsTrig  = 0;
385   if(type==0) {
386     nJetsTrig  = GetNJets(fContainerFullMC);
387   }
388   else if(type==1) {
389     nJetsTrig  = GetNJets(fContainerChargedMC);
390   }
391   else if(type==2) {
392     nJetsTrig  = GetNJets(fContainerFullMC);
393   }
394
395   for(Int_t ijt=0; ijt<nJetsTrig; ijt++) {
396
397     AliEmcalJet *jetTrigMC = static_cast<AliEmcalJet*>(GetAcceptJetFromArray(ijt, typetMC));
398     if(!jetTrigMC) continue; //jet not selected
399
400     Double_t jetTrigPtMC = GetJetPt(jetTrigMC,typetMC);
401
402     if(jetTrigPtMC<fPtMinTriggerJet)
403       continue;
404
405     if(type==1)
406       fh1TriggersCharged[0]->Fill(jetTrigPtMC);
407     if(type==2)
408       fh1TriggersFull[0]->Fill(jetTrigPtMC);
409
410     AliEmcalJet *jetTrigDet = jetTrigMC->ClosestJet();
411     if(!jetTrigDet) {
412       //trigger is lost
413       if(type==1)
414         fh1TriggersLostCharged->Fill(jetTrigPtMC);
415       if(type==2)
416         fh1TriggersLostFull->Fill(jetTrigPtMC);
417       continue;
418     }
419
420     if(type==1)
421       fh1TriggersCharged[1]->Fill(GetJetPt(jetTrigDet,typet));
422     if(type==2)
423       fh1TriggersFull[1]->Fill(GetJetPt(jetTrigDet,typet));
424
425
426     AliEmcalJet *jetAssocMC = GetLeadingJetOppositeHemisphere(type,typeaMC,jetTrigMC);
427     if(!jetAssocMC) continue;
428
429     Double_t jetAssocPtMC = GetJetPt(jetAssocMC,typeaMC);
430       
431     //Now check if jets are also there on detector level
432     AliEmcalJet *jetAssocDet = jetAssocMC->ClosestJet();
433     if(!jetAssocDet) {
434       //dijet is lost
435       if(type==1)
436         fh3AssocLostPtDeltaPhiCharged->Fill(jetTrigPtMC,jetAssocPtMC,GetDeltaPhi(jetTrigMC,jetAssocMC));
437       if(type==2)
438         fh3AssocLostPtDeltaPhiFull->Fill(jetTrigPtMC,jetAssocPtMC,GetDeltaPhi(jetTrigMC,jetAssocMC));
439       continue;
440     }
441
442     if(fDoPtBias) {
443       if(type==0 || type==1) {
444         if(GetJetPt(jetAssocDet,typea)>GetJetPt(jetTrigDet,typet))
445           continue;
446       }
447     }
448
449     FillDiJetResponse(jetTrigMC,jetAssocMC,jetTrigDet,jetAssocDet,type);
450
451
452   }//trigger jet loop
453
454 }
455
456 //________________________________________________________________________
457 void AliAnalysisTaskEmcalDiJetResponse::FillDiJetResponse(const AliEmcalJet *jetTrigMC, const AliEmcalJet *jetAssocMC, const AliEmcalJet *jetTrigDet, const AliEmcalJet *jetAssocDet, Int_t type) {
458
459   //Fill dijet response
460
461   Int_t typet = 0;
462   Int_t typea = 0;
463   Int_t typetMC = 0;
464   Int_t typeaMC = 0;
465   if(type==0) { //full-full
466     typetMC = fContainerFullMC;
467     typeaMC = fContainerFullMC;
468     typet = fContainerFull;
469     typea = fContainerFull;
470   }
471   else if(type==1) { //charged-charged
472     typetMC = fContainerChargedMC;
473     typeaMC = fContainerChargedMC;
474     typet = fContainerCharged;
475     typea = fContainerCharged;
476   }
477   else if(type==2) { //full-charged
478     typetMC = fContainerFullMC;
479     typeaMC = fContainerChargedMC;
480     typet = fContainerFull;
481     typea = fContainerCharged;
482   }
483   else {
484     AliWarning(Form("%s: type %d of dijet correlation not defined!",GetName(),type));
485     return;
486   }
487
488   Double_t jetTrigPtMC  = GetJetPt(jetTrigMC,typetMC);
489   Double_t jetAssocPtMC = GetJetPt(jetAssocMC,typeaMC);
490
491   Double_t varDet[2] = {TMath::Abs(GetJetPt(jetTrigDet,typet)*TMath::Sin(GetDeltaPhi(jetTrigDet,jetAssocDet))),(jetTrigDet->Eta()+jetAssocDet->Eta())/2.};
492   Double_t varPart[2] = {TMath::Abs(jetTrigPtMC*TMath::Sin(GetDeltaPhi(jetTrigMC,jetAssocMC))),(jetTrigMC->Eta()+jetAssocMC->Eta())/2.};
493
494   Double_t ajDet  = (GetJetPt(jetTrigDet,typet)-GetJetPt(jetAssocDet,typea))/(GetJetPt(jetTrigDet,typet)+GetJetPt(jetAssocDet,typea));
495   Double_t ajPart = (jetTrigPtMC-jetAssocPtMC)/(jetTrigPtMC+jetAssocPtMC);
496
497   //Store dijet vars: pt,trig MC; pt,trig DET; pt,ass MC; pt,ass DET; dPhi MC; dPhi Det; kT MC; kT Det;
498   Double_t diJetVars[10] = {
499     jetTrigPtMC,
500     GetJetPt(jetTrigDet,typet),
501     jetAssocPtMC,
502     GetJetPt(jetAssocDet,typea),
503     GetDeltaPhi(jetTrigMC,jetAssocMC),
504     GetDeltaPhi(jetTrigDet,jetAssocDet),
505     varPart[fnUsedResponseVar],
506     varDet[fnUsedResponseVar],
507     ajDet,
508     ajPart
509   }; 
510   
511   if(type==1)
512     fhnDiJetResponseCharged->Fill(diJetVars);
513   else if(type==2)
514     fhnDiJetResponseFullCharged->Fill(diJetVars);
515
516
517 }
518
519 //________________________________________________________________________
520 void AliAnalysisTaskEmcalDiJetResponse::FillMatchHistos() {
521   //
522   // Fill Particle-Detector level matching histos
523   //
524
525   for(int i = 0; i < GetNJets(fContainerFull);++i) {
526     AliEmcalJet *jetDet = static_cast<AliEmcalJet*>(GetAcceptJetFromArray(i, fContainerFull));
527     if(!jetDet) continue;
528
529     AliEmcalJet *jetPart = jetDet->ClosestJet();
530     if(!jetPart) continue;
531
532     Double_t matchVars[6] = {
533       jetPart->Pt(),
534       jetDet->Pt(),
535       GetDeltaPhi(jetPart->Phi(),jetDet->Phi()),
536       jetPart->Eta()-jetDet->Eta(),
537       GetDeltaR(jetPart,jetDet),
538       TMath::Min((Float_t)i+0.5,2.5)
539     };
540     fhnMatchingFull->Fill(matchVars);
541
542   }//loop over full jets
543
544   for(int i = 0; i < GetNJets(fContainerCharged);++i) {
545     AliEmcalJet *jetDet = static_cast<AliEmcalJet*>(GetAcceptJetFromArray(i, fContainerCharged));
546     if(!jetDet) continue;
547
548     AliEmcalJet *jetPart = jetDet->ClosestJet();
549     if(!jetPart) continue;
550
551     Double_t matchVars[6] = {
552       jetPart->Pt(),
553       jetDet->Pt(),
554       GetDeltaPhi(jetPart->Phi(),jetDet->Phi()),
555       jetPart->Eta()-jetDet->Eta(),
556       GetDeltaR(jetPart,jetDet),
557       TMath::Min((Float_t)i+0.5,2.5)
558     };
559     fhnMatchingCharged->Fill(matchVars);
560
561   }//loop over charged jets
562
563 }
564
565
566 //________________________________________________________________________
567 Bool_t AliAnalysisTaskEmcalDiJetResponse::RetrieveEventObjects() {
568   //
569   // retrieve event objects
570   //
571
572   if (!AliAnalysisTaskEmcalDiJetBase::RetrieveEventObjects())
573     return kFALSE;
574
575   return kTRUE;
576
577 }
578
579 //_______________________________________________________________________
580 void AliAnalysisTaskEmcalDiJetResponse::Terminate(Option_t *) 
581 {
582   // Called once at the end of the analysis.
583 }