Enlarging array
[u/mrichter/AliRoot.git] / EMCAL / AliEMCALJetFinderInputSimPrep.cxx
CommitLineData
f7d5860b 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
ee6b678f 17/* $Id$ */
f7d5860b 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"
5d12ce38 49#include "AliMC.h"
f7d5860b 50
51
52ClassImp(AliEMCALJetFinderInputSimPrep)
53
54AliEMCALJetFinderInputSimPrep::AliEMCALJetFinderInputSimPrep()
55{
44f59d68 56 // Default constructor
f7d5860b 57if (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}
66AliEMCALJetFinderInputSimPrep::~AliEMCALJetFinderInputSimPrep()
67{
68if (fDebug > 0) Info("~AliEMCALJetFinderInputSimPrep","Beginning Destructor");
69
70}
71
72void AliEMCALJetFinderInputSimPrep::Reset(AliEMCALJetFinderResetType_t resettype)
73{
44f59d68 74 // Method to reset data
f7d5860b 75if (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
111Int_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
115if (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
44f59d68 147void AliEMCALJetFinderInputSimPrep::FillHits()
f7d5860b 148{
44f59d68 149 // Fill from the hits to input object from simulation
f7d5860b 150if (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();
f7d5860b 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);
5255a1ce 186 etH = eloss*TMath::Sin(theta);
f7d5860b 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}
44f59d68 203void AliEMCALJetFinderInputSimPrep::FillTracks()
f7d5860b 204{
44f59d68 205 // Fill from particles simulating a TPC to input object from simulation
f7d5860b 206
207 if (fDebug > 1) Info("FillTracks","Beginning FillTracks");
208
209 TParticlePDG* pdgP = 0;
44f59d68 210 TParticle *mPart;
f7d5860b 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++) {
44f59d68 226 mPart = gAlice->GetMCApp()->Particle(part);
f7d5860b 227 //if (part%10) gObjectTable->Print();
44f59d68 228 pdgP = mPart->GetPDG();
f7d5860b 229
230 if (fDebug > 5) Info("FillTracks","Checking if track is a primary");
231
232 if (fFileType == kPythia) {
44f59d68 233 if (mPart->GetStatusCode() != 1) continue;
f7d5860b 234 } else if (fFileType == kHijing) {
44f59d68 235 if (mPart->GetFirstDaughter() >= 0 && mPart->GetFirstDaughter() < npart) continue;
f7d5860b 236 }
237
44f59d68 238 if (fDebug > 15) Info("FillTracks","Checking if track (eta - %f, phi - %f) is in acceptance",mPart->Eta(),mPart->Phi());
f7d5860b 239 if (fDebug > 10) Info("FillTracks","Checking if EMCAL acceptance ( %f < eta < %f, %f < phi < %f) is in acceptance",fEtaMin,fEtaMax,fPhiMin,fPhiMax);
240
44f59d68 241 if (mPart->Eta() > fEtaMax || mPart->Eta() < fEtaMin) continue;
242 if (mPart->Phi() > fPhiMax || mPart->Phi() < fPhiMin ) continue;
f7d5860b 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();
44f59d68 265 Float_t rB = 3335.6 * mPart->Pt() / bfield; // [cm] (case of |charge|=1)
f7d5860b 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 ){
44f59d68 276 Smear(mPart);/*
f7d5860b 277 TParticle *tmp = Smear(MPart);
278 fInputObject.AddTrack(Smear(MPart));
279 delete tmp;*/
280 }else{
44f59d68 281 fInputObject.AddTrack(*mPart);
f7d5860b 282 }
283 break;
284 case kEM: // All Electromagnetic particles
44f59d68 285 if (mPart->GetPdgCode() == kElectron ||
286 mPart->GetPdgCode() == kMuonPlus ||
287 mPart->GetPdgCode() == kMuonMinus ||
288 mPart->GetPdgCode() == kPositron ){
f7d5860b 289 if (fDebug > 5) Info("FillTracks","Storing electron or positron");
290 if (fSmearType == kSmear ||
291 fSmearType == kSmearEffic ){
44f59d68 292 Smear(mPart);/*
f7d5860b 293 TParticle *tmp = Smear(MPart);
294 fInputObject.AddTrack(tmp);
295 delete tmp;*/
296 }else{
44f59d68 297 fInputObject.AddTrack(*mPart);
f7d5860b 298 }
299 }
44f59d68 300 if ( mPart->GetPdgCode() == kGamma ){
301 fInputObject.AddTrack(*mPart);
f7d5860b 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 ){
44f59d68 311 Smear(mPart);/*
f7d5860b 312 TParticle *tmp = Smear(MPart);
313 fInputObject.AddTrack(tmp);
314 delete tmp;*/
315 }else{
44f59d68 316 fInputObject.AddTrack(*mPart);
f7d5860b 317 }
318 }
319 break;
320 case kNeutral: // All particles with no charge
321 if (pdgP->Charge() == 0){
44f59d68 322 fInputObject.AddTrack(*mPart);
f7d5860b 323 if (fDebug > 5) Info("FillTracks","Storing neutral");
324 }
325 break;
326 case kHadron: //All hadrons
44f59d68 327 if (mPart->GetPdgCode() != kElectron &&
328 mPart->GetPdgCode() != kPositron &&
329 mPart->GetPdgCode() != kMuonPlus &&
330 mPart->GetPdgCode() != kMuonMinus &&
331 mPart->GetPdgCode() != kGamma )
f7d5860b 332 {
333 if (fDebug > 5) Info("FillTracks","Storing hadron");
334 if (pdgP->Charge() == 0){
44f59d68 335 fInputObject.AddTrack(*mPart);
f7d5860b 336 }else{
337 if (fSmearType == kSmear ||
338 fSmearType == kSmearEffic ){
44f59d68 339 Smear(mPart);/*
f7d5860b 340 TParticle *tmp = Smear(MPart);
341 fInputObject.AddTrack(tmp);
342 delete tmp;*/
343 }else{
44f59d68 344 fInputObject.AddTrack(*mPart);
f7d5860b 345 }
346 }
347 }
348 break;
349 case kChargedHadron: // only charged hadrons
44f59d68 350 if (mPart->GetPdgCode() != kElectron &&
351 mPart->GetPdgCode() != kPositron &&
352 mPart->GetPdgCode() != kGamma &&
353 mPart->GetPdgCode() != kMuonPlus &&
354 mPart->GetPdgCode() != kMuonMinus &&
f7d5860b 355 pdgP->Charge() != 0 ){
356 if (fDebug > 5) Info("FillTracks","Storing charged hadron");
357 if (fSmearType == kSmear ||
358 fSmearType == kSmearEffic ){
44f59d68 359 Smear(mPart);/*
f7d5860b 360 TParticle *tmp = Smear(MPart);
361 fInputObject.AddTrack(tmp);
362 delete tmp;*/
363 }else{
44f59d68 364 fInputObject.AddTrack(*mPart);
f7d5860b 365 }
366 }
367 break;
368 case kNoTracks:
369 break;
370 default:
371 break;
44f59d68 372 delete mPart;
f7d5860b 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
44f59d68 381void AliEMCALJetFinderInputSimPrep::FillPartons()
f7d5860b 382{
44f59d68 383 // Fill partons to
384 // input object from simulation
f7d5860b 385if (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
44f59d68 409void AliEMCALJetFinderInputSimPrep::FillParticles()
f7d5860b 410{
44f59d68 411 // Fill particles to input object from simulation
f7d5860b 412if (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++) {
44f59d68 418 TParticle *mPart = gAlice->GetMCApp()->Particle(part);
419 pdgP = mPart->GetPDG();
f7d5860b 420
421 if (fDebug > 10) Info("FillParticles","Checking if particle is a primary");
422
423 if (fFileType == kPythia) {
44f59d68 424 if (mPart->GetStatusCode() != 1) continue;
f7d5860b 425 } else if (fFileType == kHijing) {
44f59d68 426 if (mPart->GetFirstDaughter() >= 0 && mPart->GetFirstDaughter() < npart) continue;
f7d5860b 427 }
428
429 if (fDebug > 10) Info("FillParticles","Checking if particle is in acceptance");
430
44f59d68 431 if (mPart->Eta() > fEtaMax || mPart->Eta() < fEtaMin) continue;
432 if (mPart->Phi() > fPhiMax || mPart->Phi() < fPhiMin ) continue;
f7d5860b 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");
44f59d68 448 fInputObject.AddParticle(mPart);
f7d5860b 449 break;
450 case kEM: // All Electromagnetic particles
44f59d68 451 if (mPart->GetPdgCode() == kElectron ||
452 mPart->GetPdgCode() == kPositron ||
453 mPart->GetPdgCode() == kGamma){
f7d5860b 454 if (fDebug > 5) Info("FillParticles","Storing electromagnetic particle");
44f59d68 455 fInputObject.AddParticle(mPart);
f7d5860b 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");
44f59d68 461 fInputObject.AddParticle(mPart);
f7d5860b 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");
44f59d68 467 fInputObject.AddParticle(mPart);
f7d5860b 468 }
469 break;
470 case kHadron: //All hadrons
44f59d68 471 if (mPart->GetPdgCode() != kElectron &&
472 mPart->GetPdgCode() != kPositron &&
473 mPart->GetPdgCode() != kGamma )
f7d5860b 474 {
475
476 if (fDebug > 5) Info("FillParticles","Storing hadron");
44f59d68 477 fInputObject.AddParticle(mPart);
f7d5860b 478 }
479 break;
480 case kChargedHadron: // only charged hadrons
44f59d68 481 if (mPart->GetPdgCode() != kElectron &&
482 mPart->GetPdgCode() != kPositron &&
483 mPart->GetPdgCode() != kGamma &&
f7d5860b 484 pdgP->Charge() != 0 ){
485 if (fDebug > 5) Info("FillParticles","Storing charged hadron");
44f59d68 486 fInputObject.AddParticle(mPart);
f7d5860b 487 }
488 break;
489 case kNoTracks:
490 break;
491 default:
492 break;
493 } //end of switch
494 }// end of loop over particles
495
496}
44f59d68 497void AliEMCALJetFinderInputSimPrep::FillDigits()
f7d5860b 498{
44f59d68 499 // Fill digits to input object
f7d5860b 500
501}
502
503void AliEMCALJetFinderInputSimPrep::Smear(TParticle *particle)
504{
44f59d68 505 // Smear particle momentum
f7d5860b 506if (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
523void AliEMCALJetFinderInputSimPrep::FillPartonTracks(AliEMCALParton *parton)
524{
44f59d68 525 // Populate parton tracks for input distributions
f7d5860b 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 {
5d12ce38 533 tempPart = gAlice->GetMCApp()->Particle(part);
f7d5860b 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 }
44f59d68 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];
f7d5860b 552 ntracks=0;
553 for (Int_t part = 0; part < npart; part++)
554 {
5d12ce38 555 tempPart = gAlice->GetMCApp()->Particle(part);
f7d5860b 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 {
44f59d68 571 energy[ntracks] = tempPart->Pt();
572 eta[ntracks] = tempPart->Eta();
573 phi[ntracks] = tempPart->Phi();
574 pdg[ntracks] = tempPart->GetPdgCode();
f7d5860b 575 ntracks++;
576 }
577 }
44f59d68 578 parton->SetTrackList(ntracks,energy,eta,phi,pdg);
f7d5860b 579}
580