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