EMCAL/AliEMCALJetFinderInputSimPrep.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16
  17 /* $Id$ */
  18
  19 //_________________________________________________________________________
  20 //  Class for JetFinder Input preparation from simulated data
  21 //
  22 //*-- Author: Mark Horner (LBL/UCT)
  23 //
  24
  25
  26
  27 #include "AliEMCALJetFinderInputSimPrep.h"
  28
  29 #include <TParticle.h>
  30 #include <TParticlePDG.h>
  31 #include <TPDGCode.h>
  32 #include <TFile.h>
  33 #include <TTree.h>
  34 #include <TObjectTable.h>
  35 #include <TMath.h>
  36
  37
  38 #include "AliRun.h"
  39 #include "AliMagF.h"
  40 #include "AliEMCALFast.h"
  41 #include "AliEMCAL.h"
  42 #include "AliEMCALHit.h"
  43 #include "AliEMCALGeometry.h"
  44 #include "AliEMCALGetter.h"
  45 #include "AliGenEventHeader.h"
  46 #include "AliGenPythiaEventHeader.h"
  47 #include "AliGenerator.h"
  48 #include "AliHeader.h"
  49 #include "AliMC.h"
  50
  51
  52 ClassImp(AliEMCALJetFinderInputSimPrep)
  53
  54 AliEMCALJetFinderInputSimPrep::AliEMCALJetFinderInputSimPrep()
  55 {
  56 if (fDebug > 0) Info("AliEMCALJetFinderInputSimPrep","Beginning Constructor");
  57
  58   fDebug = 0;
  59   fSmearType = kSmearEffic;
  60   fEMCALType = kHits;
  61   fTrackType = kCharged;
  62   fEfficiency = 0.90;
  63
  64 }
  65 AliEMCALJetFinderInputSimPrep::~AliEMCALJetFinderInputSimPrep()
  66 {
  67 if (fDebug > 0) Info("~AliEMCALJetFinderInputSimPrep","Beginning Destructor");
  68
  69 }
  70
  71 void AliEMCALJetFinderInputSimPrep::Reset(AliEMCALJetFinderResetType_t resettype)
  72 {
  73 if (fDebug > 1) Info("Reset","Beginning Reset");
  74         switch (resettype){
  75
  76         case kResetData:
  77                 fInputObject.Reset(resettype);
  78                 break;
  79         case kResetTracks:
  80                 fInputObject.Reset(resettype);
  81                 break;
  82         case kResetDigits:
  83                 fInputObject.Reset(resettype);
  84                 break;
  85         case kResetParameters:
  86                 fSmearType = kSmearEffic;
  87                 fEMCALType = kHits;
  88                 fTrackType = kCharged;
  89                 break;
  90         case kResetAll:
  91                 fSmearType = kSmearEffic;
  92                 fEMCALType = kHits;
  93                 fTrackType = kCharged;
  94                 fInputObject.Reset(resettype);
  95                 break;
  96         case kResetPartons:
  97           Warning("FillFromFile", "kResetPartons not implemented") ;
  98           break;
  99         case kResetParticles:
 100           Warning("FillFromFile", "kResetParticles not implemented") ;
 101           break;
 102         case kResetJets:
 103           Warning("FillFromFile", "kResetJets not implemented") ;
 104           break;
 105         }// end switch
 106
 107 }
 108
 109 Int_t AliEMCALJetFinderInputSimPrep::FillFromFile(TString *filename, AliEMCALJetFinderFileType_t filetype,Int_t EventNumber)
 110 {
 111   // gObjectTable->Print();
 112         // Test that file exists - Getter doesn't like bogus filenames
 113 if (fDebug > 1) Info("FillFromFile","Beginning FillFromFile");
 114    fFileType = filetype;
 115    TFile file(filename->Data());
 116    if (!file.IsOpen()){
 117            Error("FillFromFile","Could not open file in FillFromFile");
 118            return -1;
 119    }
 120    file.Close();
 121    /*
 122    gAlice =  static_cast<AliRun *>(file.Get("gAlice"));
 123    AliEMCAL* pEMCAL = (AliEMCAL*) gAlice->GetModule("EMCAL");
 124    if (gAlice->GetEvent(EventNumber) < 0){
 125            Error("FillFromFile","Could not open event in FillFromFile");
 126            file.Close();
 127            return -1;
 128    }*/
 129    AliEMCALGetter *gime = AliEMCALGetter::Instance(filename->Data());
 130    gime->Event(EventNumber) ;
 131
 132    if ( fEMCALType == kHits ||
 133         fEMCALType == kTimeCut )  FillHits();
 134    if ( fTrackType != kNoTracks  )  FillTracks();
 135    if ( fFileType  != kData){
 136            FillPartons();
 137 //         FillParticles();
 138    }
 139    //gime->CloseFile();
 140 //   gObjectTable->Print();
 141    delete gime;
 142    return 0;
 143 }
 144
 145 void AliEMCALJetFinderInputSimPrep::FillHits()          // Fill from the hits to input object from simulation
 146 {
 147 if (fDebug > 1) Info("FillHits","Beginning FillHits");
 148
 149 // Access hit information
 150     AliEMCALHit *mHit;
 151     AliEMCAL* pEMCAL = (AliEMCAL*) gAlice->GetModule("EMCAL");
 152     Info("AliEMCALJetFinderInputSimPrep","Title of the geometry is %s",pEMCAL->GetTitle());
 153     AliEMCALGeometry* geom =  AliEMCALGetter::Instance()->EMCALGeometry();
 154 //    Float_t samplingF = geom->GetSampling();
 155     Float_t samplingF = 11.6;
 156     Info("AliEMCALJetFinderInputSimPrep","Sampling fraction is %f",samplingF);
 157     TTree *treeH = AliEMCALGetter::Instance()->TreeH();
 158     Int_t ntracks = (Int_t) treeH->GetEntries();
 159 //
 160 //   Loop over tracks
 161 //
 162     Int_t nbytes = 0;
 163     Double_t etH = 0.0;
 164
 165     for (Int_t track=0; track<ntracks;track++) {
 166         gAlice->ResetHits();
 167         nbytes += treeH->GetEvent(track);
 168 //
 169 //  Loop over hits
 170 //
 171         for(mHit=(AliEMCALHit*) pEMCAL->FirstHit(-1);
 172             mHit;
 173             mHit=(AliEMCALHit*) pEMCAL->NextHit())
 174         {
 175                 if (fEMCALType == kTimeCut &&
 176                     (mHit->GetTime()>fTimeCut) ) continue;
 177             Float_t x      =    mHit->X();         // x-pos of hit
 178             Float_t y      =    mHit->Y();         // y-pos
 179             Float_t z      =    mHit->Z();         // z-pos
 180             Float_t eloss  =    mHit->GetEnergy(); // deposited energy
 181
 182             Float_t r      =    TMath::Sqrt(x*x+y*y);
 183             Float_t theta  =    TMath::ATan2(r,z);
 184             Float_t eta    =   -TMath::Log(TMath::Tan(theta/2.));
 185             Float_t phi    =    TMath::ATan2(y,x);
 186             etH = samplingF*eloss*TMath::Sin(theta);
 187             if (fDebug > 10) Info("FillHits","Values of hit energy %i",Int_t(1e7*etH));
 188             if (geom->TowerIndexFromEtaPhi(eta,180.0/TMath::Pi()*phi) == -1)
 189             {
 190                     if (fDebug >5)
 191                     {
 192                             Error("FillHits","Hit not inside EMCAL!");
 193                             mHit->Dump();
 194                     }
 195                     continue;
 196             }
 197             fInputObject.AddEnergyToDigit(geom->TowerIndexFromEtaPhi(eta,180.0/TMath::Pi()*phi)-1,Int_t(1e7*etH));
 198         } // Hit Loop
 199     } // Track Loop
 200
 201
 202 }
 203 void AliEMCALJetFinderInputSimPrep::FillTracks()        // Fill from particles simulating a TPC to input object from simulation
 204 {
 205
 206     if (fDebug > 1) Info("FillTracks","Beginning FillTracks");
 207
 208     TParticlePDG* pdgP = 0;
 209     TParticle *MPart;
 210     Int_t npart = (gAlice->GetHeader())->GetNprimary();
 211     Float_t bfield,rEMCAL;
 212
 213     if (fDebug > 1) Info("FillTracks","Defining the geometry");
 214
 215     AliEMCAL* pEMCAL = (AliEMCAL*) gAlice->GetModule("EMCAL");
 216     AliEMCALGeometry* geom =  AliEMCALGeometry::GetInstance(pEMCAL->GetTitle(), "");
 217     fEtaMax = geom->GetArm1EtaMax();
 218     fEtaMin = geom->GetArm1EtaMin();
 219     fPhiMax = TMath::Pi()/180.0*geom->GetArm1PhiMax();
 220     fPhiMin = TMath::Pi()/180.0*geom->GetArm1PhiMin();
 221
 222     if (fDebug > 1) Info("FillTracks","Starting particle loop");
 223
 224     for (Int_t part = 0; part < npart; part++) {
 225         MPart = gAlice->GetMCApp()->Particle(part);
 226         //if (part%10) gObjectTable->Print();
 227         pdgP = MPart->GetPDG();
 228
 229         if (fDebug > 5) Info("FillTracks","Checking if track is a primary");
 230
 231         if (fFileType == kPythia) {
 232             if (MPart->GetStatusCode() != 1) continue;
 233         } else if (fFileType == kHijing) {
 234             if (MPart->GetFirstDaughter() >= 0 && MPart->GetFirstDaughter() < npart) continue;
 235         }
 236
 237         if (fDebug > 15) Info("FillTracks","Checking if track (eta - %f, phi - %f) is in acceptance",MPart->Eta(),MPart->Phi());
 238         if (fDebug > 10) Info("FillTracks","Checking if EMCAL acceptance  ( %f < eta < %f, %f < phi < %f) is in acceptance",fEtaMin,fEtaMax,fPhiMin,fPhiMax);
 239
 240         if (MPart->Eta() > fEtaMax || MPart->Eta() < fEtaMin)    continue;
 241         if (MPart->Phi() > fPhiMax || MPart->Phi() < fPhiMin )   continue;
 242
 243 /*
 244         {kProton, kProtonBar, kElectron, kPositron,
 245          kNuE, kNuEBar, kGamma, kNeutron, kNeutronBar,
 246          kMuonPlus, kMuonMinus, kK0Long , kPiPlus, kPiMinus,
 247          kKPlus, kKMinus, kLambda0, kLambda0Bar, kK0Short,
 248          kSigmaMinus, kSigmaPlus, kSigma0, kPi0, kK0, kK0Bar,
 249          0,kNuMu,kNuMuBar
 250 */
 251
 252         if (fDebug > 15) Info("FillTracks","Checking if track is rejected");
 253
 254         if ((fSmearType == kEfficiency  ||
 255              fSmearType == kSmearEffic) &&
 256              pdgP->Charge()!=0) {
 257             if (AliEMCALFast::RandomReject(fEfficiency)) {
 258                 continue;
 259             }
 260         }
 261
 262         bfield = gAlice->Field()->SolenoidField();
 263         rEMCAL = AliEMCALGeometry::GetInstance()->GetIPDistance();
 264         Float_t rB = 3335.6 * MPart->Pt() / bfield;  // [cm]  (case of |charge|=1)
 265         if (2.*rB < rEMCAL) continue;  // track curls away
 266
 267         //if (part%10) gObjectTable->Print();
 268         switch(fTrackType)
 269         {
 270
 271            case kAll:  // All Stable particles to be included
 272                 if (fDebug > 5) Info("FillTracks","Storing track");
 273                 if (fSmearType == kSmear ||
 274                     fSmearType == kSmearEffic ){
 275                         Smear(MPart);/*
 276                         TParticle *tmp = Smear(MPart);
 277                         fInputObject.AddTrack(Smear(MPart));
 278                         delete tmp;*/
 279                 }else{
 280                         fInputObject.AddTrack(*MPart);
 281                 }
 282            break;
 283            case kEM:   // All Electromagnetic particles
 284                 if (MPart->GetPdgCode() == kElectron  ||
 285                     MPart->GetPdgCode() == kMuonPlus  ||
 286                     MPart->GetPdgCode() == kMuonMinus ||
 287                     MPart->GetPdgCode() == kPositron ){
 288                       if (fDebug > 5) Info("FillTracks","Storing electron or positron");
 289                       if (fSmearType == kSmear ||
 290                             fSmearType == kSmearEffic ){
 291                               Smear(MPart);/*
 292                               TParticle *tmp = Smear(MPart);
 293                               fInputObject.AddTrack(tmp);
 294                               delete tmp;*/
 295                       }else{
 296                           fInputObject.AddTrack(*MPart);
 297                       }
 298                 }
 299                 if ( MPart->GetPdgCode() == kGamma ){
 300                         fInputObject.AddTrack(*MPart);
 301                         if (fDebug > 5) Info("FillTracks","Storing gamma");
 302                 }
 303
 304            break;
 305            case kCharged: // All Particles with non-zero charge
 306                 if (pdgP->Charge() != 0) {
 307                         if (fDebug > 5) Info("FillTracks","Storing charged track");
 308                         if (fSmearType == kSmear ||
 309                             fSmearType == kSmearEffic ){
 310                                 Smear(MPart);/*
 311                                 TParticle *tmp = Smear(MPart);
 312                                 fInputObject.AddTrack(tmp);
 313                                 delete tmp;*/
 314                         }else{
 315                          fInputObject.AddTrack(*MPart);
 316                         }
 317                 }
 318            break;
 319            case kNeutral: // All particles with no charge
 320                 if (pdgP->Charge() == 0){
 321                         fInputObject.AddTrack(*MPart);
 322                         if (fDebug > 5) Info("FillTracks","Storing neutral");
 323                 }
 324            break;
 325            case kHadron: //All hadrons
 326                 if (MPart->GetPdgCode() != kElectron  &&
 327                     MPart->GetPdgCode() != kPositron  &&
 328                     MPart->GetPdgCode() != kMuonPlus  &&
 329                     MPart->GetPdgCode() != kMuonMinus &&
 330                     MPart->GetPdgCode() != kGamma )
 331                 {
 332                         if (fDebug > 5) Info("FillTracks","Storing hadron");
 333                         if (pdgP->Charge() == 0){
 334                             fInputObject.AddTrack(*MPart);
 335                         }else{
 336                                 if (fSmearType == kSmear ||
 337                                     fSmearType == kSmearEffic ){
 338                                         Smear(MPart);/*
 339                                         TParticle *tmp = Smear(MPart);
 340                                         fInputObject.AddTrack(tmp);
 341                                         delete tmp;*/
 342                                 }else{
 343                                     fInputObject.AddTrack(*MPart);
 344                                 }
 345                         }
 346                 }
 347            break;
 348            case kChargedHadron:  // only charged hadrons
 349                 if (MPart->GetPdgCode() != kElectron  &&
 350                     MPart->GetPdgCode() != kPositron  &&
 351                     MPart->GetPdgCode() != kGamma     &&
 352                     MPart->GetPdgCode() != kMuonPlus  &&
 353                     MPart->GetPdgCode() != kMuonMinus &&
 354                     pdgP->Charge()      != 0       ){
 355                         if (fDebug > 5) Info("FillTracks","Storing charged hadron");
 356                         if (fSmearType == kSmear ||
 357                             fSmearType == kSmearEffic ){
 358                                 Smear(MPart);/*
 359                                 TParticle *tmp = Smear(MPart);
 360                                 fInputObject.AddTrack(tmp);
 361                                 delete tmp;*/
 362                         }else{
 363                                fInputObject.AddTrack(*MPart);
 364                         }
 365                 }
 366            break;
 367            case kNoTracks:
 368            break;
 369            default:
 370            break;
 371            delete MPart;
 372         }       //end of switch
 373 //      Info("FillTracks","After Particle Storage");
 374         //if (part%10) gObjectTable->Print();
 375
 376    }    //end of particle loop
 377    //gObjectTable->Print();
 378 }
 379
 380 void AliEMCALJetFinderInputSimPrep::FillPartons()               // Fill partons to input object from simulation
 381 {
 382 if (fDebug > 1) Info("FillPartons","Beginning FillPartons");
 383
 384   AliGenEventHeader* evHeader = ((AliHeader*)(gAlice->GetHeader()))->GenEventHeader();
 385   Float_t triggerJetValues[4];
 386   AliEMCALParton tempParton;
 387
 388   if (fFileType == kPythia)
 389   {
 390         Int_t ntriggerjets = ((AliGenPythiaEventHeader*)evHeader)->NTriggerJets();
 391         if (fDebug > 1) Info("FillPartons","Number of trigger jets --> %i",ntriggerjets);
 392         for (Int_t jetcounter = 0 ; jetcounter < ntriggerjets; jetcounter++){
 393         ((AliGenPythiaEventHeader*)evHeader)->TriggerJet(jetcounter,triggerJetValues);
 394         TLorentzVector tempLParton;
 395         tempLParton.SetPxPyPzE(triggerJetValues[0],triggerJetValues[1],triggerJetValues[2],triggerJetValues[3]);
 396         tempParton.SetEnergy(tempLParton.Energy());
 397         tempParton.SetEta(tempLParton.Eta());
 398         tempParton.SetPhi(tempLParton.Phi());
 399         FillPartonTracks(&tempParton);
 400         fInputObject.AddParton(&tempParton);
 401         }
 402
 403   }
 404 }
 405
 406 void AliEMCALJetFinderInputSimPrep::FillParticles()             // Fill particles to input object from simulation
 407 {
 408 if (fDebug > 1) Info("FillParticles","Beginning FillParticles");
 409
 410     Int_t npart = (gAlice->GetHeader())->GetNprimary();
 411     TParticlePDG* pdgP = 0;
 412
 413     for (Int_t part = 0; part < npart; part++) {
 414         TParticle *MPart = gAlice->GetMCApp()->Particle(part);
 415         pdgP = MPart->GetPDG();
 416
 417         if (fDebug > 10) Info("FillParticles","Checking if particle is a primary");
 418
 419         if (fFileType == kPythia) {
 420             if (MPart->GetStatusCode() != 1) continue;
 421         } else if (fFileType == kHijing) {
 422             if (MPart->GetFirstDaughter() >= 0 && MPart->GetFirstDaughter() < npart) continue;
 423         }
 424
 425         if (fDebug > 10) Info("FillParticles","Checking if particle is in acceptance");
 426
 427         if (MPart->Eta() > fEtaMax || MPart->Eta() < fEtaMin)    continue;
 428         if (MPart->Phi() > fPhiMax || MPart->Phi() < fPhiMin )   continue;
 429
 430
 431 /*
 432         {kProton, kProtonBar, kElectron, kPositron,
 433          kNuE, kNuEBar, kGamma, kNeutron, kNeutronBar,
 434          kMuonPlus, kMuonMinus, kK0Long , kPiPlus, kPiMinus,
 435          kKPlus, kKMinus, kLambda0, kLambda0Bar, kK0Short,
 436          kSigmaMinus, kSigmaPlus, kSigma0, kPi0, kK0, kK0Bar,
 437          0,kNuMu,kNuMuBar
 438 */
 439         switch(fTrackType)
 440         {
 441
 442            case kAll:  // All Stable particles to be included
 443                 if (fDebug > 5) Info("FillParticles","Storing particle");
 444                 fInputObject.AddParticle(MPart);
 445            break;
 446            case kEM:   // All Electromagnetic particles
 447                 if (MPart->GetPdgCode() == kElectron ||
 448                     MPart->GetPdgCode() == kPositron ||
 449                     MPart->GetPdgCode() == kGamma){
 450                         if (fDebug > 5) Info("FillParticles","Storing electromagnetic particle");
 451                         fInputObject.AddParticle(MPart);
 452                 }
 453            break;
 454            case kCharged: // All Particles with non-zero charge
 455                 if (pdgP->Charge() != 0) {
 456                         if (fDebug > 5) Info("FillParticles","Storing charged particle");
 457                         fInputObject.AddParticle(MPart);
 458                 }
 459            break;
 460            case kNeutral: // All particles with no charge
 461                 if (pdgP->Charge() == 0){
 462                         if (fDebug > 5) Info("FillParticles","Storing neutral particle");
 463                         fInputObject.AddParticle(MPart);
 464                 }
 465            break;
 466            case kHadron: //All hadrons
 467                 if (MPart->GetPdgCode() != kElectron &&
 468                     MPart->GetPdgCode() != kPositron &&
 469                     MPart->GetPdgCode() != kGamma )
 470                 {
 471
 472                 if (fDebug > 5) Info("FillParticles","Storing hadron");
 473                     fInputObject.AddParticle(MPart);
 474                 }
 475            break;
 476            case kChargedHadron:  // only charged hadrons
 477                 if (MPart->GetPdgCode() != kElectron &&
 478                     MPart->GetPdgCode() != kPositron &&
 479                     MPart->GetPdgCode() != kGamma    &&
 480                     pdgP->Charge()      != 0       ){
 481                 if (fDebug > 5) Info("FillParticles","Storing charged hadron");
 482                        fInputObject.AddParticle(MPart);
 483                 }
 484            break;
 485            case kNoTracks:
 486            break;
 487            default:
 488            break;
 489         }       //end of switch
 490     }// end of loop over particles
 491
 492 }
 493 void AliEMCALJetFinderInputSimPrep::FillDigits()                // Fill digits to input object
 494 {
 495
 496 }
 497
 498 void AliEMCALJetFinderInputSimPrep::Smear(TParticle *particle)
 499 {
 500 if (fDebug > 5) Info("Smear","Beginning Smear");
 501
 502         Float_t tmp = AliEMCALFast::SmearMomentum(1,particle->P());
 503         Float_t tmpt = (tmp/particle->P()) * particle->Pt();
 504         if (fDebug > 10) Info("Smear","Smearing changed P from %f to %f",particle->Pt(),tmpt);
 505         TLorentzVector tmpvector;
 506         tmpvector.SetPtEtaPhiM(tmpt,particle->Eta(), particle->Phi(), particle->GetMass());
 507 // create a new particle with smearing momentum - and no daughter or parent information and the same
 508 // vertex
 509
 510         TParticle tmparticle(particle->GetPdgCode(), 1, 0, 0,0,0,
 511                         tmpvector.Px()  , tmpvector.Py(), tmpvector.Pz(), tmpvector.Energy(),
 512                         particle->Vx(), particle->Vy(), particle->Vz(), particle->T());
 513         fInputObject.AddTrack(tmparticle);
 514         return;
 515 }
 516
 517 void AliEMCALJetFinderInputSimPrep::FillPartonTracks(AliEMCALParton *parton)
 518 {
 519
 520         if (fDebug>1) Info("FillPartonTracks","Beginning FillPartonTracks");
 521         Int_t npart = (gAlice->GetHeader())->GetNprimary();
 522         Int_t ntracks =0;
 523         TParticlePDG *getpdg;
 524         TParticle *tempPart;
 525         for (Int_t part = 0; part < npart; part++)
 526         {
 527                 tempPart = gAlice->GetMCApp()->Particle(part);
 528                 if (tempPart->GetStatusCode() != 1) continue;
 529                 if (tempPart->Eta() > fEtaMax || tempPart->Eta() < fEtaMin ||
 530                     tempPart->Phi() > fPhiMax || tempPart->Phi() < fPhiMin ){
 531                         if (fDebug>10) Info("FillPartonTracks","Excluding particle not pointing at the EMCAL");
 532                         continue;
 533                 }
 534                 getpdg = tempPart->GetPDG();
 535                 if (getpdg->Charge() == 0.0  ) {
 536                         if (fDebug>10) Info("FillPartonTracks","Excluding a neutral particle");
 537                         continue;
 538                 }
 539                 if ( (parton->Eta() - tempPart->Eta())*(parton->Eta() - tempPart->Eta()) +
 540                 (parton->Phi() - tempPart->Phi())*(parton->Phi() - tempPart->Phi()) < 1.0 ) ntracks++;
 541         }
 542         Float_t *Energy = new Float_t[ntracks];
 543         Float_t *Eta    = new Float_t[ntracks];
 544         Float_t *Phi    = new Float_t[ntracks];
 545         Int_t   *Pdg    = new Int_t[ntracks];
 546         ntracks=0;
 547         for (Int_t part = 0; part < npart; part++)
 548         {
 549                 tempPart = gAlice->GetMCApp()->Particle(part);
 550                 if (tempPart->GetStatusCode() != 1) continue;
 551                 if (tempPart->Eta() > fEtaMax || tempPart->Eta() < fEtaMin ||
 552                     tempPart->Phi() > fPhiMax || tempPart->Phi() < fPhiMin ){
 553                         if (fDebug>10) Info("FillPartonTracks","Excluding particle not pointing at the EMCAL");
 554                         continue;
 555                 }
 556                 if (tempPart->GetStatusCode() != 1) continue;
 557                 getpdg = tempPart->GetPDG();
 558                 if (getpdg->Charge() == 0.0  ) {
 559                         if (fDebug>10) Info("FillPartonTracks","Excluding a neutral particle");
 560                         continue;
 561                 }
 562                 if (   (parton->Eta() - tempPart->Eta())*(parton->Eta() - tempPart->Eta()) +
 563                 (parton->Phi() - tempPart->Phi())*(parton->Phi() - tempPart->Phi()) < 1.0 )
 564                 {
 565                         Energy[ntracks] = tempPart->Pt();
 566                         Eta[ntracks] = tempPart->Eta();
 567                         Phi[ntracks] = tempPart->Phi();
 568                         Pdg[ntracks] = tempPart->GetPdgCode();
 569                         ntracks++;
 570                 }
 571         }
 572         parton->SetTrackList(ntracks,Energy,Eta,Phi,Pdg);
 573 }
 574