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         // Default constructor
  57 if (fDebug > 0) Info("AliEMCALJetFinderInputSimPrep","Beginning Constructor");
  58
  59   fDebug = 0;
  60   fSmearType = kSmearEffic;
  61   fEMCALType = kHits;
  62   fTrackType = kCharged;
  63   fEfficiency = 0.90;
  64
  65 }
  66 AliEMCALJetFinderInputSimPrep::~AliEMCALJetFinderInputSimPrep()
  67 {
  68 if (fDebug > 0) Info("~AliEMCALJetFinderInputSimPrep","Beginning Destructor");
  69
  70 }
  71
  72 void AliEMCALJetFinderInputSimPrep::Reset(AliEMCALJetFinderResetType_t resettype)
  73 {
  74         // Method to reset data
  75 if (fDebug > 1) Info("Reset","Beginning Reset");
  76         switch (resettype){
  77
  78         case kResetData:
  79                 fInputObject.Reset(resettype);
  80                 break;
  81         case kResetTracks:
  82                 fInputObject.Reset(resettype);
  83                 break;
  84         case kResetDigits:
  85                 fInputObject.Reset(resettype);
  86                 break;
  87         case kResetParameters:
  88                 fSmearType = kSmearEffic;
  89                 fEMCALType = kHits;
  90                 fTrackType = kCharged;
  91                 break;
  92         case kResetAll:
  93                 fSmearType = kSmearEffic;
  94                 fEMCALType = kHits;
  95                 fTrackType = kCharged;
  96                 fInputObject.Reset(resettype);
  97                 break;
  98         case kResetPartons:
  99           Warning("FillFromFile", "kResetPartons not implemented") ;
 100           break;
 101         case kResetParticles:
 102           Warning("FillFromFile", "kResetParticles not implemented") ;
 103           break;
 104         case kResetJets:
 105           Warning("FillFromFile", "kResetJets not implemented") ;
 106           break;
 107         }// end switch
 108
 109 }
 110
 111 Int_t AliEMCALJetFinderInputSimPrep::FillFromFile(TString *filename, AliEMCALJetFinderFileType_t filetype,Int_t EventNumber)
 112 {
 113   // gObjectTable->Print();
 114         // Test that file exists - Getter doesn't like bogus filenames
 115 if (fDebug > 1) Info("FillFromFile","Beginning FillFromFile");
 116    fFileType = filetype;
 117    TFile file(filename->Data());
 118    if (!file.IsOpen()){
 119            Error("FillFromFile","Could not open file in FillFromFile");
 120            return -1;
 121    }
 122    file.Close();
 123    /*
 124    gAlice =  static_cast<AliRun *>(file.Get("gAlice"));
 125    AliEMCAL* pEMCAL = (AliEMCAL*) gAlice->GetModule("EMCAL");
 126    if (gAlice->GetEvent(EventNumber) < 0){
 127            Error("FillFromFile","Could not open event in FillFromFile");
 128            file.Close();
 129            return -1;
 130    }*/
 131    AliEMCALGetter *gime = AliEMCALGetter::Instance(filename->Data());
 132    gime->Event(EventNumber) ;
 133
 134    if ( fEMCALType == kHits ||
 135         fEMCALType == kTimeCut )  FillHits();
 136    if ( fTrackType != kNoTracks  )  FillTracks();
 137    if ( fFileType  != kData){
 138            FillPartons();
 139 //         FillParticles();
 140    }
 141    //gime->CloseFile();
 142 //   gObjectTable->Print();
 143    delete gime;
 144    return 0;
 145 }
 146
 147 void AliEMCALJetFinderInputSimPrep::FillHits()
 148 {
 149         // Fill from the hits to input object from simulation
 150 if (fDebug > 1) Info("FillHits","Beginning FillHits");
 151
 152 // Access hit information
 153     AliEMCALHit *mHit;
 154     AliEMCAL* pEMCAL = (AliEMCAL*) gAlice->GetModule("EMCAL");
 155     Info("AliEMCALJetFinderInputSimPrep","Title of the geometry is %s",pEMCAL->GetTitle());
 156     AliEMCALGeometry* geom =  AliEMCALGetter::Instance()->EMCALGeometry();
 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 = 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()
 204 {
 205         // Fill from particles simulating a TPC to input object from simulation
 206
 207     if (fDebug > 1) Info("FillTracks","Beginning FillTracks");
 208
 209     TParticlePDG* pdgP = 0;
 210     TParticle *mPart;
 211     Int_t npart = (gAlice->GetHeader())->GetNprimary();
 212     Float_t bfield,rEMCAL;
 213
 214     if (fDebug > 1) Info("FillTracks","Defining the geometry");
 215
 216     AliEMCAL* pEMCAL = (AliEMCAL*) gAlice->GetModule("EMCAL");
 217     AliEMCALGeometry* geom =  AliEMCALGeometry::GetInstance(pEMCAL->GetTitle(), "");
 218     fEtaMax = geom->GetArm1EtaMax();
 219     fEtaMin = geom->GetArm1EtaMin();
 220     fPhiMax = TMath::Pi()/180.0*geom->GetArm1PhiMax();
 221     fPhiMin = TMath::Pi()/180.0*geom->GetArm1PhiMin();
 222
 223     if (fDebug > 1) Info("FillTracks","Starting particle loop");
 224
 225     for (Int_t part = 0; part < npart; part++) {
 226         mPart = gAlice->GetMCApp()->Particle(part);
 227         //if (part%10) gObjectTable->Print();
 228         pdgP = mPart->GetPDG();
 229
 230         if (fDebug > 5) Info("FillTracks","Checking if track is a primary");
 231
 232         if (fFileType == kPythia) {
 233             if (mPart->GetStatusCode() != 1) continue;
 234         } else if (fFileType == kHijing) {
 235             if (mPart->GetFirstDaughter() >= 0 && mPart->GetFirstDaughter() < npart) continue;
 236         }
 237
 238         if (fDebug > 15) Info("FillTracks","Checking if track (eta - %f, phi - %f) is in acceptance",mPart->Eta(),mPart->Phi());
 239         if (fDebug > 10) Info("FillTracks","Checking if EMCAL acceptance  ( %f < eta < %f, %f < phi < %f) is in acceptance",fEtaMin,fEtaMax,fPhiMin,fPhiMax);
 240
 241         if (mPart->Eta() > fEtaMax || mPart->Eta() < fEtaMin)    continue;
 242         if (mPart->Phi() > fPhiMax || mPart->Phi() < fPhiMin )   continue;
 243
 244 /*
 245         {kProton, kProtonBar, kElectron, kPositron,
 246          kNuE, kNuEBar, kGamma, kNeutron, kNeutronBar,
 247          kMuonPlus, kMuonMinus, kK0Long , kPiPlus, kPiMinus,
 248          kKPlus, kKMinus, kLambda0, kLambda0Bar, kK0Short,
 249          kSigmaMinus, kSigmaPlus, kSigma0, kPi0, kK0, kK0Bar,
 250          0,kNuMu,kNuMuBar
 251 */
 252
 253         if (fDebug > 15) Info("FillTracks","Checking if track is rejected");
 254
 255         if ((fSmearType == kEfficiency  ||
 256              fSmearType == kSmearEffic) &&
 257              pdgP->Charge()!=0) {
 258             if (AliEMCALFast::RandomReject(fEfficiency)) {
 259                 continue;
 260             }
 261         }
 262
 263         bfield = gAlice->Field()->SolenoidField();
 264         rEMCAL = AliEMCALGeometry::GetInstance()->GetIPDistance();
 265         Float_t rB = 3335.6 * mPart->Pt() / bfield;  // [cm]  (case of |charge|=1)
 266         if (2.*rB < rEMCAL) continue;  // track curls away
 267
 268         //if (part%10) gObjectTable->Print();
 269         switch(fTrackType)
 270         {
 271
 272            case kAll:  // All Stable particles to be included
 273                 if (fDebug > 5) Info("FillTracks","Storing track");
 274                 if (fSmearType == kSmear ||
 275                     fSmearType == kSmearEffic ){
 276                         Smear(mPart);/*
 277                         TParticle *tmp = Smear(MPart);
 278                         fInputObject.AddTrack(Smear(MPart));
 279                         delete tmp;*/
 280                 }else{
 281                         fInputObject.AddTrack(*mPart);
 282                 }
 283            break;
 284            case kEM:   // All Electromagnetic particles
 285                 if (mPart->GetPdgCode() == kElectron  ||
 286                     mPart->GetPdgCode() == kMuonPlus  ||
 287                     mPart->GetPdgCode() == kMuonMinus ||
 288                     mPart->GetPdgCode() == kPositron ){
 289                       if (fDebug > 5) Info("FillTracks","Storing electron or positron");
 290                       if (fSmearType == kSmear ||
 291                             fSmearType == kSmearEffic ){
 292                               Smear(mPart);/*
 293                               TParticle *tmp = Smear(MPart);
 294                               fInputObject.AddTrack(tmp);
 295                               delete tmp;*/
 296                       }else{
 297                           fInputObject.AddTrack(*mPart);
 298                       }
 299                 }
 300                 if ( mPart->GetPdgCode() == kGamma ){
 301                         fInputObject.AddTrack(*mPart);
 302                         if (fDebug > 5) Info("FillTracks","Storing gamma");
 303                 }
 304
 305            break;
 306            case kCharged: // All Particles with non-zero charge
 307                 if (pdgP->Charge() != 0) {
 308                         if (fDebug > 5) Info("FillTracks","Storing charged track");
 309                         if (fSmearType == kSmear ||
 310                             fSmearType == kSmearEffic ){
 311                                 Smear(mPart);/*
 312                                 TParticle *tmp = Smear(MPart);
 313                                 fInputObject.AddTrack(tmp);
 314                                 delete tmp;*/
 315                         }else{
 316                          fInputObject.AddTrack(*mPart);
 317                         }
 318                 }
 319            break;
 320            case kNeutral: // All particles with no charge
 321                 if (pdgP->Charge() == 0){
 322                         fInputObject.AddTrack(*mPart);
 323                         if (fDebug > 5) Info("FillTracks","Storing neutral");
 324                 }
 325            break;
 326            case kHadron: //All hadrons
 327                 if (mPart->GetPdgCode() != kElectron  &&
 328                     mPart->GetPdgCode() != kPositron  &&
 329                     mPart->GetPdgCode() != kMuonPlus  &&
 330                     mPart->GetPdgCode() != kMuonMinus &&
 331                     mPart->GetPdgCode() != kGamma )
 332                 {
 333                         if (fDebug > 5) Info("FillTracks","Storing hadron");
 334                         if (pdgP->Charge() == 0){
 335                             fInputObject.AddTrack(*mPart);
 336                         }else{
 337                                 if (fSmearType == kSmear ||
 338                                     fSmearType == kSmearEffic ){
 339                                         Smear(mPart);/*
 340                                         TParticle *tmp = Smear(MPart);
 341                                         fInputObject.AddTrack(tmp);
 342                                         delete tmp;*/
 343                                 }else{
 344                                     fInputObject.AddTrack(*mPart);
 345                                 }
 346                         }
 347                 }
 348            break;
 349            case kChargedHadron:  // only charged hadrons
 350                 if (mPart->GetPdgCode() != kElectron  &&
 351                     mPart->GetPdgCode() != kPositron  &&
 352                     mPart->GetPdgCode() != kGamma     &&
 353                     mPart->GetPdgCode() != kMuonPlus  &&
 354                     mPart->GetPdgCode() != kMuonMinus &&
 355                     pdgP->Charge()      != 0       ){
 356                         if (fDebug > 5) Info("FillTracks","Storing charged hadron");
 357                         if (fSmearType == kSmear ||
 358                             fSmearType == kSmearEffic ){
 359                                 Smear(mPart);/*
 360                                 TParticle *tmp = Smear(MPart);
 361                                 fInputObject.AddTrack(tmp);
 362                                 delete tmp;*/
 363                         }else{
 364                                fInputObject.AddTrack(*mPart);
 365                         }
 366                 }
 367            break;
 368            case kNoTracks:
 369            break;
 370            default:
 371            break;
 372            delete mPart;
 373         }       //end of switch
 374 //      Info("FillTracks","After Particle Storage");
 375         //if (part%10) gObjectTable->Print();
 376
 377    }    //end of particle loop
 378    //gObjectTable->Print();
 379 }
 380
 381 void AliEMCALJetFinderInputSimPrep::FillPartons()
 382 {
 383         // Fill partons to
 384         // input object from simulation
 385 if (fDebug > 1) Info("FillPartons","Beginning FillPartons");
 386
 387   AliGenEventHeader* evHeader = ((AliHeader*)(gAlice->GetHeader()))->GenEventHeader();
 388   Float_t triggerJetValues[4];
 389   AliEMCALParton tempParton;
 390
 391   if (fFileType == kPythia)
 392   {
 393         Int_t ntriggerjets = ((AliGenPythiaEventHeader*)evHeader)->NTriggerJets();
 394         if (fDebug > 1) Info("FillPartons","Number of trigger jets --> %i",ntriggerjets);
 395         for (Int_t jetcounter = 0 ; jetcounter < ntriggerjets; jetcounter++){
 396         ((AliGenPythiaEventHeader*)evHeader)->TriggerJet(jetcounter,triggerJetValues);
 397         TLorentzVector tempLParton;
 398         tempLParton.SetPxPyPzE(triggerJetValues[0],triggerJetValues[1],triggerJetValues[2],triggerJetValues[3]);
 399         tempParton.SetEnergy(tempLParton.Energy());
 400         tempParton.SetEta(tempLParton.Eta());
 401         tempParton.SetPhi(tempLParton.Phi());
 402         FillPartonTracks(&tempParton);
 403         fInputObject.AddParton(&tempParton);
 404         }
 405
 406   }
 407 }
 408
 409 void AliEMCALJetFinderInputSimPrep::FillParticles()
 410 {
 411         // Fill particles to input object from simulation
 412 if (fDebug > 1) Info("FillParticles","Beginning FillParticles");
 413
 414     Int_t npart = (gAlice->GetHeader())->GetNprimary();
 415     TParticlePDG* pdgP = 0;
 416
 417     for (Int_t part = 0; part < npart; part++) {
 418         TParticle *mPart = gAlice->GetMCApp()->Particle(part);
 419         pdgP = mPart->GetPDG();
 420
 421         if (fDebug > 10) Info("FillParticles","Checking if particle is a primary");
 422
 423         if (fFileType == kPythia) {
 424             if (mPart->GetStatusCode() != 1) continue;
 425         } else if (fFileType == kHijing) {
 426             if (mPart->GetFirstDaughter() >= 0 && mPart->GetFirstDaughter() < npart) continue;
 427         }
 428
 429         if (fDebug > 10) Info("FillParticles","Checking if particle is in acceptance");
 430
 431         if (mPart->Eta() > fEtaMax || mPart->Eta() < fEtaMin)    continue;
 432         if (mPart->Phi() > fPhiMax || mPart->Phi() < fPhiMin )   continue;
 433
 434
 435 /*
 436         {kProton, kProtonBar, kElectron, kPositron,
 437          kNuE, kNuEBar, kGamma, kNeutron, kNeutronBar,
 438          kMuonPlus, kMuonMinus, kK0Long , kPiPlus, kPiMinus,
 439          kKPlus, kKMinus, kLambda0, kLambda0Bar, kK0Short,
 440          kSigmaMinus, kSigmaPlus, kSigma0, kPi0, kK0, kK0Bar,
 441          0,kNuMu,kNuMuBar
 442 */
 443         switch(fTrackType)
 444         {
 445
 446            case kAll:  // All Stable particles to be included
 447                 if (fDebug > 5) Info("FillParticles","Storing particle");
 448                 fInputObject.AddParticle(mPart);
 449            break;
 450            case kEM:   // All Electromagnetic particles
 451                 if (mPart->GetPdgCode() == kElectron ||
 452                     mPart->GetPdgCode() == kPositron ||
 453                     mPart->GetPdgCode() == kGamma){
 454                         if (fDebug > 5) Info("FillParticles","Storing electromagnetic particle");
 455                         fInputObject.AddParticle(mPart);
 456                 }
 457            break;
 458            case kCharged: // All Particles with non-zero charge
 459                 if (pdgP->Charge() != 0) {
 460                         if (fDebug > 5) Info("FillParticles","Storing charged particle");
 461                         fInputObject.AddParticle(mPart);
 462                 }
 463            break;
 464            case kNeutral: // All particles with no charge
 465                 if (pdgP->Charge() == 0){
 466                         if (fDebug > 5) Info("FillParticles","Storing neutral particle");
 467                         fInputObject.AddParticle(mPart);
 468                 }
 469            break;
 470            case kHadron: //All hadrons
 471                 if (mPart->GetPdgCode() != kElectron &&
 472                     mPart->GetPdgCode() != kPositron &&
 473                     mPart->GetPdgCode() != kGamma )
 474                 {
 475
 476                 if (fDebug > 5) Info("FillParticles","Storing hadron");
 477                     fInputObject.AddParticle(mPart);
 478                 }
 479            break;
 480            case kChargedHadron:  // only charged hadrons
 481                 if (mPart->GetPdgCode() != kElectron &&
 482                     mPart->GetPdgCode() != kPositron &&
 483                     mPart->GetPdgCode() != kGamma    &&
 484                     pdgP->Charge()      != 0       ){
 485                 if (fDebug > 5) Info("FillParticles","Storing charged hadron");
 486                        fInputObject.AddParticle(mPart);
 487                 }
 488            break;
 489            case kNoTracks:
 490            break;
 491            default:
 492            break;
 493         }       //end of switch
 494     }// end of loop over particles
 495
 496 }
 497 void AliEMCALJetFinderInputSimPrep::FillDigits()
 498 {
 499         // Fill digits to input object
 500
 501 }
 502
 503 void AliEMCALJetFinderInputSimPrep::Smear(TParticle *particle)
 504 {
 505         // Smear particle momentum
 506 if (fDebug > 5) Info("Smear","Beginning Smear");
 507
 508         Float_t tmp = AliEMCALFast::SmearMomentum(1,particle->P());
 509         Float_t tmpt = (tmp/particle->P()) * particle->Pt();
 510         if (fDebug > 10) Info("Smear","Smearing changed P from %f to %f",particle->Pt(),tmpt);
 511         TLorentzVector tmpvector;
 512         tmpvector.SetPtEtaPhiM(tmpt,particle->Eta(), particle->Phi(), particle->GetMass());
 513 // create a new particle with smearing momentum - and no daughter or parent information and the same
 514 // vertex
 515
 516         TParticle tmparticle(particle->GetPdgCode(), 1, 0, 0,0,0,
 517                         tmpvector.Px()  , tmpvector.Py(), tmpvector.Pz(), tmpvector.Energy(),
 518                         particle->Vx(), particle->Vy(), particle->Vz(), particle->T());
 519         fInputObject.AddTrack(tmparticle);
 520         return;
 521 }
 522
 523 void AliEMCALJetFinderInputSimPrep::FillPartonTracks(AliEMCALParton *parton)
 524 {
 525         // Populate parton tracks for input distributions
 526         if (fDebug>1) Info("FillPartonTracks","Beginning FillPartonTracks");
 527         Int_t npart = (gAlice->GetHeader())->GetNprimary();
 528         Int_t ntracks =0;
 529         TParticlePDG *getpdg;
 530         TParticle *tempPart;
 531         for (Int_t part = 0; part < npart; part++)
 532         {
 533                 tempPart = gAlice->GetMCApp()->Particle(part);
 534                 if (tempPart->GetStatusCode() != 1) continue;
 535                 if (tempPart->Eta() > fEtaMax || tempPart->Eta() < fEtaMin ||
 536                     tempPart->Phi() > fPhiMax || tempPart->Phi() < fPhiMin ){
 537                         if (fDebug>10) Info("FillPartonTracks","Excluding particle not pointing at the EMCAL");
 538                         continue;
 539                 }
 540                 getpdg = tempPart->GetPDG();
 541                 if (getpdg->Charge() == 0.0  ) {
 542                         if (fDebug>10) Info("FillPartonTracks","Excluding a neutral particle");
 543                         continue;
 544                 }
 545                 if ( (parton->Eta() - tempPart->Eta())*(parton->Eta() - tempPart->Eta()) +
 546                 (parton->Phi() - tempPart->Phi())*(parton->Phi() - tempPart->Phi()) < 1.0 ) ntracks++;
 547         }
 548         Float_t *energy = new Float_t[ntracks];
 549         Float_t *eta    = new Float_t[ntracks];
 550         Float_t *phi    = new Float_t[ntracks];
 551         Int_t   *pdg    = new Int_t[ntracks];
 552         ntracks=0;
 553         for (Int_t part = 0; part < npart; part++)
 554         {
 555                 tempPart = gAlice->GetMCApp()->Particle(part);
 556                 if (tempPart->GetStatusCode() != 1) continue;
 557                 if (tempPart->Eta() > fEtaMax || tempPart->Eta() < fEtaMin ||
 558                     tempPart->Phi() > fPhiMax || tempPart->Phi() < fPhiMin ){
 559                         if (fDebug>10) Info("FillPartonTracks","Excluding particle not pointing at the EMCAL");
 560                         continue;
 561                 }
 562                 if (tempPart->GetStatusCode() != 1) continue;
 563                 getpdg = tempPart->GetPDG();
 564                 if (getpdg->Charge() == 0.0  ) {
 565                         if (fDebug>10) Info("FillPartonTracks","Excluding a neutral particle");
 566                         continue;
 567                 }
 568                 if (   (parton->Eta() - tempPart->Eta())*(parton->Eta() - tempPart->Eta()) +
 569                 (parton->Phi() - tempPart->Phi())*(parton->Phi() - tempPart->Phi()) < 1.0 )
 570                 {
 571                         energy[ntracks] = tempPart->Pt();
 572                         eta[ntracks] = tempPart->Eta();
 573                         phi[ntracks] = tempPart->Phi();
 574                         pdg[ntracks] = tempPart->GetPdgCode();
 575                         ntracks++;
 576                 }
 577         }
 578         parton->SetTrackList(ntracks,energy,eta,phi,pdg);
 579 }
 580