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