]> git.uio.no Git - u/mrichter/AliRoot.git/blame - EMCAL/jetfinder/AliEMCALJetFinderInputSimPrep.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
The present commit corresponds to an important change in the way the
[u/mrichter/AliRoot.git] / EMCAL / jetfinder / AliEMCALJetFinderInputSimPrep.cxx
CommitLineData
45a58699 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 "AliRunLoader.h"
40#include "AliStack.h"
41#include "AliMagF.h"
42#include "AliEMCALFast.h"
43#include "AliEMCAL.h"
44#include "AliEMCALHit.h"
45#include "AliEMCALGeometry.h"
46#include "AliEMCALLoader.h"
47#include "AliGenEventHeader.h"
48#include "AliGenPythiaEventHeader.h"
49#include "AliGenerator.h"
50#include "AliHeader.h"
51#include "AliMC.h"
52
53
54
55ClassImp(AliEMCALJetFinderInputSimPrep)
56
18a21c7c 57AliEMCALJetFinderInputSimPrep::AliEMCALJetFinderInputSimPrep()
58 : fEMCALType(kHits),fSmearType(kSmearEffic),fTrackType(kCharged),fFileType(kPythia),fEfficiency(0.90),
59 fTimeCut(0.),fEtaMax(0.),fEtaMin(0.),fPhiMax(0.),fPhiMin(0.)
45a58699 60{
61 // Default constructor
62if (fDebug > 0) Info("AliEMCALJetFinderInputSimPrep","Beginning Constructor");
63
64 fDebug = 0;
65 fSmearType = kSmearEffic;
66 fEMCALType = kHits;
67 fTrackType = kCharged;
68 fEfficiency = 0.90;
69
70}
71AliEMCALJetFinderInputSimPrep::~AliEMCALJetFinderInputSimPrep()
72{
73if (fDebug > 0) Info("~AliEMCALJetFinderInputSimPrep","Beginning Destructor");
74
75}
76
77void AliEMCALJetFinderInputSimPrep::Reset(AliEMCALJetFinderResetType_t resettype)
78{
79 // Method to reset data
80 if (fDebug > 1) Info("Reset","Beginning Reset");
81 switch (resettype){
82
83 case kResetData:
84 fInputObject.Reset(resettype);
85 break;
86 case kResetTracks:
87 fInputObject.Reset(resettype);
88 break;
89 case kResetDigits:
90 fInputObject.Reset(resettype);
91 break;
92 case kResetParameters:
93 fSmearType = kSmearEffic;
94 fEMCALType = kHits;
95 fTrackType = kCharged;
96 break;
97 case kResetAll:
98 fSmearType = kSmearEffic;
99 fEMCALType = kHits;
100 fTrackType = kCharged;
101 fInputObject.Reset(resettype);
102 break;
103 case kResetPartons:
104 Warning("FillFromFile", "kResetPartons not implemented") ;
105 break;
106 case kResetParticles:
107 Warning("FillFromFile", "kResetParticles not implemented") ;
108 break;
109 case kResetJets:
110 Warning("FillFromFile", "kResetJets not implemented") ;
111 break;
112 }// end switch
113
114}
115
116Int_t AliEMCALJetFinderInputSimPrep::FillFromFile(TString *filename, AliEMCALJetFinderFileType_t filetype,Int_t EventNumber,TString data)
117{
118 if (fDebug > 1) Info("FillFromFile","Beginning FillFromFile");
119 fFileType = filetype;
120 AliRunLoader *rl=AliRunLoader::Open(filename->Data());
121 if (fDebug > 1) Info("FillFromFile","Opened file %s",filename->Data());
122 if (data.Contains("S"))
123 {
124 rl->LoadKinematics();
125 rl->LoadSDigits();
126 rl->GetEvent(EventNumber);
127 }else
128 {
129 rl->LoadKinematics();
130 rl->LoadHits();
131 rl->GetEvent(EventNumber);
132 }
133 if (fDebug > 1) Info("FillFromFile","Got event %i with option \"XH\"",EventNumber);
134
135 if ( fEMCALType == kHits ||
136 fEMCALType == kTimeCut )
137 {
138 if (data.Contains("S"))
139 {
140 FillSDigits();
141 }else
142 {
143 FillHits();
144 }
145 }
146 if ( fTrackType != kNoTracks ) FillTracks();
147 if ( fFileType != kData){
148 FillPartons();
149 }
150 return 0;
151}
152
153void AliEMCALJetFinderInputSimPrep::FillHits()
154{
155 // Fill from the hits to input object from simulation
156if (fDebug > 1) Info("FillHits","Beginning FillHits");
157
158// Access hit information
159// AliEMCAL* pEMCAL = (AliEMCAL*) gAlice->GetModule("EMCAL");
160// Info("AliEMCALJetFinderInputSimPrep","Title of the geometry is %s",pEMCAL->GetTitle());
dc7da436 161// AliEMCALGeometry* geom = AliEMCALGeometry::GetInstance();
45a58699 162 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(AliRunLoader::GetRunLoader()->GetDetectorLoader("EMCAL"));
163
164 // TTree *treeH = AliEMCALGetter::Instance()->TreeH();
165// Int_t ntracks = (Int_t) treeH->GetEntries();
166//
167// Loop over tracks
168//
169 Double_t etH = 0.0;
170
171/* for (Int_t track=0; track<ntracks;track++) {
172 gAlice->ResetHits();
173 nbytes += treeH->GetEvent(track);*/
174//
175// Loop over hits
176//
177 for(Int_t i =0; i < emcalLoader->Hits()->GetEntries() ;i++)
178 {
179 const AliEMCALHit *mHit = emcalLoader->Hit(i);
180 if (fEMCALType == kTimeCut &&
181 (mHit->GetTime()>fTimeCut) ) continue;
182 Float_t x = mHit->X(); // x-pos of hit
183 Float_t y = mHit->Y(); // y-pos
184 Float_t z = mHit->Z(); // z-pos
185 Float_t eloss = mHit->GetEnergy(); // deposited energy
186
187 Float_t r = TMath::Sqrt(x*x+y*y);
188 Float_t theta = TMath::ATan2(r,z);
dc7da436 189 //Float_t eta = -TMath::Log(TMath::Tan(theta/2.));
190 //Float_t phi = TMath::ATan2(y,x);
45a58699 191 etH = eloss*TMath::Sin(theta);
192 if (fDebug > 10) Info("FillHits","Values of hit energy %i",Int_t(1e7*etH));
dc7da436 193 /* have to be tune for TRD1; May 31,06
194 if (geom->TowerIndexFromEtaPhi(eta,180.0/TMath::Pi()*phi) == -1 )
45a58699 195 {
196 if (fDebug >5)
197 {
198 Error("FillHits","Hit not inside EMCAL!");
199 mHit->Dump();
200 }
201 continue;
202 }
203 fInputObject.AddEnergyToDigit(geom->TowerIndexFromEtaPhi(eta,180.0/TMath::Pi()*phi)-1,Int_t(1e7*etH));
dc7da436 204 */
45a58699 205 } // Hit Loop
206// } // Track Loop
207
208
209}
210void AliEMCALJetFinderInputSimPrep::FillTracks()
211{
212 // Fill from particles simulating a TPC to input object from simulation
213 if (fDebug > 1) Info("FillTracks","Beginning FillTracks");
214
215 AliRunLoader *rl = AliRunLoader::GetRunLoader();
216 TParticlePDG* pdgP = 0;
217 TParticle *mPart;
218 Int_t npart = rl->Stack()->GetNprimary();
219 if (fDebug > 1) Info("FillTracks","Header says there are %i primaries",npart);
220 Float_t bfield,rEMCAL;
221
222 if (fDebug > 1) Info("FillTracks","Defining the geometry");
223
224 AliEMCAL* pEMCAL = (AliEMCAL*) gAlice->GetModule("EMCAL");
225 AliEMCALGeometry* geom = AliEMCALGeometry::GetInstance();
226 if (geom == 0 && pEMCAL)
227 geom = AliEMCALGeometry::GetInstance(pEMCAL->GetTitle(), "");
228 if (geom == 0)
229 geom = AliEMCALGeometry::GetInstance("EMCAL_55_25","");//","SHISH_77_TRD1_2X2_FINAL_110DEG");
230 fEtaMax = geom->GetArm1EtaMax();
231 fEtaMin = geom->GetArm1EtaMin();
232 fPhiMax = TMath::Pi()/180.0*geom->GetArm1PhiMax();
233 fPhiMin = TMath::Pi()/180.0*geom->GetArm1PhiMin();
234
235 if (fDebug > 1) Info("FillTracks","Starting particle loop");
236
237 if (fDebug > 1) Info("FillTracks","Particle loop of %i",npart);
238 for (Int_t part = 0; part < npart; part++) {
239 mPart = rl->Stack()->Particle(part);
240 pdgP = mPart->GetPDG();
241
242 if (fDebug > 5) Info("FillTracks","Checking if track is a primary");
243
244 if (fFileType == kPythia) {
245 if (mPart->GetStatusCode() != 1) continue;
246 } else if (fFileType == kHijing) {
247 if (mPart->GetFirstDaughter() >= 0 && mPart->GetFirstDaughter() < npart) continue;
248 }
249
250 if (fDebug > 15) Info("FillTracks","Checking if track (eta - %f, phi - %f) is in acceptance",mPart->Eta(),mPart->Phi());
251 if (fDebug > 10) Info("FillTracks","Checking if EMCAL acceptance ( %f < eta < %f, %f < phi < %f) is in acceptance",fEtaMin,fEtaMax,fPhiMin,fPhiMax);
252
253 if ((!fPythiaComparison)&&(mPart->Eta() > fEtaMax || mPart->Eta() < fEtaMin)) continue;
254 if ((!fPythiaComparison)&&(mPart->Phi() > fPhiMax || mPart->Phi() < fPhiMin)) continue;
255
256/*
257 {kProton, kProtonBar, kElectron, kPositron,
258 kNuE, kNuEBar, kGamma, kNeutron, kNeutronBar,
259 kMuonPlus, kMuonMinus, kK0Long , kPiPlus, kPiMinus,
260 kKPlus, kKMinus, kLambda0, kLambda0Bar, kK0Short,
261 kSigmaMinus, kSigmaPlus, kSigma0, kPi0, kK0, kK0Bar,
262 0,kNuMu,kNuMuBar
263*/
264
265 if (fDebug > 15) Info("FillTracks","Checking if track is rejected");
266
267 if ((fSmearType == kEfficiency ||
268 fSmearType == kSmearEffic) &&
269 pdgP->Charge()!=0) {
270 if (AliEMCALFast::RandomReject(fEfficiency)) {
271 continue;
272 }
273 }
274
f7a1cc68 275 if (gAlice && ((AliMagF*)TGeoGlobalMagField::Instance()->GetField()))
276 bfield = ((AliMagF*)TGeoGlobalMagField::Instance()->GetField())->SolenoidField();
45a58699 277 else
278 bfield = 0.4;
279 rEMCAL = AliEMCALGeometry::GetInstance()->GetIPDistance();
280 Float_t rB = 3335.6 * mPart->Pt() / bfield; // [cm] (case of |charge|=1)
281 if (2.*rB < rEMCAL) continue; // track curls away
282
283 //if (part%10) gObjectTable->Print();
284 switch(fTrackType)
285 {
286
287 case kAllP: // All Stable particles to be included
288 if (fDebug > 5) Info("FillTracks","Storing track");
289 if (fSmearType == kSmear ||
290 fSmearType == kSmearEffic ){
291 Smear(mPart);/*
292 TParticle *tmp = Smear(MPart);
293 fInputObject.AddTrack(Smear(MPart));
294 delete tmp;*/
295 }else{
296 fInputObject.AddTrack(*mPart);
297 }
298 break;
299 case kEM: // All Electromagnetic particles
300 if (mPart->GetPdgCode() == kElectron ||
301 mPart->GetPdgCode() == kMuonPlus ||
302 mPart->GetPdgCode() == kMuonMinus ||
303 mPart->GetPdgCode() == kPositron ){
304 if (fDebug > 5) Info("FillTracks","Storing electron or positron");
305 if (fSmearType == kSmear ||
306 fSmearType == kSmearEffic ){
307 Smear(mPart);/*
308 TParticle *tmp = Smear(MPart);
309 fInputObject.AddTrack(tmp);
310 delete tmp;*/
311 }else{
312 fInputObject.AddTrack(*mPart);
313 }
314 }
315 if ( mPart->GetPdgCode() == kGamma ){
316 fInputObject.AddTrack(*mPart);
317 if (fDebug > 5) Info("FillTracks","Storing gamma");
318 }
319
320 break;
321 case kCharged: // All Particles with non-zero charge
322 if (pdgP->Charge() != 0) {
323 if (fDebug > 5) Info("FillTracks","Storing charged track");
324 if (fSmearType == kSmear ||
325 fSmearType == kSmearEffic ){
326 Smear(mPart);/*
327 TParticle *tmp = Smear(MPart);
328 fInputObject.AddTrack(tmp);
329 delete tmp;*/
330 }else{
331 fInputObject.AddTrack(*mPart);
332 }
333 }
334 break;
335 case kNeutral: // All particles with no charge
336 if (pdgP->Charge() == 0){
337 fInputObject.AddTrack(*mPart);
338 if (fDebug > 5) Info("FillTracks","Storing neutral");
339 }
340 break;
341 case kHadron: //All hadrons
342 if (mPart->GetPdgCode() != kElectron &&
343 mPart->GetPdgCode() != kPositron &&
344 mPart->GetPdgCode() != kMuonPlus &&
345 mPart->GetPdgCode() != kMuonMinus &&
346 mPart->GetPdgCode() != kGamma )
347 {
348 if (fDebug > 5) Info("FillTracks","Storing hadron");
349 if (pdgP->Charge() == 0){
350 fInputObject.AddTrack(*mPart);
351 }else{
352 if (fSmearType == kSmear ||
353 fSmearType == kSmearEffic ){
354 Smear(mPart);/*
355 TParticle *tmp = Smear(MPart);
356 fInputObject.AddTrack(tmp);
357 delete tmp;*/
358 }else{
359 fInputObject.AddTrack(*mPart);
360 }
361 }
362 }
363 break;
364 case kChargedHadron: // only charged hadrons
365 if (mPart->GetPdgCode() != kElectron &&
366 mPart->GetPdgCode() != kPositron &&
367 mPart->GetPdgCode() != kGamma &&
368 mPart->GetPdgCode() != kMuonPlus &&
369 mPart->GetPdgCode() != kMuonMinus &&
370 pdgP->Charge() != 0 ){
371 if (fDebug > 5) Info("FillTracks","Storing charged hadron");
372 if (fSmearType == kSmear ||
373 fSmearType == kSmearEffic ){
374 Smear(mPart);/*
375 TParticle *tmp = Smear(MPart);
376 fInputObject.AddTrack(tmp);
377 delete tmp;*/
378 }else{
379 fInputObject.AddTrack(*mPart);
380 }
381 }
382 break;
383 case kNoTracks:
384 break;
385 case kEMChargedPi0:
386 if (pdgP->Charge() != 0 || mPart->GetPdgCode() == kPi0 ||
387 mPart->GetPdgCode() == kGamma )
388 {
389 if (fDebug > 5) Info("FillTracks","Storing charged track");
390 if (fSmearType == kSmear ||
391 fSmearType == kSmearEffic ){
392 Smear(mPart);/*
393 TParticle *tmp = Smear(MPart);
394 fInputObject.AddTrack(tmp);
395 delete tmp;*/
396 }else{
397 fInputObject.AddTrack(*mPart);
398 }
399 }
400 break;
401 case kNoNeutronNeutrinoKlong:
402 if ( mPart->GetPdgCode() != kNeutron &&
403 mPart->GetPdgCode() != kNeutronBar &&
404 mPart->GetPdgCode() != kK0Long &&
405 mPart->GetPdgCode() != kNuE &&
406 mPart->GetPdgCode() != kNuEBar &&
407 mPart->GetPdgCode() != kNuMu &&
408 mPart->GetPdgCode() != kNuMuBar &&
409 mPart->GetPdgCode() != kNuTau &&
410 mPart->GetPdgCode() != kNuTauBar )
411 {
412 if (fDebug > 5) Info("FillTracks","Storing charged track");
413 if (fSmearType == kSmear ||
414 fSmearType == kSmearEffic ){
415 Smear(mPart);/*
416 TParticle *tmp = Smear(MPart);
417 fInputObject.AddTrack(tmp);
418 delete tmp;*/
419 }else{
420 fInputObject.AddTrack(*mPart);
421 }
422 }
423 break;
424 default:
425 break;
426 // delete mPart;
427 } //end of switch
428// Info("FillTracks","After Particle Storage");
429 //if (part%10) gObjectTable->Print();
430
431 } //end of particle loop
432 //gObjectTable->Print();
433}
434
435void AliEMCALJetFinderInputSimPrep::FillPartons()
436{
437 // Fill partons to
438 // input object from simulation
439if (fDebug > 1) Info("FillPartons","Beginning FillPartons");
440
441 AliRunLoader *rl = AliRunLoader::GetRunLoader();
442 AliGenEventHeader* evHeader = rl->GetHeader()->GenEventHeader();
443 Float_t triggerJetValues[4];
444 AliEMCALParton tempParton;
445
446 if (fFileType == kPythia)
447 {
448 Int_t ntriggerjets = ((AliGenPythiaEventHeader*)evHeader)->NTriggerJets();
449 if (fDebug > 1) Info("FillPartons","Number of trigger jets --> %i",ntriggerjets);
450 for (Int_t jetcounter = 0 ; jetcounter < ntriggerjets; jetcounter++){
451 ((AliGenPythiaEventHeader*)evHeader)->TriggerJet(jetcounter,triggerJetValues);
452 TLorentzVector tempLParton;
453 tempLParton.SetPxPyPzE(triggerJetValues[0],triggerJetValues[1],triggerJetValues[2],triggerJetValues[3]);
454 tempParton.SetEnergy(tempLParton.Energy());
455 tempParton.SetEta(tempLParton.Eta());
456 tempParton.SetPhi(tempLParton.Phi());
457 FillPartonTracks(&tempParton);
458 fInputObject.AddParton(&tempParton);
459 }
460
461 }
462}
463
464void AliEMCALJetFinderInputSimPrep::FillParticles()
465{
466 // Fill particles to input object from simulation
467if (fDebug > 1) Info("FillParticles","Beginning FillParticles");
468
469 Int_t npart = (gAlice->GetHeader())->GetNprimary();
470 TParticlePDG* pdgP = 0;
471
472 for (Int_t part = 0; part < npart; part++) {
473 TParticle *mPart = gAlice->GetMCApp()->Particle(part);
474 pdgP = mPart->GetPDG();
475
476 if (fDebug > 10) Info("FillParticles","Checking if particle is a primary");
477
478 if (fFileType == kPythia) {
479 if (mPart->GetStatusCode() != 1) continue;
480 } else if (fFileType == kHijing) {
481 if (mPart->GetFirstDaughter() >= 0 && mPart->GetFirstDaughter() < npart) continue;
482 }
483
484 if (fDebug > 10) Info("FillParticles","Checking if particle is in acceptance");
485
486 if ((!fPythiaComparison)&&(mPart->Eta() > fEtaMax || mPart->Eta() < fEtaMin)) continue;
487 if ((!fPythiaComparison)&&(mPart->Phi() > fPhiMax || mPart->Phi() < fPhiMin)) continue;
488
489
490/*
491 {kProton, kProtonBar, kElectron, kPositron,
492 kNuE, kNuEBar, kGamma, kNeutron, kNeutronBar,
493 kMuonPlus, kMuonMinus, kK0Long , kPiPlus, kPiMinus,
494 kKPlus, kKMinus, kLambda0, kLambda0Bar, kK0Short,
495 kSigmaMinus, kSigmaPlus, kSigma0, kPi0, kK0, kK0Bar,
496 0,kNuMu,kNuMuBar
497*/
498 switch(fTrackType)
499 {
500
501 case kAllP: // All Stable particles to be included
502 if (fDebug > 5) Info("FillParticles","Storing particle");
503 fInputObject.AddParticle(mPart);
504 break;
505 case kEM: // All Electromagnetic particles
506 if (mPart->GetPdgCode() == kElectron ||
507 mPart->GetPdgCode() == kPositron ||
508 mPart->GetPdgCode() == kGamma){
509 if (fDebug > 5) Info("FillParticles","Storing electromagnetic particle");
510 fInputObject.AddParticle(mPart);
511 }
512 break;
513 case kCharged: // All Particles with non-zero charge
514 if (pdgP->Charge() != 0) {
515 if (fDebug > 5) Info("FillParticles","Storing charged particle");
516 fInputObject.AddParticle(mPart);
517 }
518 break;
519 case kNeutral: // All particles with no charge
520 if (pdgP->Charge() == 0){
521 if (fDebug > 5) Info("FillParticles","Storing neutral particle");
522 fInputObject.AddParticle(mPart);
523 }
524 break;
525 case kHadron: //All hadrons
526 if (mPart->GetPdgCode() != kElectron &&
527 mPart->GetPdgCode() != kPositron &&
528 mPart->GetPdgCode() != kGamma )
529 {
530
531 if (fDebug > 5) Info("FillParticles","Storing hadron");
532 fInputObject.AddParticle(mPart);
533 }
534 break;
535 case kChargedHadron: // only charged hadrons
536 if (mPart->GetPdgCode() != kElectron &&
537 mPart->GetPdgCode() != kPositron &&
538 mPart->GetPdgCode() != kGamma &&
539 pdgP->Charge() != 0 ){
540 if (fDebug > 5) Info("FillParticles","Storing charged hadron");
541 fInputObject.AddParticle(mPart);
542 }
543 break;
544 case kNoTracks:
545 break;
546 case kEMChargedPi0:
547 if (pdgP->Charge() != 0 || mPart->GetPdgCode() == kPi0 ||
548 mPart->GetPdgCode() == kGamma )
549 {
550 if (fDebug > 5) Info("FillTracks","Storing kEMChargedPi0 track");
551 if (fSmearType == kSmear ||
552 fSmearType == kSmearEffic ){
553 Smear(mPart);/*
554 TParticle *tmp = Smear(MPart);
555 fInputObject.AddTrack(tmp);
556 delete tmp;*/
557 }else{
558 fInputObject.AddTrack(*mPart);
559 }
560 }
561 break;
562 case kNoNeutronNeutrinoKlong:
563 if ( mPart->GetPdgCode() != kNeutron &&
564 mPart->GetPdgCode() != kNeutronBar &&
565 mPart->GetPdgCode() != kK0Long &&
566 mPart->GetPdgCode() != kNuE &&
567 mPart->GetPdgCode() != kNuEBar &&
568 mPart->GetPdgCode() != kNuMu &&
569 mPart->GetPdgCode() != kNuMuBar &&
570 mPart->GetPdgCode() != kNuTau &&
571 mPart->GetPdgCode() != kNuTauBar )
572 {
573 if (fDebug > 5) Info("FillTracks","Storing kNoNeutronNeutrinoKlong track");
574 if (fSmearType == kSmear ||
575 fSmearType == kSmearEffic ){
576 Smear(mPart);/*
577 TParticle *tmp = Smear(MPart);
578 fInputObject.AddTrack(tmp);
579 delete tmp;*/
580 }else{
581 fInputObject.AddTrack(*mPart);
582 }
583 }
584 break;
585 default:
586 break;
587 } //end of switch
588 }// end of loop over particles
589
590}
591void AliEMCALJetFinderInputSimPrep::FillDigits()
592{
593 // Fill digits to input object
594
595}
596void AliEMCALJetFinderInputSimPrep::FillSDigits()
597{
598Info("FillSDigits","Beginning FillSDigits");
599 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(AliRunLoader::GetRunLoader()->GetDetectorLoader("EMCAL"));
600
601 // Fill digits to input object
602 for (Int_t towerid=0; towerid < emcalLoader->SDigits()->GetEntries(); towerid++)
603 {
604 fInputObject.AddEnergyToDigit(emcalLoader->SDigit(towerid)->GetId(), emcalLoader->SDigit(towerid)->GetAmp());
605 }
606
607}
608
609void AliEMCALJetFinderInputSimPrep::Smear(TParticle *particle)
610{
611 // Smear particle momentum
612if (fDebug > 5) Info("Smear","Beginning Smear");
613
614 Float_t tmp = AliEMCALFast::SmearMomentum(1,particle->P());
615 Float_t tmpt = (tmp/particle->P()) * particle->Pt();
616 if (fDebug > 10) Info("Smear","Smearing changed P from %f to %f",particle->Pt(),tmpt);
617 TLorentzVector tmpvector;
618 tmpvector.SetPtEtaPhiM(tmpt,particle->Eta(), particle->Phi(), particle->GetMass());
619// create a new particle with smearing momentum - and no daughter or parent information and the same
620// vertex
621
622 TParticle tmparticle(particle->GetPdgCode(), 1, 0, 0,0,0,
623 tmpvector.Px() , tmpvector.Py(), tmpvector.Pz(), tmpvector.Energy(),
624 particle->Vx(), particle->Vy(), particle->Vz(), particle->T());
625 fInputObject.AddTrack(tmparticle);
626 return;
627}
628
629void AliEMCALJetFinderInputSimPrep::FillPartonTracks(AliEMCALParton *parton)
630{
631 // Populate parton tracks for input distributions
632 if (fDebug>1) Info("FillPartonTracks","Beginning FillPartonTracks");
633 AliRunLoader *rl = AliRunLoader::GetRunLoader();
634
635 Int_t npart = rl->Stack()->GetNprimary();
636 Int_t ntracks =0;
637 TParticlePDG *getpdg;
638 TParticle *tempPart;
639 for (Int_t part = 0; part < npart; part++)
640 {
641 tempPart = rl->Stack()->Particle(part);
642 if (tempPart->GetStatusCode() != 1) continue;
643 if ((!fPythiaComparison)&&(tempPart->Eta() > fEtaMax || tempPart->Eta() < fEtaMin ||
644 tempPart->Phi() > fPhiMax || tempPart->Phi() < fPhiMin) ){
645 if (fDebug>10) Info("FillPartonTracks","Excluding particle not pointing at the EMCAL");
646 continue;
647 }
648 getpdg = tempPart->GetPDG();
649 if (getpdg->Charge() == 0.0 ) {
650 if (fDebug>10) Info("FillPartonTracks","Excluding a neutral particle");
651 continue;
652 }
653 if ( (parton->Eta() - tempPart->Eta())*(parton->Eta() - tempPart->Eta()) +
654 (parton->Phi() - tempPart->Phi())*(parton->Phi() - tempPart->Phi()) < 1.0 ) ntracks++;
655 }
656 Float_t *energy = new Float_t[ntracks];
657 Float_t *eta = new Float_t[ntracks];
658 Float_t *phi = new Float_t[ntracks];
659 Int_t *pdg = new Int_t[ntracks];
660 ntracks=0;
661 for (Int_t part = 0; part < npart; part++)
662 {
663 tempPart = rl->Stack()->Particle(part);
664 if (tempPart->GetStatusCode() != 1) continue;
665 if ((!fPythiaComparison)&&(tempPart->Eta() > fEtaMax || tempPart->Eta() < fEtaMin ||
666 tempPart->Phi() > fPhiMax || tempPart->Phi() < fPhiMin) ){
667 if (fDebug>10) Info("FillPartonTracks","Excluding particle not pointing at the EMCAL");
668 continue;
669 }
670 if (tempPart->GetStatusCode() != 1) continue;
671 getpdg = tempPart->GetPDG();
672 if (getpdg->Charge() == 0.0 ) {
673 if (fDebug>10) Info("FillPartonTracks","Excluding a neutral particle");
674 continue;
675 }
676 if ( (parton->Eta() - tempPart->Eta())*(parton->Eta() - tempPart->Eta()) +
677 (parton->Phi() - tempPart->Phi())*(parton->Phi() - tempPart->Phi()) < 1.0 )
678 {
679 energy[ntracks] = tempPart->Pt();
680 eta[ntracks] = tempPart->Eta();
681 phi[ntracks] = tempPart->Phi();
682 pdg[ntracks] = tempPart->GetPdgCode();
683 ntracks++;
684 }
685 }
686 parton->SetTrackList(ntracks,energy,eta,phi,pdg);
687 delete[] energy;
688 delete[] eta;
689 delete[] phi;
690 delete[] pdg;
691}
692