PWGJE/EMCALJetTasks/AliJetContainer.cxx

   1 // $Id$
   2 //
   3 // Container with name, TClonesArray and cuts for jets
   4 //
   5 // Author: M. Verweij
   6
   7 #include <TClonesArray.h>
   8
   9 #include "AliEmcalJet.h"
  10 #include "AliVEvent.h"
  11 #include "AliLog.h"
  12 #include "AliEMCALGeometry.h"
  13 #include "AliParticleContainer.h"
  14 #include "AliClusterContainer.h"
  15 #include "AliLocalRhoParameter.h"
  16 #include "AliPythiaInfo.h"
  17
  18 #include "AliJetContainer.h"
  19
  20 ClassImp(AliJetContainer)
  21
  22 //________________________________________________________________________
  23 AliJetContainer::AliJetContainer():
  24   AliEmcalContainer("AliJetContainer"),
  25   fJetAcceptanceType(kUser),
  26   fJetRadius(0),
  27   fRhoName(),
  28   fLocalRhoName(),
  29   fRhoMassName(),
  30   fPythiaInfoName(),
  31   fFlavourSelection(0),
  32   fPtBiasJetTrack(0),
  33   fPtBiasJetClus(0),
  34   fJetPtCut(1),
  35   fJetAreaCut(-1),
  36   fAreaEmcCut(0),
  37   fJetMinEta(-0.9),
  38   fJetMaxEta(0.9),
  39   fJetMinPhi(-10),
  40   fJetMaxPhi(10),
  41   fMaxClusterPt(1000),
  42   fMaxTrackPt(100),
  43   fZLeadingEmcCut(10.),
  44   fZLeadingChCut(10.),
  45   fNEFMinCut(-10.),
  46   fNEFMaxCut(10.),
  47   fLeadingHadronType(0),
  48   fNLeadingJets(1),
  49   fJetBitMap(0),
  50   fJetTrigger(0),
  51   fTagStatus(-1),
  52   fParticleContainer(0),
  53   fClusterContainer(0),
  54   fRho(0),
  55   fLocalRho(0),
  56   fRhoMass(0),
  57   fPythiaInfo(0),
  58   fGeom(0),
  59   fRunNumber(0)
  60 {
  61   // Default constructor.
  62
  63   fClassName = "AliEmcalJet";
  64 }
  65
  66 //________________________________________________________________________
  67 AliJetContainer::AliJetContainer(const char *name):
  68   AliEmcalContainer(name),
  69   fJetAcceptanceType(kUser),
  70   fJetRadius(0),
  71   fRhoName(),
  72   fLocalRhoName(),
  73   fRhoMassName(),
  74   fPythiaInfoName(),
  75   fFlavourSelection(0),
  76   fPtBiasJetTrack(0),
  77   fPtBiasJetClus(0),
  78   fJetPtCut(1),
  79   fJetAreaCut(-1),
  80   fAreaEmcCut(0),
  81   fJetMinEta(-0.9),
  82   fJetMaxEta(0.9),
  83   fJetMinPhi(-10),
  84   fJetMaxPhi(10),
  85   fMaxClusterPt(1000),
  86   fMaxTrackPt(100),
  87   fZLeadingEmcCut(10.),
  88   fZLeadingChCut(10.),
  89   fNEFMinCut(-10.),
  90   fNEFMaxCut(10.),
  91   fLeadingHadronType(0),
  92   fNLeadingJets(1),
  93   fJetBitMap(0),
  94   fJetTrigger(0),
  95   fTagStatus(-1),
  96   fParticleContainer(0),
  97   fClusterContainer(0),
  98   fRho(0),
  99   fLocalRho(0),
 100   fRhoMass(0),
 101   fPythiaInfo(0),
 102   fGeom(0),
 103   fRunNumber(0)
 104 {
 105   // Standard constructor.
 106
 107   fClassName = "AliEmcalJet";
 108 }
 109
 110 //________________________________________________________________________
 111 void AliJetContainer::SetArray(AliVEvent *event)
 112 {
 113   // Set jet array
 114
 115   AliEmcalContainer::SetArray(event);
 116
 117   if(fJetAcceptanceType==kTPC) {
 118     AliDebug(2,Form("%s: set TPC acceptance cuts",GetName()));
 119     SetJetEtaPhiTPC();
 120   }
 121   else if(fJetAcceptanceType==kEMCAL) {
 122     AliDebug(2,Form("%s: set EMCAL acceptance cuts",GetName()));
 123     SetJetEtaPhiEMCAL();
 124  }
 125 }
 126
 127
 128 //________________________________________________________________________
 129 void AliJetContainer::SetEMCALGeometry() {
 130   fGeom = AliEMCALGeometry::GetInstance();
 131   if (!fGeom) {
 132     AliError(Form("%s: Can not create geometry", GetName()));
 133     return;
 134   }
 135 }
 136
 137 //________________________________________________________________________
 138 void AliJetContainer::LoadRho(AliVEvent *event)
 139 {
 140   // Load rho
 141
 142   if (!fRhoName.IsNull() && !fRho) {
 143     fRho = dynamic_cast<AliRhoParameter*>(event->FindListObject(fRhoName));
 144     if (!fRho) {
 145       AliError(Form("%s: Could not retrieve rho %s!", GetName(), fRhoName.Data()));
 146       return;
 147     }
 148   }
 149 }
 150
 151 //________________________________________________________________________
 152 void AliJetContainer::LoadLocalRho(AliVEvent *event)
 153 {
 154   // Load local rho
 155
 156   if (!fLocalRhoName.IsNull() && !fLocalRho) {
 157     fLocalRho = dynamic_cast<AliLocalRhoParameter*>(event->FindListObject(fLocalRhoName));
 158     if (!fLocalRho) {
 159       AliError(Form("%s: Could not retrieve rho %s!", GetName(), fLocalRhoName.Data()));
 160       return;
 161     }
 162   }
 163 }
 164
 165 //________________________________________________________________________
 166 void AliJetContainer::LoadRhoMass(AliVEvent *event)
 167 {
 168   // Load rho
 169
 170   if (!fRhoMassName.IsNull() && !fRhoMass) {
 171     fRhoMass = dynamic_cast<AliRhoParameter*>(event->FindListObject(fRhoMassName));
 172     if (!fRhoMass) {
 173       AliError(Form("%s: Could not retrieve rho_mass %s!", GetName(), fRhoMassName.Data()));
 174       return;
 175     }
 176   }
 177 }
 178 //________________________________________________________________________
 179 void AliJetContainer::LoadPythiaInfo(AliVEvent *event)
 180 {
 181     // Load parton info
 182
 183     if (!fPythiaInfoName.IsNull() && !fPythiaInfo) {
 184         fPythiaInfo = dynamic_cast<AliPythiaInfo*>(event->FindListObject(fPythiaInfoName));
 185         if (!fPythiaInfo) {
 186            AliError(Form("%s: Could not retrieve parton infos! %s!", GetName(), fPythiaInfoName.Data()));
 187         return;
 188         }
 189     }
 190 }
 191
 192
 193
 194 //________________________________________________________________________
 195 AliEmcalJet* AliJetContainer::GetLeadingJet(const char* opt)
 196 {
 197   // Get the leading jet; if opt contains "rho" the sorting is according to pt-A*rho
 198
 199   TString option(opt);
 200   option.ToLower();
 201
 202   Int_t tempID = fCurrentID;
 203
 204   AliEmcalJet *jetMax = GetNextAcceptJet(0);
 205   AliEmcalJet *jet = 0;
 206
 207   if (option.Contains("rho")) {
 208     while ((jet = GetNextAcceptJet())) {
 209       if ( (jet->Pt()-jet->Area()*GetRhoVal()) > (jetMax->Pt()-jetMax->Area()*GetRhoVal()) )
 210         jetMax = jet;
 211     }
 212   }
 213   else {
 214     while ((jet = GetNextAcceptJet())) {
 215       if (jet->Pt() > jetMax->Pt()) jetMax = jet;
 216     }
 217   }
 218
 219   fCurrentID = tempID;
 220
 221   return jetMax;
 222 }
 223
 224 //________________________________________________________________________
 225 AliEmcalJet* AliJetContainer::GetJet(Int_t i) const {
 226
 227   //Get i^th jet in array
 228
 229   if(i<0 || i>fClArray->GetEntriesFast()) return 0;
 230   AliEmcalJet *jet = static_cast<AliEmcalJet*>(fClArray->At(i));
 231   return jet;
 232
 233 }
 234
 235 //________________________________________________________________________
 236 AliEmcalJet* AliJetContainer::GetAcceptJet(Int_t i) {
 237
 238   //Only return jet if is accepted
 239
 240   AliEmcalJet *jet = GetJet(i);
 241   if(!AcceptJet(jet)) return 0;
 242
 243   return jet;
 244 }
 245
 246 //________________________________________________________________________
 247 AliEmcalJet* AliJetContainer::GetJetWithLabel(Int_t lab) const {
 248
 249   //Get particle with label lab in array
 250
 251   Int_t i = GetIndexFromLabel(lab);
 252   return GetJet(i);
 253 }
 254
 255 //________________________________________________________________________
 256 AliEmcalJet* AliJetContainer::GetAcceptJetWithLabel(Int_t lab) {
 257
 258   //Get particle with label lab in array
 259
 260   Int_t i = GetIndexFromLabel(lab);
 261   return GetAcceptJet(i);
 262 }
 263
 264 //________________________________________________________________________
 265 AliEmcalJet* AliJetContainer::GetNextAcceptJet(Int_t i) {
 266
 267   //Get next accepted jet; if i >= 0 (re)start counter from i; return 0 if no accepted jet could be found
 268
 269   if (i>=0) fCurrentID = i;
 270
 271   const Int_t njets = GetNEntries();
 272   AliEmcalJet *jet = 0;
 273   while (fCurrentID < njets && !jet) {
 274     jet = GetAcceptJet(fCurrentID);
 275     fCurrentID++;
 276   }
 277
 278   return jet;
 279 }
 280
 281 //________________________________________________________________________
 282 AliEmcalJet* AliJetContainer::GetNextJet(Int_t i) {
 283
 284   //Get next jet; if i >= 0 (re)start counter from i; return 0 if no jet could be found
 285
 286   if (i>=0) fCurrentID = i;
 287
 288   const Int_t njets = GetNEntries();
 289   AliEmcalJet *jet = 0;
 290   while (fCurrentID < njets && !jet) {
 291     jet = GetJet(fCurrentID);
 292     fCurrentID++;
 293   }
 294
 295   return jet;
 296 }
 297
 298 //________________________________________________________________________
 299 Double_t AliJetContainer::GetJetPtCorr(Int_t i) const {
 300   AliEmcalJet *jet = GetJet(i);
 301
 302   return jet->Pt() - fRho->GetVal()*jet->Area();
 303 }
 304
 305 //________________________________________________________________________
 306 Double_t AliJetContainer::GetJetPtCorrLocal(Int_t i) const {
 307   AliEmcalJet *jet = GetJet(i);
 308
 309   return jet->Pt() - fLocalRho->GetLocalVal(jet->Phi(), fJetRadius)*jet->Area();
 310 }
 311
 312 //________________________________________________________________________
 313 void AliJetContainer::GetMomentum(TLorentzVector &mom, Int_t i) const
 314 {
 315   //Get momentum of the i^th jet in array
 316
 317   AliEmcalJet *jet = GetJet(i);
 318   if(jet) jet->GetMom(mom);
 319 }
 320
 321 //________________________________________________________________________
 322 Bool_t AliJetContainer::AcceptBiasJet(const AliEmcalJet *jet)
 323 {
 324   // Accept jet with a bias.
 325
 326   if (fLeadingHadronType == 0) {
 327     if (jet->MaxTrackPt() < fPtBiasJetTrack) return kFALSE;
 328   }
 329   else if (fLeadingHadronType == 1) {
 330     if (jet->MaxClusterPt() < fPtBiasJetClus) return kFALSE;
 331   }
 332   else {
 333     if (jet->MaxTrackPt() < fPtBiasJetTrack && jet->MaxClusterPt() < fPtBiasJetClus) return kFALSE;
 334   }
 335
 336   return kTRUE;
 337 }
 338
 339 //________________________________________________________________________
 340 Bool_t AliJetContainer::AcceptJet(const AliEmcalJet *jet)
 341 {
 342    // Return true if jet is accepted.
 343
 344   fRejectionReason = 0;
 345
 346   if (!jet) {
 347     AliDebug(11,"No jet found");
 348     fRejectionReason |= kNullObject;
 349     return kFALSE;
 350   }
 351
 352   if (jet->Pt() <= fJetPtCut) {
 353     AliDebug(11,Form("Cut rejecting jet: JetPtCut %.1f",fJetPtCut));
 354     fRejectionReason |= kPtCut;
 355     return kFALSE;
 356   }
 357
 358   Double_t jetPhi = jet->Phi();
 359   Double_t jetEta = jet->Eta();
 360
 361   // if limits are given in (-pi, pi) range
 362   if (fJetMinPhi < 0) jetPhi -= TMath::Pi() * 2;
 363
 364   if (jetEta < fJetMinEta || jetEta > fJetMaxEta || jetPhi < fJetMinPhi || jetPhi > fJetMaxPhi) {
 365     AliDebug(11,"Cut rejecting jet: Acceptance");
 366     fRejectionReason |= kAcceptanceCut;
 367     return kFALSE;
 368   }
 369
 370   if (jet->TestBits(fJetBitMap) != (Int_t)fJetBitMap) {
 371     AliDebug(11,"Cut rejecting jet: Bit map");
 372     fRejectionReason |= kBitMapCut;
 373     return kFALSE;
 374   }
 375
 376   if (jet->Area() <= fJetAreaCut)  {
 377     AliDebug(11,"Cut rejecting jet: Area");
 378     fRejectionReason |= kAreaCut;
 379     return kFALSE;
 380   }
 381
 382   if (jet->AreaEmc() < fAreaEmcCut) {
 383     AliDebug(11,"Cut rejecting jet: AreaEmc");
 384     fRejectionReason |= kAreaEmcCut;
 385     return kFALSE;
 386   }
 387
 388   if (fZLeadingChCut < 1 && GetZLeadingCharged(jet) > fZLeadingChCut) {
 389     AliDebug(11,"Cut rejecting jet: ZLeading");
 390     fRejectionReason |= kZLeadingChCut;
 391     return kFALSE;
 392   }
 393
 394   if (fZLeadingEmcCut < 1 && GetZLeadingEmc(jet) > fZLeadingEmcCut) {
 395     AliDebug(11,"Cut rejecting jet: ZLeadEmc");
 396     fRejectionReason |= kZLeadingEmcCut;
 397     return kFALSE;
 398   }
 399
 400   if (jet->NEF() < fNEFMinCut || jet->NEF() > fNEFMaxCut) {
 401     AliDebug(11,"Cut rejecting jet: NEF");
 402     fRejectionReason |= kNEFCut;
 403     return kFALSE;
 404   }
 405
 406   if (!AcceptBiasJet(jet)) {
 407     AliDebug(11,"Cut rejecting jet: Bias");
 408     fRejectionReason |= kMinLeadPtCut;
 409     return kFALSE;
 410   }
 411
 412   if (jet->MaxTrackPt() > fMaxTrackPt) {
 413     AliDebug(11,"Cut rejecting jet: MaxTrackPt");
 414     fRejectionReason |= kMaxTrackPtCut;
 415     return kFALSE;
 416
 417   }
 418
 419   if (jet->MaxClusterPt() > fMaxClusterPt) {
 420     AliDebug(11,"Cut rejecting jet: MaxClusPt");
 421     fRejectionReason |= kMaxClusterPtCut;
 422     return kFALSE;
 423   }
 424
 425   if (fFlavourSelection != 0 && !jet->TestFlavourTag(fFlavourSelection)) {
 426     AliDebug(11,"Cut rejecting jet: Flavour");
 427     fRejectionReason |= kFlavourCut;
 428     return kFALSE;
 429   }
 430
 431   if(fTagStatus>-1 && jet->GetTagStatus()!=fTagStatus) {
 432     AliDebug(11,"Cut rejecting jet: tag status");
 433     fRejectionReason |= kTagStatus;
 434     return kFALSE;
 435   }
 436
 437   return kTRUE;
 438 }
 439
 440 //________________________________________________________________________
 441 Double_t AliJetContainer::GetLeadingHadronPt(const AliEmcalJet *jet) const
 442 {
 443   if (fLeadingHadronType == 0)       // charged leading hadron
 444     return jet->MaxTrackPt();
 445   else if (fLeadingHadronType == 1)  // neutral leading hadron
 446     return jet->MaxClusterPt();
 447   else                               // charged or neutral
 448     return jet->MaxPartPt();
 449 }
 450
 451 //________________________________________________________________________
 452 void AliJetContainer::GetLeadingHadronMomentum(TLorentzVector &mom, const AliEmcalJet *jet) const
 453 {
 454   Double_t maxClusterPt = 0;
 455   Double_t maxClusterEta = 0;
 456   Double_t maxClusterPhi = 0;
 457
 458   Double_t maxTrackPt = 0;
 459   Double_t maxTrackEta = 0;
 460   Double_t maxTrackPhi = 0;
 461
 462   if (fClusterContainer && fClusterContainer->GetArray() && (fLeadingHadronType == 1 || fLeadingHadronType == 2)) {
 463     AliVCluster *cluster = jet->GetLeadingCluster(fClusterContainer->GetArray());
 464     if (cluster) {
 465       TLorentzVector nPart;
 466       cluster->GetMomentum(nPart, const_cast<Double_t*>(fVertex));
 467
 468       maxClusterEta = nPart.Eta();
 469       maxClusterPhi = nPart.Phi();
 470       maxClusterPt = nPart.Pt();
 471     }
 472   }
 473
 474   if (fParticleContainer && fParticleContainer->GetArray() && (fLeadingHadronType == 0 || fLeadingHadronType == 2)) {
 475     AliVParticle *track = jet->GetLeadingTrack(fParticleContainer->GetArray());
 476     if (track) {
 477       maxTrackEta = track->Eta();
 478       maxTrackPhi = track->Phi();
 479       maxTrackPt = track->Pt();
 480     }
 481   }
 482
 483   if (maxTrackPt > maxClusterPt)
 484     mom.SetPtEtaPhiM(maxTrackPt,maxTrackEta,maxTrackPhi,0.139);
 485   else
 486     mom.SetPtEtaPhiM(maxClusterPt,maxClusterEta,maxClusterPhi,0.139);
 487 }
 488
 489 //________________________________________________________________________
 490 Double_t AliJetContainer::GetZLeadingEmc(const AliEmcalJet *jet) const
 491 {
 492
 493   if (fClusterContainer && fClusterContainer->GetArray()) {
 494     TLorentzVector mom;
 495
 496     AliVCluster *cluster = jet->GetLeadingCluster(fClusterContainer->GetArray());
 497     if (cluster) {
 498       cluster->GetMomentum(mom, const_cast<Double_t*>(fVertex));
 499
 500       return GetZ(jet,mom);
 501     }
 502     else
 503       return -1;
 504   }
 505   else
 506     return -1;
 507 }
 508
 509 //________________________________________________________________________
 510 Double_t AliJetContainer::GetZLeadingCharged(const AliEmcalJet *jet) const
 511 {
 512
 513   if (fParticleContainer && fParticleContainer->GetArray() ) {
 514     TLorentzVector mom;
 515
 516     AliVParticle *track = jet->GetLeadingTrack(fParticleContainer->GetArray());
 517     if (track) {
 518       mom.SetPtEtaPhiM(track->Pt(),track->Eta(),track->Phi(),0.139);
 519
 520       return GetZ(jet,mom);
 521     }
 522     else
 523       return -1;
 524   }
 525   else
 526     return -1;
 527 }
 528
 529 //________________________________________________________________________
 530 Double_t AliJetContainer::GetZ(const AliEmcalJet *jet, TLorentzVector mom) const
 531 {
 532
 533   Double_t pJetSq = jet->Px()*jet->Px() + jet->Py()*jet->Py() + jet->Pz()*jet->Pz();
 534
 535   if(pJetSq>1e-6)
 536     return (mom.Px()*jet->Px() + mom.Py()*jet->Py() + mom.Pz()*jet->Pz())/pJetSq;
 537   else {
 538     AliWarning(Form("%s: strange, pjet*pjet seems to be zero pJetSq: %f",GetName(), pJetSq));
 539     return -1;
 540   }
 541
 542 }
 543
 544 //________________________________________________________________________
 545 void AliJetContainer::SetJetEtaPhiEMCAL()
 546 {
 547   //Set default cuts for full jets
 548
 549   if(!fGeom) SetEMCALGeometry();
 550   if(fGeom) {
 551     SetJetEtaLimits(fGeom->GetArm1EtaMin() + fJetRadius, fGeom->GetArm1EtaMax() - fJetRadius);
 552
 553     if(fRunNumber>=177295 && fRunNumber<=197470) //small SM masked in 2012 and 2013
 554       SetJetPhiLimits(1.405+fJetRadius,3.135-fJetRadius);
 555     else
 556       SetJetPhiLimits(fGeom->GetArm1PhiMin() * TMath::DegToRad() + fJetRadius, fGeom->GetArm1PhiMax() * TMath::DegToRad() - fJetRadius);
 557   }
 558   else {
 559     AliWarning("Could not get instance of AliEMCALGeometry. Using manual settings for EMCAL year 2011!!");
 560     SetJetEtaLimits(-0.7+fJetRadius,0.7-fJetRadius);
 561     SetJetPhiLimits(1.405+fJetRadius,3.135-fJetRadius);
 562   }
 563 }
 564
 565 //________________________________________________________________________
 566 void AliJetContainer::SetJetEtaPhiTPC()
 567 {
 568   //Set default cuts for charged jets
 569
 570   SetJetEtaLimits(-0.9+fJetRadius, 0.9-fJetRadius);
 571   SetJetPhiLimits(-10, 10);
 572 }
 573
 574 //________________________________________________________________________
 575 void AliJetContainer::PrintCuts()
 576 {
 577   // Reset cuts to default values
 578
 579   Printf("PtBiasJetTrack: %f",fPtBiasJetTrack);
 580   Printf("PtBiasJetClus: %f",fPtBiasJetClus);
 581   Printf("JetPtCut: %f", fJetPtCut);
 582   Printf("JetAreaCut: %f",fJetAreaCut);
 583   Printf("AreaEmcCut: %f",fAreaEmcCut);
 584   Printf("JetMinEta: %f", fJetMinEta);
 585   Printf("JetMaxEta: %f", fJetMaxEta);
 586   Printf("JetMinPhi: %f", fJetMinPhi);
 587   Printf("JetMaxPhi: %f", fJetMaxPhi);
 588   Printf("MaxClusterPt: %f",fMaxClusterPt);
 589   Printf("MaxTrackPt: %f",fMaxTrackPt);
 590   Printf("LeadingHadronType: %d",fLeadingHadronType);
 591   Printf("ZLeadingEmcCut: %f",fZLeadingEmcCut);
 592   Printf("ZLeadingChCut: %f",fZLeadingChCut);
 593
 594 }
 595
 596 //________________________________________________________________________
 597 void AliJetContainer::ResetCuts()
 598 {
 599   // Reset cuts to default values
 600
 601   fPtBiasJetTrack = 0;
 602   fPtBiasJetClus = 0;
 603   fJetPtCut = 1;
 604   fJetAreaCut = -1;
 605   fAreaEmcCut = 0;
 606   fJetMinEta = -0.9;
 607   fJetMaxEta = 0.9;
 608   fJetMinPhi = -10;
 609   fJetMaxPhi = 10;
 610   fMaxClusterPt = 1000;
 611   fMaxTrackPt = 100;
 612   fLeadingHadronType = 0;
 613   fZLeadingEmcCut = 10.;
 614   fZLeadingChCut = 10.;
 615 }
 616
 617 //________________________________________________________________________
 618 void AliJetContainer::SetClassName(const char *clname)
 619 {
 620   // Set the class name
 621
 622   TClass cls(clname);
 623   if (cls.InheritsFrom("AliEmcalJet")) fClassName = clname;
 624   else AliError(Form("Unable to set class name %s for a AliJetContainer, it must inherits from AliEmcalJet!",clname));
 625 }
 626
 627 //________________________________________________________________________
 628 Double_t AliJetContainer::GetFractionSharedPt(const AliEmcalJet *jet1) const
 629 {
 630   //
 631   // Get fraction of shared pT between matched jets
 632   // Uses ClosestJet() jet pT as baseline: fraction = \Sum_{const,jet1} pT,const,i / pT,jet,closest
 633   // Only works if tracks array of both jets is the same
 634   //
 635
 636   AliEmcalJet *jet2 = jet1->ClosestJet();
 637   if(!jet2) return -1;
 638
 639   Double_t fraction = 0.;
 640   Double_t jetPt2 = jet2->Pt();
 641
 642   if(jetPt2>0) {
 643     Double_t sumPt = 0.;
 644     AliVParticle *vpf = 0x0;
 645     Int_t iFound = 0;
 646     for(Int_t icc=0; icc<jet2->GetNumberOfTracks(); icc++) {
 647       Int_t idx = (Int_t)jet2->TrackAt(icc);
 648       iFound = 0;
 649       for(Int_t icf=0; icf<jet1->GetNumberOfTracks(); icf++) {
 650         if(idx == jet1->TrackAt(icf) && iFound==0 ) {
 651           iFound=1;
 652           vpf = static_cast<AliVParticle*>(jet1->TrackAt(icf, fParticleContainer->GetArray()));
 653           if(vpf) sumPt += vpf->Pt();
 654           continue;
 655         }
 656       }
 657     }
 658     fraction = sumPt/jetPt2;
 659   } else
 660     fraction = -1;
 661
 662   return fraction;
 663 }
 664