PWGJE/EMCALJetTasks/AliEmcalJetTask.cxx

   1 // $Id$
   2 //
   3 // Emcal jet finder task.
   4 //
   5 // Authors: C.Loizides, S.Aiola
   6
   7 #include <vector>
   8 #include "AliEmcalJetTask.h"
   9
  10 #include <TChain.h>
  11 #include <TClonesArray.h>
  12 #include <TList.h>
  13 #include <TLorentzVector.h>
  14 #include <TMath.h>
  15
  16 #include "AliAnalysisManager.h"
  17 #include "AliCentrality.h"
  18 #include "AliEMCALGeometry.h"
  19 #include "AliESDEvent.h"
  20 #include "AliEmcalJet.h"
  21 #include "AliEmcalParticle.h"
  22 #include "AliFJWrapper.h"
  23 #include "AliMCEvent.h"
  24 #include "AliVCluster.h"
  25 #include "AliVEvent.h"
  26 #include "AliVParticle.h"
  27
  28 ClassImp(AliEmcalJetTask)
  29
  30 //________________________________________________________________________
  31 AliEmcalJetTask::AliEmcalJetTask() :
  32   AliAnalysisTaskSE("AliEmcalJetTask"),
  33   fTracksName("Tracks"),
  34   fCaloName("CaloClusters"),
  35   fJetsName("Jets"),
  36   fJetType(kNone),
  37   fConstSel(kAllJets),
  38   fMCConstSel(kAllJets),
  39   fMarkConst(kFALSE),
  40   fRadius(0.4),
  41   fMinJetTrackPt(0.15),
  42   fMinJetClusPt(0.15),
  43   fPhiMin(0),
  44   fPhiMax(TMath::TwoPi()),
  45   fEtaMin(-1),
  46   fEtaMax(1),
  47   fMinJetArea(0.01),
  48   fMinJetPt(1.0),
  49   fJetPhiMin(-10),
  50   fJetPhiMax(10),
  51   fJetEtaMin(-1),
  52   fJetEtaMax(1),
  53   fGhostArea(0.01),
  54   fIsInit(0),
  55   fIsPSelSet(0),
  56   fIsMcPart(0),
  57   fIsEmcPart(0),
  58   fJets(0),
  59   fEvent(0),
  60   fTracks(0),
  61   fClus(0)
  62 {
  63   // Default constructor.
  64 }
  65
  66 //________________________________________________________________________
  67 AliEmcalJetTask::AliEmcalJetTask(const char *name) :
  68   AliAnalysisTaskSE(name),
  69   fTracksName("Tracks"),
  70   fCaloName("CaloClusters"),
  71   fJetsName("Jets"),
  72   fJetType(kAKT|kFullJet|kRX1Jet),
  73   fConstSel(kAllJets),
  74   fMCConstSel(kAllJets),
  75   fMarkConst(kFALSE),
  76   fRadius(0.4),
  77   fMinJetTrackPt(0.15),
  78   fMinJetClusPt(0.15),
  79   fPhiMin(0),
  80   fPhiMax(TMath::TwoPi()),
  81   fEtaMin(-1),
  82   fEtaMax(1),
  83   fMinJetArea(0.01),
  84   fMinJetPt(1.0),
  85   fJetPhiMin(-10),
  86   fJetPhiMax(10),
  87   fJetEtaMin(-1),
  88   fJetEtaMax(1),
  89   fGhostArea(0.01),
  90   fIsInit(0),
  91   fIsPSelSet(0),
  92   fIsMcPart(0),
  93   fIsEmcPart(0),
  94   fJets(0),
  95   fEvent(0),
  96   fTracks(0),
  97   fClus(0)
  98 {
  99   // Standard constructor.
 100
 101   fBranchNames="ESD:AliESDRun.,AliESDHeader.,PrimaryVertex.";
 102 }
 103
 104 //________________________________________________________________________
 105 AliEmcalJetTask::~AliEmcalJetTask()
 106 {
 107   // Destructor
 108 }
 109
 110 //________________________________________________________________________
 111 void AliEmcalJetTask::UserCreateOutputObjects()
 112 {
 113   // Create user objects.
 114
 115   fJets = new TClonesArray("AliEmcalJet");
 116   fJets->SetName(fJetsName);
 117 }
 118
 119 //________________________________________________________________________
 120 void AliEmcalJetTask::UserExec(Option_t *)
 121 {
 122   // Main loop, called for each event.
 123
 124   if (!fIsInit) {
 125     if (!DoInit())
 126       return;
 127     fIsInit = kTRUE;
 128   }
 129
 130   FindJets();
 131 }
 132
 133 //________________________________________________________________________
 134 void AliEmcalJetTask::Terminate(Option_t *)
 135 {
 136   // Called once at the end of the analysis.
 137 }
 138
 139 //________________________________________________________________________
 140 void AliEmcalJetTask::FindJets()
 141 {
 142   // Find jets.
 143
 144   if (!fTracks && !fClus)
 145     return;
 146
 147   TString name("kt");
 148   fastjet::JetAlgorithm jalgo(fastjet::kt_algorithm);
 149   if ((fJetType & kAKT) != 0) {
 150     name  = "antikt";
 151     jalgo = fastjet::antikt_algorithm;
 152     AliDebug(1,"Using AKT algorithm");
 153   }
 154   else {
 155     AliDebug(1,"Using KT algorithm");
 156   }
 157
 158   if ((fJetType & kR020Jet) != 0)
 159     fRadius = 0.2;
 160   else if ((fJetType & kR030Jet) != 0)
 161     fRadius = 0.3;
 162   else if ((fJetType & kR040Jet) != 0)
 163     fRadius = 0.4;
 164
 165   // setup fj wrapper
 166   AliFJWrapper fjw(name, name);
 167   fjw.SetAreaType(fastjet::active_area_explicit_ghosts);
 168   fjw.SetGhostArea(fGhostArea);
 169   fjw.SetR(fRadius);
 170   fjw.SetAlgorithm(jalgo);
 171   fjw.SetMaxRap(1);
 172   fjw.Clear();
 173
 174   // get primary vertex
 175   Double_t vertex[3] = {0, 0, 0};
 176   fEvent->GetPrimaryVertex()->GetXYZ(vertex);
 177
 178   AliDebug(2,Form("Jet type = %d", fJetType));
 179
 180   if ((fIsMcPart || ((fJetType & kFullJet) != 0) || ((fJetType & kChargedJet) != 0)) && fTracks) {
 181     const Int_t Ntracks = fTracks->GetEntries();
 182     for (Int_t iTracks = 0; iTracks < Ntracks; ++iTracks) {
 183       AliVParticle *t = static_cast<AliVParticle*>(fTracks->At(iTracks));
 184       if (!t)
 185         continue;
 186       if (fIsMcPart) {
 187         if (((fJetType & kChargedJet) != 0) && (t->Charge() == 0)) {
 188           AliDebug(2,Form("Skipping track %d because it is neutral.", iTracks));
 189           continue;
 190         }
 191         if (((fJetType & kNeutralJet) != 0) && (t->Charge() != 0)) {
 192           AliDebug(2,Form("Skipping track %d because it is charged.", iTracks));
 193           continue;
 194         }
 195       }
 196       if (fIsMcPart || t->GetLabel() != 0) {
 197         if (fMCConstSel == kNone) {
 198           AliDebug(2,Form("Skipping track %d because bit mask is 0.", iTracks));
 199           continue;
 200         }
 201         if (fMCConstSel != kAllJets) {
 202           if (t->TestBits(fMCConstSel) != (Int_t)fMCConstSel) {
 203             AliDebug(2,Form("Skipping track %d because it does not match the bit mask (%d, %d)", iTracks, fMCConstSel, t->TestBits(TObject::kBitMask)));
 204             continue;
 205           }
 206           else {
 207             AliDebug(2,Form("Track %d matches the bit mask (%d, %d)", iTracks, fMCConstSel, t->TestBits(TObject::kBitMask)));
 208           }
 209         }
 210       }
 211       else {
 212         if (fConstSel == kNone) {
 213           AliDebug(2,Form("Skipping track %d because bit mask is 0.", iTracks));
 214           continue;
 215         }
 216         if (fConstSel != kAllJets) {
 217           if (t->TestBits(fConstSel) != (Int_t)fConstSel) {
 218             AliDebug(2,Form("Skipping track %d because it does not match the bit mask (%d, %d)", iTracks, fConstSel, t->TestBits(TObject::kBitMask)));
 219             continue;
 220           }
 221           else {
 222             AliDebug(2,Form("Track %d matches the bit mask (%d, %d)", iTracks, fConstSel, t->TestBits(TObject::kBitMask)));
 223           }
 224         }
 225       }
 226       if (t->Pt() < fMinJetTrackPt)
 227         continue;
 228       Double_t eta = t->Eta();
 229       Double_t phi = t->Phi();
 230       if ((eta<fEtaMin) || (eta>fEtaMax) ||
 231           (phi<fPhiMin) || (phi>fPhiMax))
 232         continue;
 233
 234       // offset of 100 for consistency with cluster ids
 235       AliDebug(2,Form("Track %d accepted (label = %d, pt = %f)", iTracks, t->GetLabel(), t->Pt()));
 236       fjw.AddInputVector(t->Px(), t->Py(), t->Pz(), t->P(), iTracks + 100);
 237     }
 238   }
 239
 240   if ((((fJetType & kFullJet) != 0) || ((fJetType & kNeutralJet) != 0)) && fClus) {
 241     const Int_t Nclus = fClus->GetEntries();
 242     for (Int_t iClus = 0; iClus < Nclus; ++iClus) {
 243       AliVCluster *c = 0;
 244       Double_t cEta=0,cPhi=0,cPt=0;
 245       Double_t cPx=0,cPy=0,cPz=0;
 246       if (fIsEmcPart) {
 247         AliEmcalParticle *ep = static_cast<AliEmcalParticle*>(fClus->At(iClus));
 248         if (!ep)
 249           continue;
 250
 251         c = ep->GetCluster();
 252         if (!c)
 253           continue;
 254
 255         if (c->GetLabel() > 0) {
 256           if (fMCConstSel == kNone) {
 257             AliDebug(2,Form("Skipping cluster %d because bit mask is 0.", iClus));
 258             continue;
 259           }
 260           if (fMCConstSel != kAllJets) {
 261             if (ep->TestBits(fMCConstSel) != (Int_t)fMCConstSel) {
 262               AliDebug(2,Form("Skipping cluster %d because it does not match the bit mask (%d, %d)", iClus, fMCConstSel, ep->TestBits(TObject::kBitMask)));
 263               continue;
 264             }
 265             else {
 266               AliDebug(2,Form("Cluster %d matches the bit mask (%d, %d)", iClus, fMCConstSel, ep->TestBits(TObject::kBitMask)));
 267             }
 268           }
 269         }
 270         else {
 271           if (fConstSel == kNone) {
 272             AliDebug(2,Form("Skipping cluster %d because bit mask is 0.", iClus));
 273             continue;
 274           }
 275           if (fConstSel != kAllJets) {
 276             if (ep->TestBits(fConstSel) != (Int_t)fConstSel) {
 277               AliDebug(2,Form("Skipping cluster %d because it does not match the bit mask (%d, %d)", iClus, fConstSel, ep->TestBits(TObject::kBitMask)));
 278               continue;
 279             }
 280             else {
 281               AliDebug(2,Form("Cluster %d matches the bit mask (%d, %d)", iClus, fConstSel, ep->TestBits(TObject::kBitMask)));
 282             }
 283           }
 284         }
 285
 286         cEta = ep->Eta();
 287         cPhi = ep->Phi();
 288         cPt  = ep->Pt();
 289         cPx  = ep->Px();
 290         cPy  = ep->Py();
 291         cPz  = ep->Pz();
 292       } else {
 293         c = static_cast<AliVCluster*>(fClus->At(iClus));
 294         if (!c)
 295           continue;
 296
 297         if (c->GetLabel() > 0) {
 298           if (fMCConstSel == kNone) {
 299             AliDebug(2,Form("Skipping cluster %d because bit mask is 0.", iClus));
 300             continue;
 301           }
 302           if (fMCConstSel != kAllJets) {
 303             if (c->TestBits(fMCConstSel) != (Int_t)fMCConstSel) {
 304               AliDebug(2,Form("Skipping cluster %d because it does not match the bit mask (%d, %d)", iClus, fMCConstSel, c->TestBits(TObject::kBitMask)));
 305               continue;
 306             }
 307             else {
 308               AliDebug(2,Form("Cluster %d matches the bit mask (%d, %d)", iClus, fMCConstSel, c->TestBits(TObject::kBitMask)));
 309             }
 310           }
 311         }
 312         else {
 313           if (fConstSel == kNone) {
 314             AliDebug(2,Form("Skipping cluster %d because bit mask is 0.", iClus));
 315             continue;
 316           }
 317           if (fConstSel != kAllJets) {
 318             if (c->TestBits(fConstSel) != (Int_t)fConstSel) {
 319               AliDebug(2,Form("Skipping cluster %d because it does not match the bit mask (%d, %d)", iClus, fConstSel, c->TestBits(TObject::kBitMask)));
 320               continue;
 321             }
 322             else {
 323               AliDebug(2,Form("Cluster %d matches the bit mask (%d, %d)", iClus, fConstSel, c->TestBits(TObject::kBitMask)));
 324             }
 325           }
 326         }
 327
 328         TLorentzVector nP;
 329         c->GetMomentum(nP, vertex);
 330         cEta = nP.Eta();
 331         cPhi = nP.Phi();
 332         cPt  = nP.Pt();
 333         cPx  = nP.Px();
 334         cPy  = nP.Py();
 335         cPz  = nP.Pz();
 336       }
 337       if (!c->IsEMCAL())
 338         continue;
 339       if (cPt < fMinJetClusPt)
 340         continue;
 341       if ((cEta<fEtaMin) || (cEta>fEtaMax) ||
 342           (cPhi<fPhiMin) || (cPhi>fPhiMax))
 343         continue;
 344       // offset of 100 to skip ghost particles uid = -1
 345       AliDebug(2,Form("Cluster %d accepted (label = %d)", iClus, c->GetLabel()));
 346       fjw.AddInputVector(cPx, cPy, cPz, TMath::Sqrt(cPx*cPx+cPy*cPy+cPz*cPz), -iClus - 100);
 347     }
 348   }
 349
 350   // run jet finder
 351   fjw.Run();
 352
 353   // get geometry
 354   AliEMCALGeometry *geom = AliEMCALGeometry::GetInstance();
 355   if (!geom) {
 356     AliFatal(Form("%s: Can not create geometry", GetName()));
 357     return;
 358   }
 359
 360   // loop over fastjet jets
 361   std::vector<fastjet::PseudoJet> jets_incl = fjw.GetInclusiveJets();
 362   // sort jets according to jet pt
 363   static Int_t indexes[9999] = {-1};
 364   GetSortedArray(indexes, jets_incl);
 365
 366   AliDebug(1,Form("%d jets found", (Int_t)jets_incl.size()));
 367   for (UInt_t ijet=0, jetCount=0; ijet<jets_incl.size(); ++ijet) {
 368     Int_t ij = indexes[ijet];
 369     AliDebug(3,Form("Jet pt = %f, area = %f", jets_incl[ij].perp(), fjw.GetJetArea(ij)));
 370
 371     if (jets_incl[ij].perp()<fMinJetPt)
 372       continue;
 373     if (fjw.GetJetArea(ij)<fMinJetArea)
 374       continue;
 375     if ((jets_incl[ij].eta()<fJetEtaMin) || (jets_incl[ij].eta()>fJetEtaMax) ||
 376         (jets_incl[ij].phi()<fJetPhiMin) || (jets_incl[ij].phi()>fJetPhiMax))
 377       continue;
 378
 379     AliEmcalJet *jet = new ((*fJets)[jetCount])
 380       AliEmcalJet(jets_incl[ij].perp(), jets_incl[ij].eta(), jets_incl[ij].phi(), jets_incl[ij].m());
 381
 382     // loop over constituents
 383     std::vector<fastjet::PseudoJet> constituents(fjw.GetJetConstituents(ij));
 384     jet->SetNumberOfTracks(constituents.size());
 385     jet->SetNumberOfClusters(constituents.size());
 386
 387     Int_t nt            = 0;
 388     Int_t nc            = 0;
 389     Double_t neutralE   = 0;
 390     Double_t maxCh      = 0;
 391     Double_t maxNe      = 0;
 392     Int_t gall          = 0;
 393     Int_t gemc          = 0;
 394     Int_t cemc          = 0;
 395     Int_t ncharged      = 0;
 396     Int_t nneutral      = 0;
 397     Double_t mcpt       = 0;
 398     Double_t emcpt      = 0;
 399
 400     for(UInt_t ic = 0; ic < constituents.size(); ++ic) {
 401       Int_t uid = constituents[ic].user_index();
 402       AliDebug(3,Form("Processing constituent %d", uid));
 403       if ((uid == -1) /*&& (constituents[ic].kt2() < 1e-25)*/) { //ghost particle
 404         ++gall;
 405         Double_t gphi = constituents[ic].phi();
 406         if (gphi<0)
 407           gphi += TMath::TwoPi();
 408         gphi *= TMath::RadToDeg();
 409         Double_t geta = constituents[ic].eta();
 410         if ((gphi > geom->GetArm1PhiMin()) && (gphi < geom->GetArm1PhiMax()) &&
 411             (geta > geom->GetArm1EtaMin()) && (geta < geom->GetArm1EtaMax()))
 412           ++gemc;
 413       } else if ((uid > 0) && fTracks) { // track constituent
 414         Int_t tid = uid - 100;
 415         AliVParticle *t = static_cast<AliVParticle*>(fTracks->At(tid));
 416         if (!t) {
 417           AliError(Form("Could not find track %d",tid));
 418           continue;
 419         }
 420         if (jetCount < fMarkConst) {
 421           AliDebug(2,Form("Marking track %d with bit map %d", tid, fJetType));
 422           t->SetBit(fJetType);
 423         }
 424         Double_t cEta = t->Eta();
 425         Double_t cPhi = t->Phi();
 426         Double_t cPt  = t->Pt();
 427         Double_t cP   = t->P();
 428         if (t->Charge() == 0) {
 429           neutralE += cP;
 430           ++nneutral;
 431           if (cPt > maxNe)
 432             maxNe = cPt;
 433         } else {
 434           ++ncharged;
 435           if (cPt > maxCh)
 436             maxCh = cPt;
 437         }
 438         if (fIsMcPart || t->GetLabel() != 0) // check if MC particle
 439           mcpt += cPt;
 440
 441         if (cPhi<0)
 442           cPhi += TMath::TwoPi();
 443         cPhi *= TMath::RadToDeg();
 444         if ((cPhi > geom->GetArm1PhiMin()) && (cPhi < geom->GetArm1PhiMax()) &&
 445             (cEta > geom->GetArm1EtaMin()) && (cEta < geom->GetArm1EtaMax())) {
 446           emcpt += cPt;
 447           ++cemc;
 448         }
 449
 450         jet->AddTrackAt(tid, nt);
 451         ++nt;
 452       } else if (fClus) { // cluster constituent
 453         Int_t cid = -uid - 100;
 454         AliVCluster *c = 0;
 455         Double_t cEta=0,cPhi=0,cPt=0,cP=0;
 456         if (fIsEmcPart) {
 457           AliEmcalParticle *ep = static_cast<AliEmcalParticle*>(fClus->At(cid));
 458           if (!ep)
 459             continue;
 460           c = ep->GetCluster();
 461           if (!c)
 462             continue;
 463           if (jetCount < fMarkConst)
 464             ep->SetBit(fJetType);
 465           cEta = ep->Eta();
 466           cPhi = ep->Phi();
 467           cPt  = ep->Pt();
 468           cP   = ep->P();
 469         } else {
 470           c = static_cast<AliVCluster*>(fClus->At(cid));
 471           if (!c)
 472             continue;
 473           if (jetCount < fMarkConst)
 474             c->SetBit(fJetType);
 475           TLorentzVector nP;
 476           c->GetMomentum(nP, vertex);
 477           cEta = nP.Eta();
 478           cPhi = nP.Phi();
 479           cPt  = nP.Pt();
 480           cP   = nP.P();
 481         }
 482
 483         neutralE += cP;
 484         if (cPt > maxNe)
 485           maxNe = cPt;
 486
 487         if (c->GetLabel() > 0) // MC particle
 488           mcpt += cPt;
 489
 490         if (cPhi<0)
 491           cPhi += TMath::TwoPi();
 492         cPhi *= TMath::RadToDeg();
 493         if ((cPhi > geom->GetArm1PhiMin()) && (cPhi < geom->GetArm1PhiMax()) &&
 494             (cEta > geom->GetArm1EtaMin()) && (cEta < geom->GetArm1EtaMax())) {
 495           emcpt += cPt;
 496           ++cemc;
 497         }
 498
 499         jet->AddClusterAt(cid, nc);
 500         ++nc;
 501         ++nneutral;
 502       } else {
 503         AliError(Form("%s: No logical way to end up here.", GetName()));
 504         continue;
 505       }
 506     } /* end of constituent loop */
 507
 508     jet->SetNumberOfTracks(nt);
 509     jet->SetNumberOfClusters(nc);
 510     jet->SortConstituents();
 511     jet->SetMaxNeutralPt(maxNe);
 512     jet->SetMaxChargedPt(maxCh);
 513     jet->SetNEF(neutralE / jet->E());
 514     fastjet::PseudoJet area(fjw.GetJetAreaVector(ij));
 515     jet->SetArea(area.perp());
 516     jet->SetAreaEta(area.eta());
 517     jet->SetAreaPhi(area.phi());
 518     jet->SetNumberOfCharged(ncharged);
 519     jet->SetNumberOfNeutrals(nneutral);
 520     jet->SetMCPt(mcpt);
 521     jet->SetNEmc(cemc);
 522     jet->SetPtEmc(emcpt);
 523
 524     if (gall > 0)
 525       jet->SetAreaEmc(fjw.GetJetArea(ij) * gemc / gall);
 526     else
 527       jet->SetAreaEmc(-1);
 528     if ((jet->Phi() > geom->GetArm1PhiMin() * TMath::DegToRad()) &&
 529         (jet->Phi() < geom->GetArm1PhiMax() * TMath::DegToRad()) &&
 530         (jet->Eta() > geom->GetArm1EtaMin()) &&
 531         (jet->Eta() < geom->GetArm1EtaMax()))
 532       jet->SetAxisInEmcal(kTRUE);
 533
 534     AliDebug(2,Form("Added jet n. %d, pt = %f, area = %f, constituents = %d", jetCount, jet->Pt(), jet->Area(), (Int_t)constituents.size()));
 535     jetCount++;
 536   }
 537   //fJets->Sort();
 538 }
 539
 540 //________________________________________________________________________
 541 Bool_t AliEmcalJetTask::GetSortedArray(Int_t indexes[], std::vector<fastjet::PseudoJet> array) const
 542 {
 543   // Get the leading jets.
 544
 545   static Float_t pt[9999] = {0};
 546
 547   const Int_t n = (Int_t)array.size();
 548
 549   if (n < 1)
 550     return kFALSE;
 551
 552   for (Int_t i = 0; i < n; i++)
 553     pt[i] = array[i].perp();
 554
 555   TMath::Sort(n, pt, indexes);
 556
 557   return kTRUE;
 558 }
 559
 560 //________________________________________________________________________
 561 Bool_t AliEmcalJetTask::DoInit()
 562 {
 563   // Init. Return true if successful.
 564
 565   // get input collections
 566   AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager();
 567
 568   // get the event
 569   fEvent = InputEvent();
 570   if (!fEvent) {
 571     AliError(Form("%s: Could not retrieve event! Returning", GetName()));
 572     return 0;
 573   }
 574
 575   // add jets to event if not yet there
 576   fJets->Delete();
 577   if (!(fEvent->FindListObject(fJetsName)))
 578     fEvent->AddObject(fJets);
 579   else {
 580     AliError(Form("%s: Object with name %s already in event! Returning", GetName(), fJetsName.Data()));
 581     return 0;
 582   }
 583
 584   if (fTracksName == "Tracks")
 585     am->LoadBranch("Tracks");
 586   if (!fTracks && !fTracksName.IsNull()) {
 587     fTracks = dynamic_cast<TClonesArray*>(fEvent->FindListObject(fTracksName));
 588     if (!fTracks) {
 589       AliError(Form("%s: Pointer to tracks %s == 0", GetName(), fTracksName.Data()));
 590       return 0;
 591     }
 592   }
 593   if (fTracks) {
 594     TClass cls(fTracks->GetClass()->GetName());
 595     if (cls.InheritsFrom("AliMCParticle") || cls.InheritsFrom("AliAODMCParticle"))
 596       fIsMcPart = 1;
 597   }
 598
 599   if (fCaloName == "CaloClusters")
 600     am->LoadBranch("CaloClusters");
 601   if (!fClus && !fCaloName.IsNull()) {
 602     fClus = dynamic_cast<TClonesArray*>(fEvent->FindListObject(fCaloName));
 603     if (!fClus) {
 604       AliError(Form("%s: Pointer to clus %s == 0", GetName(), fCaloName.Data()));
 605       return 0;
 606     }
 607   }
 608   if (fClus) {
 609     TClass cls(fClus->GetClass()->GetName());
 610     if (cls.InheritsFrom("AliEmcalParticle"))
 611       fIsEmcPart = 1;
 612   }
 613
 614   return 1;
 615 }