PWGJE/EMCALJetTasks/AliAnalysisTaskEmcalJet.cxx

   1 // $Id$
   2 //
   3 // Emcal jet analysis base task.
   4 //
   5 // Author: S.Aiola, M. Verweij
   6
   7 #include "AliAnalysisTaskEmcalJet.h"
   8
   9 #include <TClonesArray.h>
  10 #include <TList.h>
  11
  12 #include "AliAnalysisManager.h"
  13 #include "AliCentrality.h"
  14 #include "AliEMCALGeometry.h"
  15 #include "AliESDEvent.h"
  16 #include "AliEmcalJet.h"
  17 #include "AliLog.h"
  18 #include "AliRhoParameter.h"
  19 #include "AliLocalRhoParameter.h"
  20 #include "AliVCluster.h"
  21 #include "AliVEventHandler.h"
  22 #include "AliVParticle.h"
  23 #include "AliJetContainer.h"
  24
  25 ClassImp(AliAnalysisTaskEmcalJet)
  26
  27 //________________________________________________________________________
  28 AliAnalysisTaskEmcalJet::AliAnalysisTaskEmcalJet() :
  29   AliAnalysisTaskEmcal("AliAnalysisTaskEmcalJet"),
  30   fRhoName(),
  31   fLocalRhoName(),
  32   fJetCollArray(),
  33   fJets(0),
  34   fRho(0),
  35   fLocalRho(0),
  36   fRhoVal(0)
  37 {
  38   // Default constructor.
  39
  40   fJetCollArray.SetOwner(kTRUE);
  41 }
  42
  43 //________________________________________________________________________
  44 AliAnalysisTaskEmcalJet::AliAnalysisTaskEmcalJet(const char *name, Bool_t histo) :
  45   AliAnalysisTaskEmcal(name, histo),
  46   fRhoName(),
  47   fLocalRhoName(),
  48   fJetCollArray(),
  49   fJets(0),
  50   fRho(0),
  51   fLocalRho(0),
  52   fRhoVal(0)
  53 {
  54   // Standard constructor.
  55
  56   fJetCollArray.SetOwner(kTRUE);
  57 }
  58
  59 //________________________________________________________________________
  60 AliAnalysisTaskEmcalJet::~AliAnalysisTaskEmcalJet()
  61 {
  62   // Destructor
  63 }
  64
  65 //________________________________________________________________________
  66 Bool_t AliAnalysisTaskEmcalJet::AcceptBiasJet(AliEmcalJet *jet, Int_t c)
  67 {
  68   // Accept jet with a bias.
  69
  70   AliJetContainer *cont = GetJetContainer(c);
  71   if(!cont) {
  72     AliError(Form("%s:Container %d not found",GetName(),c));
  73     return 0;
  74   }
  75
  76   return cont->AcceptBiasJet(jet);
  77 }
  78
  79 //________________________________________________________________________
  80 Double_t AliAnalysisTaskEmcalJet::GetLeadingHadronPt(AliEmcalJet *jet, Int_t c)
  81 {
  82
  83   AliJetContainer *cont = GetJetContainer(c);
  84   if(!cont) {
  85     AliError(Form("%s:Container %d not found",GetName(),c));
  86     return 0;
  87   }
  88
  89   return cont->GetLeadingHadronPt(jet);
  90 }
  91
  92 //________________________________________________________________________
  93 Bool_t AliAnalysisTaskEmcalJet::AcceptJet(AliEmcalJet *jet, Int_t c)
  94 {
  95   // Return true if jet is accepted.
  96   if (!jet)
  97     return kFALSE;
  98
  99   AliJetContainer *cont = GetJetContainer(c);
 100   if(!cont) {
 101     AliError(Form("%s:Container %d not found",GetName(),c));
 102     return 0;
 103   }
 104
 105   return cont->AcceptJet(jet);
 106 }
 107
 108 //________________________________________________________________________
 109 AliRhoParameter *AliAnalysisTaskEmcalJet::GetRhoFromEvent(const char *name)
 110 {
 111   // Get rho from event.
 112
 113   AliRhoParameter *rho = 0;
 114   TString sname(name);
 115   if (!sname.IsNull()) {
 116     rho = dynamic_cast<AliRhoParameter*>(InputEvent()->FindListObject(sname));
 117     if (!rho) {
 118       AliWarning(Form("%s: Could not retrieve rho with name %s!", GetName(), name));
 119       return 0;
 120     }
 121   }
 122   return rho;
 123 }
 124
 125 //________________________________________________________________________
 126 AliLocalRhoParameter *AliAnalysisTaskEmcalJet::GetLocalRhoFromEvent(const char *name)
 127 {
 128   // Get local rho from event.
 129   AliLocalRhoParameter *rho = 0;
 130   TString sname(name);
 131   if (!sname.IsNull()) {
 132     rho = dynamic_cast<AliLocalRhoParameter*>(InputEvent()->FindListObject(sname));
 133     if (!rho) {
 134       AliWarning(Form("%s: Could not retrieve local rho with name %s!", GetName(), name));
 135       return 0;
 136     }
 137   }
 138   return rho;
 139 }
 140
 141
 142 //________________________________________________________________________
 143 void AliAnalysisTaskEmcalJet::ExecOnce()
 144 {
 145   // Init the analysis.
 146
 147   AliAnalysisTaskEmcal::ExecOnce();
 148
 149   if (!fRhoName.IsNull() && !fRho) { // get rho from the event
 150     fRho = dynamic_cast<AliRhoParameter*>(InputEvent()->FindListObject(fRhoName));
 151     if (!fRho) {
 152       AliError(Form("%s: Could not retrieve rho %s!", GetName(), fRhoName.Data()));
 153       fInitialized = kFALSE;
 154       return;
 155     }
 156   }
 157
 158   if (!fLocalRhoName.IsNull() && !fLocalRho) {
 159     fLocalRho = dynamic_cast<AliLocalRhoParameter*>(InputEvent()->FindListObject(fLocalRhoName));
 160     if (!fLocalRho) {
 161       AliError(Form("%s: Could not retrieve local rho %s!", GetName(), fLocalRhoName.Data()));
 162       fInitialized = kFALSE;
 163       return;
 164     }
 165   }
 166
 167   //Load all requested jet branches - each container knows name already
 168   if(fJetCollArray.GetEntriesFast()==0) {
 169     AliWarning("There are no jet collections");
 170     return;
 171   }
 172
 173   for(Int_t i =0; i<fJetCollArray.GetEntriesFast(); i++) {
 174     AliJetContainer *cont = static_cast<AliJetContainer*>(fJetCollArray.At(i));
 175     cont->SetRunNumber(InputEvent()->GetRunNumber());
 176     cont->SetEMCALGeometry();
 177     cont->SetArray(InputEvent());
 178     cont->LoadRho(InputEvent());
 179     cont->LoadPartonsInfo(InputEvent());
 180   }
 181
 182   //Get Jets, cuts and rho for first jet container
 183   AliJetContainer *cont = GetJetContainer(0);
 184
 185   if (!cont->GetArrayName().IsNull()) {
 186     fJets = cont->GetArray();
 187     if(!fJets && fJetCollArray.GetEntriesFast()>0) {
 188       AliError(Form("%s: Could not retrieve first jet branch!", GetName()));
 189       fInitialized = kFALSE;
 190       return;
 191     }
 192   }
 193
 194   if (!fRho) { // if rho name is not provided, tries to use the rho object of the first jet branch
 195     fRhoName = cont->GetRhoName();
 196     fRho = cont->GetRhoParameter();
 197   }
 198 }
 199
 200 //________________________________________________________________________
 201 Bool_t AliAnalysisTaskEmcalJet::IsJetCluster(AliEmcalJet* jet, Int_t iclus, Bool_t sorted) const
 202 {
 203   // Return true if cluster is in jet.
 204
 205   for (Int_t i = 0; i < jet->GetNumberOfClusters(); ++i) {
 206     Int_t ijetclus = jet->ClusterAt(i);
 207     if (sorted && ijetclus > iclus)
 208       return kFALSE;
 209     if (ijetclus == iclus)
 210       return kTRUE;
 211   }
 212   return kFALSE;
 213 }
 214
 215 //________________________________________________________________________
 216 Bool_t AliAnalysisTaskEmcalJet::IsJetTrack(AliEmcalJet* jet, Int_t itrack, Bool_t sorted) const
 217 {
 218   // Return true if track is in jet.
 219
 220   for (Int_t i = 0; i < jet->GetNumberOfTracks(); ++i) {
 221     Int_t ijettrack = jet->TrackAt(i);
 222     if (sorted && ijettrack > itrack)
 223       return kFALSE;
 224     if (ijettrack == itrack)
 225       return kTRUE;
 226   }
 227   return kFALSE;
 228 }
 229
 230 //________________________________________________________________________
 231 Bool_t AliAnalysisTaskEmcalJet::RetrieveEventObjects()
 232 {
 233   // Retrieve objects from event.
 234
 235   if (!AliAnalysisTaskEmcal::RetrieveEventObjects())
 236     return kFALSE;
 237
 238   if (fRho)
 239     fRhoVal = fRho->GetVal();
 240
 241   return kTRUE;
 242 }
 243
 244 //________________________________________________________________________
 245 AliJetContainer* AliAnalysisTaskEmcalJet::AddJetContainer(const char *n, TString defaultCutType, Float_t jetRadius) {
 246
 247   // Add particle container
 248   // will be called in AddTask macro
 249
 250   TString tmp = TString(n);
 251   if(tmp.IsNull()) return 0;
 252
 253   AliJetContainer *cont = 0x0;
 254   cont = new AliJetContainer();
 255   cont->SetArrayName(n);
 256   cont->SetJetRadius(jetRadius);
 257
 258   defaultCutType.ToUpper();
 259
 260   if(!defaultCutType.IsNull() && !defaultCutType.EqualTo("USER")) {
 261     if(defaultCutType.EqualTo("TPC"))
 262       cont->SetJetAcceptanceType(AliJetContainer::kTPC);
 263     else if(defaultCutType.EqualTo("EMCAL"))
 264       cont->SetJetAcceptanceType(AliJetContainer::kEMCAL);
 265     else
 266       AliWarning(Form("%s: default cut type %s not recognized. Not setting cuts.",GetName(),defaultCutType.Data()));
 267   } else
 268     cont->SetJetAcceptanceType(AliJetContainer::kUser);
 269
 270   fJetCollArray.Add(cont);
 271
 272   return cont;
 273 }
 274
 275 //________________________________________________________________________
 276 AliJetContainer* AliAnalysisTaskEmcalJet::GetJetContainer(Int_t i) const{
 277   // Get i^th jet container
 278
 279   if(i<0 || i>=fJetCollArray.GetEntriesFast()) return 0;
 280   AliJetContainer *cont = static_cast<AliJetContainer*>(fJetCollArray.At(i));
 281   return cont;
 282 }
 283
 284 //________________________________________________________________________
 285 AliJetContainer* AliAnalysisTaskEmcalJet::GetJetContainer(const char* name) const{
 286   // Get the jet container with name
 287
 288   AliJetContainer *cont = static_cast<AliJetContainer*>(fJetCollArray.FindObject(name));
 289   return cont;
 290 }
 291
 292 //________________________________________________________________________
 293 void AliAnalysisTaskEmcalJet::SetJetAcceptanceType(UInt_t t, Int_t c)
 294 {
 295   // Set acceptance cuts
 296   AliJetContainer *cont = GetJetContainer(c);
 297   if (cont) {
 298     cont->SetJetAcceptanceType((AliJetContainer::JetAcceptanceType)t);
 299   }
 300   else {
 301     AliError(Form("%s in SetJetAcceptanceType(...): container %d not found!",GetName(),c));
 302   }
 303 }
 304
 305 //________________________________________________________________________
 306 void AliAnalysisTaskEmcalJet::SetJetAcceptanceType(TString cutType, Int_t c) {
 307   //set acceptance cuts
 308
 309   AliJetContainer *cont = GetJetContainer(c);
 310   if (!cont) {
 311     AliError(Form("%s in SetJetAcceptanceType(...): container %d not found",GetName(),c));
 312     return;
 313   }
 314
 315   cutType.ToUpper();
 316
 317   if(!cutType.IsNull() && !cutType.EqualTo("USER")) {
 318     if(cutType.EqualTo("TPC"))
 319      cont->SetJetAcceptanceType(AliJetContainer::kTPC);
 320     else if(cutType.EqualTo("EMCAL"))
 321       cont->SetJetAcceptanceType(AliJetContainer::kEMCAL);
 322     else
 323       AliWarning(Form("%s: default cut type %s not recognized. Not setting cuts.",GetName(),cutType.Data()));
 324   } else
 325     cont->SetJetAcceptanceType(AliJetContainer::kUser);
 326 }
 327
 328 //________________________________________________________________________
 329 void AliAnalysisTaskEmcalJet::SetRhoName(const char *n, Int_t c)
 330 {
 331   if (c >= 0) {
 332     AliJetContainer *cont = GetJetContainer(c);
 333     if (cont) cont->SetRhoName(n);
 334     else AliError(Form("%s in SetRhoName(...): container %d not found",GetName(),c));
 335   }
 336   else {
 337     fRhoName = n;
 338   }
 339 }
 340
 341 //________________________________________________________________________
 342 void AliAnalysisTaskEmcalJet::SetJetEtaLimits(Float_t min, Float_t max, Int_t c)
 343 {
 344   AliJetContainer *cont = GetJetContainer(c);
 345   if (cont) cont->SetJetEtaLimits(min,max);
 346   else AliError(Form("%s in SetJetEtaLimits(...): container %d not found",GetName(),c));
 347 }
 348
 349 //________________________________________________________________________
 350 void AliAnalysisTaskEmcalJet::SetJetPhiLimits(Float_t min, Float_t max, Int_t c)
 351 {
 352   AliJetContainer *cont = GetJetContainer(c);
 353   if (cont) cont->SetJetPhiLimits(min,max);
 354   else AliError(Form("%s in SetJetPhiLimits(...): container %d not found",GetName(),c));
 355 }
 356
 357 //________________________________________________________________________
 358 void AliAnalysisTaskEmcalJet::SetJetAreaCut(Float_t cut, Int_t c)
 359 {
 360   AliJetContainer *cont = GetJetContainer(c);
 361   if (cont) cont->SetJetAreaCut(cut);
 362   else AliError(Form("%s in SetJetAreaCut(...): container %d not found",GetName(),c));
 363 }
 364
 365 //________________________________________________________________________
 366 void AliAnalysisTaskEmcalJet::SetPercAreaCut(Float_t p, Int_t c)
 367 {
 368   AliJetContainer *cont = GetJetContainer(c);
 369   if (cont) cont->SetPercAreaCut(p);
 370   else AliError(Form("%s in SetPercAreaCut(...): container %d not found",GetName(),c));
 371 }
 372
 373 //________________________________________________________________________
 374 void AliAnalysisTaskEmcalJet::SetZLeadingCut(Float_t zemc, Float_t zch, Int_t c)
 375 {
 376   AliJetContainer *cont = GetJetContainer(c);
 377   if (cont) cont->SetZLeadingCut(zemc,zch);
 378   else AliError(Form("%s in SetZLeadingCut(...): container %d not found",GetName(),c));
 379 }
 380
 381 //________________________________________________________________________
 382 void AliAnalysisTaskEmcalJet::SetNEFCut(Float_t min, Float_t max, Int_t c)
 383 {
 384   AliJetContainer *cont = GetJetContainer(c);
 385   if (cont) cont->SetNEFCut(min,max);
 386   else AliError(Form("%s in SetNEFCut(...): container %d not found",GetName(),c));
 387 }
 388
 389 //________________________________________________________________________
 390 void AliAnalysisTaskEmcalJet::SetAreaEmcCut(Double_t a, Int_t c)
 391 {
 392   AliJetContainer *cont = GetJetContainer(c);
 393   if (cont) cont->SetAreaEmcCut(a);
 394   else AliError(Form("%s in SetAreaEmcCut(...): container %d not found",GetName(),c));
 395 }
 396
 397 //________________________________________________________________________
 398 void AliAnalysisTaskEmcalJet::SetJetPtCut(Float_t cut, Int_t c)
 399 {
 400   AliJetContainer *cont = GetJetContainer(c);
 401   if (cont) cont->SetJetPtCut(cut);
 402   else AliError(Form("%s in SetJetPtCut(...): container %d not found",GetName(),c));
 403 }
 404
 405 //________________________________________________________________________
 406 void AliAnalysisTaskEmcalJet::SetJetRadius(Float_t r, Int_t c)
 407 {
 408   AliJetContainer *cont = GetJetContainer(c);
 409   if (cont) cont->SetJetRadius(r);
 410   else AliError(Form("%s in SetJetRadius(...): container %d not found",GetName(),c));
 411 }
 412
 413 //________________________________________________________________________
 414 void AliAnalysisTaskEmcalJet::SetMaxClusterPt(Float_t cut, Int_t c)
 415 {
 416   AliJetContainer *cont = GetJetContainer(c);
 417   if (cont) cont->SetMaxClusterPt(cut);
 418   else AliError(Form("%s in SetMaxClusterPt(...): container %d not found",GetName(),c));
 419 }
 420
 421 //________________________________________________________________________
 422 void AliAnalysisTaskEmcalJet::SetMaxTrackPt(Float_t cut, Int_t c)
 423 {
 424   AliJetContainer *cont = GetJetContainer(c);
 425   if (cont) cont->SetMaxTrackPt(cut);
 426   else AliError(Form("%s in SetMaxTrackPt(...): container %d not found",GetName(),c));
 427 }
 428
 429 //________________________________________________________________________
 430 void AliAnalysisTaskEmcalJet::SetPtBiasJetClus(Float_t cut, Int_t c)
 431 {
 432   AliJetContainer *cont = GetJetContainer(c);
 433   if (cont) cont->SetPtBiasJetClus(cut);
 434   else AliError(Form("%s in SetPtBiasJetClus(...): container %d not found",GetName(),c));
 435 }
 436
 437 //________________________________________________________________________
 438 void AliAnalysisTaskEmcalJet::SetPtBiasJetTrack(Float_t cut, Int_t c)
 439 {
 440   AliJetContainer *cont = GetJetContainer(c);
 441   if (cont) cont->SetPtBiasJetTrack(cut);
 442   else AliError(Form("%s in SetPtBiasJetTrack(...): container %d not found",GetName(),c));
 443 }
 444
 445 //________________________________________________________________________
 446 void AliAnalysisTaskEmcalJet::SetLeadingHadronType(Int_t t, Int_t c)
 447 {
 448   AliJetContainer *cont = GetJetContainer(c);
 449   if (cont) cont->SetLeadingHadronType(t);
 450   else AliError(Form("%s in SetLeadingHadronType(...): container %d not found",GetName(),c));
 451 }
 452
 453 //________________________________________________________________________
 454 void AliAnalysisTaskEmcalJet::SetNLeadingJets(Int_t t, Int_t c)
 455 {
 456   AliJetContainer *cont = GetJetContainer(c);
 457   if (cont) cont->SetNLeadingJets(t);
 458   else AliError(Form("%s in SetNLeadingJets(...): container %d not found",GetName(),c));
 459 }
 460
 461 //________________________________________________________________________
 462 void AliAnalysisTaskEmcalJet::SetJetBitMap(UInt_t m, Int_t c)
 463 {
 464   AliJetContainer *cont = GetJetContainer(c);
 465   if (cont) cont->SetJetBitMap(m);
 466   else AliError(Form("%s in SetJetBitMap(...): container %d not found",GetName(),c));
 467 }
 468
 469 //________________________________________________________________________
 470 void AliAnalysisTaskEmcalJet::SetIsParticleLevel(Bool_t b, Int_t c)
 471 {
 472   AliJetContainer *cont = GetJetContainer(c);
 473   if (cont) cont->SetIsParticleLevel(b);
 474   else AliError(Form("%s in SetIsParticleLevel(...): container %d not found",GetName(),c));
 475 }
 476
 477 //________________________________________________________________________
 478 const TString& AliAnalysisTaskEmcalJet::GetRhoName(Int_t c) const
 479 {
 480   if (c >= 0) {
 481     AliJetContainer *cont = GetJetContainer(c);
 482     if (cont) return cont->GetRhoName();
 483     else { AliError(Form("%s in GetRhoName(...): container %d not found. Returning fRhoName...",GetName(),c)); return fRhoName; }
 484   }
 485   else {
 486     return fRhoName;
 487   }
 488 }
 489
 490 //________________________________________________________________________
 491 TClonesArray* AliAnalysisTaskEmcalJet::GetJetArray(Int_t i) const {
 492   // Get i^th TClonesArray with AliEmcalJet
 493
 494   AliJetContainer *cont = GetJetContainer(i);
 495   if(!cont) {
 496     AliError(Form("%s:Container %d not found",GetName(),i));
 497     return 0;
 498   }
 499   return cont->GetArray();
 500 }
 501
 502 //________________________________________________________________________
 503 Double_t AliAnalysisTaskEmcalJet::GetJetRadius(Int_t i) const {
 504   // Get jet radius from jet container i
 505
 506   AliJetContainer *cont = GetJetContainer(i);
 507   if(!cont) {
 508     AliError(Form("%s:Container %d not found",GetName(),i));
 509     return 0;
 510   }
 511
 512   return cont->GetJetRadius();
 513 }
 514
 515 //________________________________________________________________________
 516 AliEmcalJet* AliAnalysisTaskEmcalJet::GetJetFromArray(Int_t j, Int_t c) const {
 517   // Get jet j if accepted from  container c
 518   // If jet not accepted return 0
 519
 520   AliJetContainer *cont = GetJetContainer(c);
 521   if(!cont) {
 522     AliError(Form("%s:Container %d not found",GetName(),c));
 523     return 0;
 524   }
 525   AliEmcalJet *jet = cont->GetJet(j);
 526
 527   return jet;
 528 }
 529
 530 //________________________________________________________________________
 531 AliEmcalJet* AliAnalysisTaskEmcalJet::GetAcceptJetFromArray(Int_t j, Int_t c) const {
 532   // Get jet j if accepted from  container c
 533   // If jet not accepted return 0
 534
 535   AliJetContainer *cont = GetJetContainer(c);
 536   if(!cont) {
 537     AliError(Form("%s:Container %d not found",GetName(),c));
 538     return 0;
 539   }
 540   AliEmcalJet *jet = cont->GetAcceptJet(j);
 541
 542   return jet;
 543 }
 544
 545 //________________________________________________________________________
 546 Int_t AliAnalysisTaskEmcalJet::GetNJets(Int_t i) const {
 547   // Get number of entries in jet array i
 548
 549   AliJetContainer *cont = GetJetContainer(i);
 550   if(!cont) {
 551     AliError(Form("%s:Container %d not found",GetName(),i));
 552     return 0;
 553   }
 554   return cont->GetNJets();
 555
 556 }
 557
 558 //________________________________________________________________________
 559 Double_t AliAnalysisTaskEmcalJet::GetRhoVal(Int_t i) const {
 560   // Get rho corresponding to jet array i
 561
 562   AliJetContainer *cont = GetJetContainer(i);
 563   if(!cont) {
 564     AliError(Form("%s:Container %d not found",GetName(),i));
 565     return 0;
 566   }
 567   return cont->GetRhoVal();
 568 }