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 }