JETAN/AliAnalysisTaskJetBackgroundSubtract.cxx

   1 // **************************************
   2 // Task used for the correction of determiantion of reconstructed jet spectra
   3 // Compares input (gen) and output (rec) jets
   4 // *******************************************
   5
   6
   7 /**************************************************************************
   8  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   9  *                                                                        *
  10  * Author: The ALICE Off-line Project.                                    *
  11  * Contributors are mentioned in the code where appropriate.              *
  12  *                                                                        *
  13  * Permission to use, copy, modify and distribute this software and its   *
  14  * documentation strictly for non-commercial purposes is hereby granted   *
  15  * without fee, provided that the above copyright notice appears in all   *
  16  * copies and that both the copyright notice and this permission notice   *
  17  * appear in the supporting documentation. The authors make no claims     *
  18  * about the suitability of this software for any purpose. It is          *
  19  * provided "as is" without express or implied warranty.                  *
  20  **************************************************************************/
  21
  22
  23 #include <TROOT.h>
  24 #include <TH1F.h>
  25 #include <TH2F.h>
  26 #include <THnSparse.h>
  27 #include <TSystem.h>
  28 #include <TInterpreter.h>
  29 #include <TList.h>
  30 #include <TLorentzVector.h>
  31 #include <TRefArray.h>
  32 #include  "TDatabasePDG.h"
  33
  34 #include "AliAnalysisTaskJetBackgroundSubtract.h"
  35 #include "AliAnalysisManager.h"
  36 #include "AliAODHandler.h"
  37 #include "AliAODTrack.h"
  38 #include "AliAODJet.h"
  39 #include "AliAODEvent.h"
  40 #include "AliInputEventHandler.h"
  41 #include "AliAODJetEventBackground.h"
  42
  43
  44 ClassImp(AliAnalysisTaskJetBackgroundSubtract)
  45
  46 AliAnalysisTaskJetBackgroundSubtract::~AliAnalysisTaskJetBackgroundSubtract(){
  47   delete fJBArray;
  48   delete fOutJetArrayList;
  49   delete fInJetArrayList;
  50 }
  51
  52 AliAnalysisTaskJetBackgroundSubtract::AliAnalysisTaskJetBackgroundSubtract():
  53   AliAnalysisTaskSE(),
  54   fAODOut(0x0),
  55   fAODIn(0x0),
  56   fAODExtension(0x0),
  57   fJBArray(new TObjArray()),
  58   fBackgroundBranch(""),
  59   fNonStdFile(""),
  60   fReplaceString1("B0"),
  61   fReplaceString2("B%d"),
  62   fSubtraction(kRhoRecalc),
  63   fInJetArrayList(0x0),
  64   fOutJetArrayList(0x0),
  65   fHistList(0x0)
  66 {
  67
  68 }
  69
  70 AliAnalysisTaskJetBackgroundSubtract::AliAnalysisTaskJetBackgroundSubtract(const char* name):
  71
  72   AliAnalysisTaskSE(name),
  73   fAODOut(0x0),
  74   fAODIn(0x0),
  75   fAODExtension(0x0),
  76   fJBArray(new TObjArray()),
  77   fBackgroundBranch(""),
  78   fNonStdFile(""),
  79   fReplaceString1("B0"),
  80   fReplaceString2("B%d"),
  81   fSubtraction(kRhoRecalc),
  82   fInJetArrayList(0x0),
  83   fOutJetArrayList(0x0),
  84   fHistList(0x0)
  85 {
  86  DefineOutput(1, TList::Class());
  87 }
  88
  89
  90
  91 Bool_t AliAnalysisTaskJetBackgroundSubtract::Notify()
  92 {
  93   //
  94   fAODIn = dynamic_cast<AliAODEvent*>(InputEvent());
  95
  96   ResetOutJets();
  97
  98   // Now we also have the Input Event Available! Fillvthe pointers in the list
  99
 100   fInJetArrayList->Clear();
 101   fOutJetArrayList->Clear();
 102
 103   for(int iJB = 0;iJB<fJBArray->GetEntries();iJB++){
 104     TObjString *ostr = (TObjString*)fJBArray->At(iJB);
 105
 106
 107     TClonesArray* jarray = 0;
 108     if(!jarray&&fAODOut){
 109       jarray = (TClonesArray*)(fAODOut->FindListObject(ostr->GetString().Data()));
 110     }
 111     if(!jarray&&fAODExtension){
 112       jarray = (TClonesArray*)(fAODExtension->GetAOD()->FindListObject(ostr->GetString().Data()));
 113     }
 114     if(!jarray&&fAODIn){
 115       jarray = (TClonesArray*)(fAODIn->FindListObject(ostr->GetString().Data()));
 116     }
 117
 118     if(!jarray){
 119       if(fDebug){
 120         Printf("%s:%d Input jet branch %s not found",(char*)__FILE__,__LINE__,ostr->GetString().Data());
 121       }
 122       continue;
 123     }
 124
 125     TString newName(ostr->GetString().Data());
 126     newName.ReplaceAll(fReplaceString1.Data(),Form(fReplaceString2.Data(),fSubtraction));
 127     TClonesArray* jarrayOut = 0;
 128     if(!jarrayOut&&fAODOut){
 129       jarrayOut = (TClonesArray*)(fAODOut->FindListObject(newName.Data()));
 130     }
 131     if(!jarrayOut&&fAODExtension){
 132       jarrayOut = (TClonesArray*)(fAODExtension->GetAOD()->FindListObject(newName.Data()));
 133     }
 134
 135     if(!jarrayOut){
 136       if(fDebug){
 137         Printf("%s:%d Output jet branch %s not found",(char*)__FILE__,__LINE__,newName.Data());
 138         PrintAODContents();
 139       }
 140       continue;
 141     }
 142     if(jarrayOut&&jarray){
 143       fOutJetArrayList->Add(jarrayOut);
 144       fInJetArrayList->Add(jarray);
 145     }
 146   }
 147   return kTRUE;
 148 }
 149
 150 void AliAnalysisTaskJetBackgroundSubtract::UserCreateOutputObjects()
 151 {
 152
 153   //
 154   // Create the output container
 155   //
 156   // Connect the AOD
 157
 158   if (fDebug > 1) printf("AnalysisTaskJetBackgroundSubtract::UserCreateOutputObjects() \n");
 159   if(fNonStdFile.Length()!=0){
 160
 161     // case that we have an AOD extension we need to fetch the jets from the extended output
 162     // we identifay the extension aod event by looking for the branchname
 163     AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
 164     fAODExtension = aodH->GetExtension(fNonStdFile.Data());
 165
 166     if(!fAODExtension){
 167       if(fDebug>1)Printf("AODExtension found for %s",fNonStdFile.Data());
 168     }
 169   }
 170   fAODOut = AODEvent();
 171
 172   // usually we do not have the input already here
 173
 174   if(!fInJetArrayList)fInJetArrayList =new TList();
 175   if(!fOutJetArrayList)fOutJetArrayList =new TList();
 176
 177   for(int iJB = 0;iJB<fJBArray->GetEntries();iJB++){
 178     TObjString *ostr = (TObjString*)fJBArray->At(iJB);
 179     TString newName(ostr->GetString().Data());
 180     if(!newName.Contains(fReplaceString1.Data())){
 181       Printf("%s:%d cannot replace string %s in %s",(char*)__FILE__,__LINE__,fReplaceString1.Data(),
 182              newName.Data());
 183       continue;
 184     }
 185
 186     // add a new branch to the output for the background subtracted jets take the names from
 187     // the input jets and replace the background flag names
 188     TClonesArray *tca = new TClonesArray("AliAODJet", 0);
 189     newName.ReplaceAll(fReplaceString1.Data(),Form(fReplaceString2.Data(),fSubtraction));
 190     if(fDebug){
 191       Printf("%s:%d created branch \n %s from \n %s",(char*)__FILE__,__LINE__,newName.Data(),
 192              ostr->GetString().Data());
 193     }
 194     tca->SetName(newName.Data());
 195     AddAODBranch("TClonesArray",&tca,fNonStdFile.Data());
 196   }
 197
 198
 199   if(!fHistList)fHistList = new TList();
 200   fHistList->SetOwner();
 201
 202   for(int iJB = 0;iJB<fJBArray->GetEntries();iJB++){
 203     TObjString *ostr = (TObjString*)fJBArray->At(iJB);
 204     TString oldName(ostr->GetString().Data());
 205     TString newName(ostr->GetString().Data());
 206     if(!newName.Contains(fReplaceString1.Data())){
 207       Printf("%s:%d cannot replace string %s in %s",(char*)__FILE__,__LINE__,fReplaceString1.Data(),
 208              newName.Data());
 209       continue;
 210     }
 211     newName.ReplaceAll(fReplaceString1.Data(),Form(fReplaceString2.Data(),fSubtraction));
 212     TH2F *hTmp = new TH2F(Form("h2PtInPtOut_%d",iJB),Form(";%s p_{T}; %s p_{T}",oldName.Data(),newName.Data()),200,0,200.,200,0.,200.);
 213     fHistList->Add(hTmp);
 214   }
 215
 216   Bool_t oldStatus = TH1::AddDirectoryStatus();
 217   TH1::AddDirectory(kFALSE);
 218
 219   //
 220   //  Histogram booking, add som control histograms here
 221   //
 222
 223   // =========== Switch on Sumw2 for all histos ===========
 224   for (Int_t i=0; i<fHistList->GetEntries(); ++i) {
 225     TH1 *h1 = dynamic_cast<TH1*>(fHistList->At(i));
 226     if (h1){
 227       h1->Sumw2();
 228       continue;
 229     }
 230     THnSparse *hn = dynamic_cast<THnSparse*>(fHistList->At(i));
 231     if(hn)hn->Sumw2();
 232   }
 233   TH1::AddDirectory(oldStatus);
 234
 235   if(fBackgroundBranch.Length()==0)
 236     AliError(Form("%s:%d No BackgroundBranch defined",(char*)__FILE__,__LINE__));
 237   if(fJBArray->GetEntries()==0)
 238     AliError(Form("%s:%d No Jet Branches defined defined",(char*)__FILE__,__LINE__));
 239 }
 240
 241 void AliAnalysisTaskJetBackgroundSubtract::Init()
 242 {
 243   //
 244   // Initialization
 245   //
 246   if (fDebug > 1) printf("AnalysisTaskJetBackgroundSubtract::Init() \n");
 247 }
 248
 249 void AliAnalysisTaskJetBackgroundSubtract::UserExec(Option_t */*option*/)
 250 {
 251
 252   if (fDebug > 1) printf("AnalysisTaskJetBackgroundSubtract::UserExec() \n");
 253   ResetOutJets();
 254   if(fBackgroundBranch.Length()==0||fJBArray->GetEntries()==0){
 255     if(fDebug)Printf("%s:%d No background subtraction done",(char*)__FILE__,__LINE__);
 256     PostData(1,fHistList);
 257   }
 258   if(fJBArray->GetEntries()!=fInJetArrayList->GetEntries()){
 259     if(fDebug)Printf("%s:%d different Array  sizes %d %d %d",(char*)__FILE__,__LINE__,fJBArray->GetEntries(),fInJetArrayList->GetEntries(),fOutJetArrayList->GetEntries());
 260     PostData(1,fHistList);
 261   }
 262
 263
 264
 265   AliAODJetEventBackground*  evBkg = 0;
 266   TClonesArray*              bkgClusters = 0;
 267   TClonesArray*              bkgClustersRC = 0;
 268   TString bkgClusterName(fBackgroundBranch.Data());
 269   bkgClusterName.ReplaceAll(Form("%s_",AliAODJetEventBackground::StdBranchName()),"");
 270   TString bkgClusterRCName(Form("%s%s",bkgClusterName.Data(),"RandomCone"));
 271
 272   if(!evBkg&&!bkgClusters&&!bkgClustersRC&&fAODOut){
 273     evBkg = (AliAODJetEventBackground*)(fAODOut->FindListObject(fBackgroundBranch.Data()));
 274     bkgClusters = (TClonesArray*)(fAODOut->FindListObject(bkgClusterName.Data()));
 275     bkgClustersRC = (TClonesArray*)(fAODOut->FindListObject(bkgClusterRCName.Data()));
 276
 277     if(fDebug&&bkgClusters)Printf("%s:%d Background cluster branch %s found",(char*)__FILE__,__LINE__,bkgClusterName.Data());
 278     if(fDebug&&bkgClustersRC)Printf("%s:%d Background cluster RC branch %s found",(char*)__FILE__,__LINE__,bkgClusterRCName.Data());
 279     if(fDebug&&evBkg)Printf("%s:%d Backgroundbranch %s found",(char*)__FILE__,__LINE__,fBackgroundBranch.Data());
 280   }
 281   if(!evBkg&&!bkgClusters&&!bkgClustersRC&&fAODExtension){
 282     evBkg = (AliAODJetEventBackground*)(fAODExtension->GetAOD()->FindListObject(fBackgroundBranch.Data()));
 283     bkgClusters = (TClonesArray*)(fAODExtension->GetAOD()->FindListObject(bkgClusterName.Data()));
 284     bkgClustersRC = (TClonesArray*)(fAODExtension->GetAOD()->FindListObject(bkgClusterRCName.Data()));
 285     if(fDebug&&bkgClusters)Printf("%s:%d Background cluster branch %s found",(char*)__FILE__,__LINE__,bkgClusterName.Data());
 286     if(fDebug&&bkgClustersRC)Printf("%s:%d Background cluster RC branch %s found",(char*)__FILE__,__LINE__,bkgClusterRCName.Data());
 287
 288     if(fDebug&&evBkg)Printf("%s:%d Backgroundbranch %s found",(char*)__FILE__,__LINE__,fBackgroundBranch.Data());
 289   }
 290
 291   if(!evBkg&&!bkgClusters&&!bkgClustersRC&&fAODIn){
 292     evBkg = (AliAODJetEventBackground*)(fAODIn->FindListObject(fBackgroundBranch.Data()));
 293     bkgClusters = (TClonesArray*)(fAODIn->FindListObject(bkgClusterName.Data()));
 294     bkgClustersRC = (TClonesArray*)(fAODIn->FindListObject(bkgClusterRCName.Data()));
 295
 296     if(fDebug&&bkgClusters)Printf("%s:%d Background cluster branch %s found",(char*)__FILE__,__LINE__,bkgClusterName.Data());
 297     if(fDebug&&bkgClustersRC)Printf("%s:%d Background cluster RC branch %s found",(char*)__FILE__,__LINE__,bkgClusterRCName.Data());
 298     if(fDebug&&evBkg)Printf("%s:%d Backgroundbranch %s found",(char*)__FILE__,__LINE__,fBackgroundBranch.Data());
 299   }
 300
 301   if(!evBkg&&(fSubtraction==kArea||fSubtraction==kRhoRecalc)){
 302     if(fDebug){
 303       Printf("%s:%d Backgroundbranch %s not found",(char*)__FILE__,__LINE__,fBackgroundBranch.Data());
 304       PrintAODContents();
 305     }
 306     PostData(1,fHistList);
 307     return;
 308   }
 309
 310   if(!bkgClusters&&(fSubtraction==kRhoRecalc)){
 311     if(fDebug){
 312       Printf("%s:%d Background cluster branch %s not found",(char*)__FILE__,__LINE__,bkgClusterName.Data());
 313       PrintAODContents();
 314     }
 315     PostData(1,fHistList);
 316     return;
 317   }
 318
 319   if(!bkgClustersRC&&(fSubtraction==kRhoRC)){
 320     if(fDebug){
 321       Printf("%s:%d Background cluster RC branch %s not found",(char*)__FILE__,__LINE__,bkgClusterRCName.Data());
 322       PrintAODContents();
 323     }
 324     PostData(1,fHistList);
 325     return;
 326   }
 327   // LOOP over all jet branches and subtract the background
 328
 329    Float_t rho = 0;
 330    Double_t meanarea = 0;
 331    if(fSubtraction==kArea)rho =  evBkg->GetBackground(2);
 332    if(fSubtraction==kRhoRecalc){
 333      meanarea=evBkg->GetMeanarea(2);
 334      rho =RecalcRho(bkgClusters,meanarea);
 335    }
 336    if(fSubtraction==kRhoRC)rho =RhoRC(bkgClustersRC);
 337
 338   for(int iJB = 0;iJB<fInJetArrayList->GetEntries();iJB++){
 339     TClonesArray* jarray = (TClonesArray*)fInJetArrayList->At(iJB);
 340     TClonesArray* jarrayOut = (TClonesArray*)fOutJetArrayList->At(iJB);
 341
 342     if(!jarray||!jarrayOut){
 343     Printf("%s:%d Array not found %d: %p %p",(char*)__FILE__,__LINE__,iJB,jarray->GetName(),jarrayOut->GetName());
 344     continue;
 345       }
 346     TH2F* h2PtInOut = (TH2F*)fHistList->FindObject(Form("h2PtInPtOut_%d",iJB));
 347     // loop over all jets
 348     Int_t nOut = 0;
 349
 350
 351     for(int i = 0;i < jarray->GetEntriesFast();i++){
 352       AliAODJet *jet = (AliAODJet*)jarray->At(i);
 353       AliAODJet tmpNewJet(*jet);
 354       Bool_t bAdd = false;
 355
 356
 357       if(fSubtraction==kArea){
 358         Double_t background = rho * jet->EffectiveAreaCharged();
 359         Float_t ptSub = jet->Pt() - background;
 360         if(fDebug>2){
 361           Printf("%s:%d Jet %d %3.3f %3.3f",(char*)__FILE__,__LINE__,i,jet->Pt(),ptSub);
 362         }
 363         if(ptSub<0){
 364           // optionally rescale it and keep??
 365           bAdd = RescaleJetMomentum(&tmpNewJet,0.1);
 366           if(h2PtInOut)h2PtInOut->Fill(jet->Pt(),0.1);
 367         }
 368         else{
 369           bAdd = RescaleJetMomentum(&tmpNewJet,ptSub);
 370           if(h2PtInOut)h2PtInOut->Fill(jet->Pt(),ptSub);
 371         }
 372         // add background estimates to the new jet object
 373         // allows to recover old p_T and rho...
 374         tmpNewJet.SetBgEnergy(background,0);
 375       }// kAREA
 376       else if(fSubtraction==kRhoRecalc){
 377         Double_t background = rho * jet->EffectiveAreaCharged();
 378         Float_t ptSub = jet->Pt() - background;
 379         if(fDebug>2){
 380           Printf("%s:%d Jet %d %3.3f %3.3f %3.3f %3.3f",(char*)__FILE__,__LINE__,i,jet->Pt(),ptSub,background,rho);}
 381         if(ptSub<0){
 382           // optionally rescale it and keep??
 383           bAdd = RescaleJetMomentum(&tmpNewJet,0.1);
 384           if(h2PtInOut)h2PtInOut->Fill(jet->Pt(),0.1);
 385         }
 386         else{
 387           bAdd = RescaleJetMomentum(&tmpNewJet,ptSub);
 388           if(h2PtInOut)h2PtInOut->Fill(jet->Pt(),ptSub);
 389         }
 390         // add background estimates to the new jet object
 391         // allows to recover old p_T and rho...
 392         tmpNewJet.SetBgEnergy(background,0);
 393
 394       }//kRhoRecalc
 395        else if(fSubtraction==kRhoRC){
 396         Double_t background = rho * jet->EffectiveAreaCharged();
 397         Float_t ptSub = jet->Pt() - background;
 398         if(fDebug>2){   Printf("%s:%d Jet %d %3.3f %3.3f %3.3f %3.3f",(char*)__FILE__,__LINE__,i,jet->Pt(),ptSub,background,rho);}
 399         if(ptSub<0){
 400           // optionally rescale it and keep??
 401           bAdd = RescaleJetMomentum(&tmpNewJet,0.1);
 402           if(h2PtInOut)h2PtInOut->Fill(jet->Pt(),0.1);
 403         }
 404         else{
 405           bAdd = RescaleJetMomentum(&tmpNewJet,ptSub);
 406           if(h2PtInOut)h2PtInOut->Fill(jet->Pt(),ptSub);
 407         }
 408         // add background estimates to the new jet object
 409         // allows to recover old p_T and rho...
 410         tmpNewJet.SetBgEnergy(background,0);
 411
 412        }//kRhoRC
 413
 414       if(bAdd){
 415         AliAODJet *newJet = new ((*jarrayOut)[nOut++]) AliAODJet(tmpNewJet);
 416         // what about track references, clear for now...
 417         newJet->GetRefTracks()->Clear();
 418       }
 419
 420
 421
 422
 423
 424
 425     }
 426
 427
 428
 429     // subtract the background
 430
 431
 432     // remove jets??
 433
 434     // sort jets...
 435
 436   }
 437   PostData(1, fHistList);
 438 }
 439
 440 void AliAnalysisTaskJetBackgroundSubtract::Terminate(Option_t */*option*/)
 441 {
 442   // Terminate analysis
 443   //
 444   if (fDebug > 1) printf("AnalysisJetBackgroundSubtract: Terminate() \n");
 445 }
 446
 447 Bool_t AliAnalysisTaskJetBackgroundSubtract::RescaleJetMomentum(AliAODJet *jet,Float_t pT){
 448   // keep the direction and the jet mass
 449   if(pT<=0)return kFALSE;
 450   Double_t pTold = jet->Pt();
 451   Double_t scale  = pT/pTold;
 452   Double_t mass  = jet->M();
 453   Double_t pNew = jet->P() * scale;
 454   jet->SetPxPyPzE(scale*jet->Px(),scale*jet->Py(),scale*jet->Pz(),TMath::Sqrt(mass*mass+pNew*pNew));
 455   return kTRUE;
 456 }
 457
 458 Double_t AliAnalysisTaskJetBackgroundSubtract::RecalcRho(TClonesArray* bkgClusters,Double_t meanarea){
 459
 460        Double_t ptarea=0.;
 461        Int_t count=0;
 462        Double_t rho=0.;
 463        const Double_t Rlimit2=0.8*0.8;  //2*jet radius.
 464        TClonesArray* jarray=0;
 465
 466         for(int iJB = 0;iJB<fInJetArrayList->GetEntries();iJB++){
 467         TObjString *ostr = (TObjString*)fInJetArrayList->At(iJB);
 468         TString jetref=ostr->GetString().Data();
 469         if(jetref.Contains("ANTIKT04")){
 470           jarray = (TClonesArray*)fInJetArrayList->At(iJB);}}
 471         if(jarray->GetEntries()>=2){
 472           AliAODJet *first = (AliAODJet*)(jarray->At(0));
 473           AliAODJet *second= (AliAODJet*)(jarray->At(1));
 474           for(Int_t k=0;k<bkgClusters->GetEntriesFast();k++){
 475             AliAODJet *clus = (AliAODJet*)(bkgClusters->At(k));
 476             if(TMath::Abs(clus->Eta())>0.5) continue;
 477             if((clus->EffectiveAreaCharged())<0.1*meanarea) continue;
 478             Double_t distance1=(first->Eta()-clus->Eta())*(first->Eta()-clus->Eta())+
 479               (first->Phi()-clus->Phi())*(first->Phi()-clus->Phi());
 480             Double_t distance2= (second->Eta()-clus->Eta())*(second->Eta()-clus->Eta())+
 481               (second->Phi()-clus->Phi())*(second->Phi()-clus->Phi());
 482             if((distance1<Rlimit2)||(distance2<Rlimit2)) continue;
 483             ptarea=ptarea+clus->Pt()/clus->EffectiveAreaCharged();
 484             count=count+1;}
 485           if(count!=0) rho=ptarea/count;
 486         }
 487         return rho;
 488 }
 489
 490    Double_t AliAnalysisTaskJetBackgroundSubtract::RhoRC(TClonesArray* bkgClustersRC){
 491
 492        Double_t ptarea=0.;
 493        Int_t count=0;
 494        Double_t rho=0.;
 495        const Double_t Rlimit2=0.8*0.8;  //2*jet radius.
 496        TClonesArray* jarray=0;
 497         for(int iJB = 0;iJB<fInJetArrayList->GetEntries();iJB++){
 498         TObjString *ostr = (TObjString*)fInJetArrayList->At(iJB);
 499         TString jetref=ostr->GetString().Data();
 500         if(jetref.Contains("ANTIKT04")){
 501          jarray = (TClonesArray*)fInJetArrayList->At(iJB);}}
 502          if(jarray->GetEntries()>=2){
 503            AliAODJet *first = (AliAODJet*)(jarray->At(0));
 504            AliAODJet *second=(AliAODJet*)(jarray->At(1));
 505          for(Int_t k=0;k<bkgClustersRC->GetEntriesFast();k++){
 506            AliAODJet *clus = (AliAODJet*)(bkgClustersRC->At(k));
 507            if(TMath::Abs(clus->Eta())>0.5) continue;
 508            Double_t distance1=(first->Eta()-clus->Eta())*(first->Eta()-clus->Eta())+
 509              (first->Phi()-clus->Phi())*(first->Phi()-clus->Phi());
 510            Double_t distance2= (second->Eta()-clus->Eta())*(second->Eta()-clus->Eta())+
 511              (second->Phi()-clus->Phi())*(second->Phi()-clus->Phi());
 512            if((distance1<Rlimit2)||(distance2<Rlimit2)) continue;
 513            ptarea=ptarea+clus->Pt()/clus->EffectiveAreaCharged();
 514            count=count+1;}
 515          if(count!=0) rho=ptarea/count;  }
 516          return rho;
 517 }
 518
 519
 520
 521
 522
 523
 524
 525
 526
 527 void AliAnalysisTaskJetBackgroundSubtract::ResetOutJets(){
 528   if(!fOutJetArrayList)return;
 529   for(int iJB = 0;iJB<fOutJetArrayList->GetEntries();iJB++){
 530     TClonesArray* jarray = (TClonesArray*)fOutJetArrayList->At(iJB);
 531     if(jarray)jarray->Delete();
 532   }
 533 }
 534
 535
 536 void AliAnalysisTaskJetBackgroundSubtract::PrintAODContents(){
 537   if(fAODIn){
 538     Printf("%s:%d >>>>>> Input",(char*)__FILE__,__LINE__);
 539     fAODIn->Print();
 540   }
 541   if(fAODExtension){
 542     Printf("%s:%d >>>>>> Extenstion",(char*)__FILE__,__LINE__);
 543     fAODExtension->GetAOD()->Print();
 544   }
 545   if(fAODOut){
 546     Printf("%s:%d >>>>>> Output",(char*)__FILE__,__LINE__);
 547     fAODOut->Print();
 548   }
 549 }