o add options for other gases and Epsilons
[u/mrichter/AliRoot.git] / TPC / Upgrade / AliToyMCEventGeneratorSimple.cxx
CommitLineData
526ddf0e 1#include <iostream>
a1a695e5 2
526ddf0e 3#include <TDatabasePDG.h>
4#include <TRandom.h>
5#include <TF1.h>
32438f4e 6#include <TStopwatch.h>
f356075c 7#include <TFile.h>
8#include <TTree.h>
609d6d39 9#include <TRandom3.h>
a1a695e5 10
0403120d 11#include <AliESDtrackCuts.h>
12#include <AliESDtrack.h>
13#include <AliTPCLaserTrack.h>
14
de0014b7 15#include "AliToyMCEvent.h"
526ddf0e 16
a1a695e5 17#include "AliToyMCEventGeneratorSimple.h"
18
cdf56507 19/*
20
21
22*/
f356075c 23
de0014b7 24ClassImp(AliToyMCEventGeneratorSimple);
526ddf0e 25
26
de0014b7 27AliToyMCEventGeneratorSimple::AliToyMCEventGeneratorSimple()
28 :AliToyMCEventGenerator()
526ddf0e 29 ,fVertexMean(0.)
30 ,fVertexSigma(7.)
32438f4e 31 ,fNtracks(20)
f356075c 32 ,fInputFileNameESD("")
33 ,fESDCuts(0x0)
29b7f480 34 ,fInputIndex(0)
35 ,fESDEvent(0x0)
36 ,fESDTree(0x0)
37 ,fInputFile(0x0)
1e62e876 38 ,fHPt(0x0)
39 ,fHEta(0x0)
40 ,fHMult(0x0)
41 ,fHistosSet(kFALSE)
42 ,fParamFile(0x0)
526ddf0e 43{
44 //
45
46}
47//________________________________________________________________
de0014b7 48AliToyMCEventGeneratorSimple::AliToyMCEventGeneratorSimple(const AliToyMCEventGeneratorSimple &gen)
49 :AliToyMCEventGenerator(gen)
526ddf0e 50 ,fVertexMean(gen.fVertexMean)
51 ,fVertexSigma(gen.fVertexSigma)
32438f4e 52 ,fNtracks(gen.fNtracks)
f356075c 53 ,fInputFileNameESD(gen.fInputFileNameESD)
54 ,fESDCuts(gen.fESDCuts)
29b7f480 55 ,fInputIndex(gen.fInputIndex)
56 ,fESDEvent(gen.fESDEvent)
57 ,fESDTree(gen.fESDTree)
58 ,fInputFile(gen.fInputFile)
1e62e876 59 ,fHPt(gen.fHPt)
60 ,fHEta(gen.fHEta)
61 ,fHMult(gen.fHMult)
62 ,fHistosSet(gen.fHistosSet)
63 ,fParamFile(0x0)
526ddf0e 64{
65}
66//________________________________________________________________
de0014b7 67AliToyMCEventGeneratorSimple::~AliToyMCEventGeneratorSimple()
526ddf0e 68{
69}
70 //________________________________________________________________
de0014b7 71AliToyMCEventGeneratorSimple& AliToyMCEventGeneratorSimple::operator = (const AliToyMCEventGeneratorSimple &gen)
526ddf0e 72{
73 //assignment operator
74 if (&gen == this) return *this;
de0014b7 75 new (this) AliToyMCEventGeneratorSimple(gen);
526ddf0e 76
77 return *this;
78}
79//________________________________________________________________
1f195c60 80void AliToyMCEventGeneratorSimple::SetParametersToyGen(const Char_t* parfilename/*="files/params.root*/, Double_t vertexMean/*=0*/, Double_t vertexSigma/*=7.*/) {
526ddf0e 81 fVertexMean = vertexMean;
82 fVertexSigma = vertexSigma;
1e62e876 83 fParamFile = new TFile(parfilename, "read");
84 fHPt = (TH1F*) fParamFile->Get("hPt");
85 fHEta = (TH1F*) fParamFile->Get("hEta");
86 fHMult = (TH1I*) fParamFile->Get("hMult") ;
87 fHistosSet = kTRUE;
88
89
526ddf0e 90}
91//________________________________________________________________
cd8ed0ac 92AliToyMCEvent* AliToyMCEventGeneratorSimple::Generate(Double_t time)
93{
94 //
95 //
96 //
97
de0014b7 98 AliToyMCEvent *retEvent = new AliToyMCEvent();
526ddf0e 99 retEvent->SetT0(time);
100 retEvent->SetX(0);
101 retEvent->SetX(0);
e83fd282 102 Double_t vertexZ=0;
103 //limit vertex to +- 10cm
104 while ( TMath::Abs(vertexZ=gRandom->Gaus(fVertexMean,fVertexSigma))>10. );
105 retEvent->SetZ(vertexZ);
526ddf0e 106
e83fd282 107 Double_t etaCuts=0.9;
526ddf0e 108 Double_t mass = TDatabasePDG::Instance()->GetParticle("pi+")->Mass();
e83fd282 109 static TF1 fpt("fpt",Form("x*(1+(sqrt(x*x+%f^2)-%f)/([0]*[1]))^(-[0])",mass,mass),0.3,100);
0403120d 110 if (fpt.GetParameter(0)<1){
111 printf("Set Parameters\n");
112 fpt.SetParameters(7.24,0.120);
113 fpt.SetNpx(200);
114 }
32438f4e 115// Int_t nTracks = 400; //TODO: draw from experim dist
116// Int_t nTracks = 20; //TODO: draw from experim dist
1e62e876 117 Int_t nTracksLocal = fNtracks;
118 if(fHistosSet) {
119 nTracksLocal = fHMult->GetRandom();
120 if(nTracksLocal > fNtracks) nTracksLocal = fNtracks;
121 }
122 std::cout << "Generating " << nTracksLocal << " tracks " << std::endl;
223d9e38 123
1e62e876 124 for(Int_t iTrack=0; iTrack<nTracksLocal; iTrack++){
526ddf0e 125 Double_t phi = gRandom->Uniform(0.0, 2*TMath::Pi());
1e62e876 126 Double_t eta = 0.;
127 Double_t pt = 0.;
128
129 if(fHistosSet) {//draw from histograms if set
130 eta = fHEta->GetRandom();
131 pt = fHPt->GetRandom();
132 }
133 else {
134 eta = gRandom->Uniform(-etaCuts, etaCuts);
135 pt = fpt.GetRandom(); // momentum for f1
136 }
526ddf0e 137 Short_t sign=1;
138 if(gRandom->Rndm() < 0.5){
139 sign =1;
140 }else{
141 sign=-1;
142 }
143
144 Double_t theta = 2*TMath::ATan(TMath::Exp(-eta))-TMath::Pi()/2.;
145 Double_t pxyz[3];
146 pxyz[0]=pt*TMath::Cos(phi);
147 pxyz[1]=pt*TMath::Sin(phi);
148 pxyz[2]=pt*TMath::Tan(theta);
149 Double_t vertex[3]={0,0,retEvent->GetZ()};
150 Double_t cv[21]={0};
0403120d 151 AliToyMCTrack *tempTrack = retEvent->AddTrack(vertex,pxyz,cv,sign);
7d26e3ce 152 // use unique ID for track number
153 // this will be used in DistortTrack track to set the cluster label
223d9e38 154 // in one simulation the track id should be unique for performance studies
155 tempTrack->SetUniqueID(fCurrentTrack++);
0403120d 156 DistortTrack(*tempTrack, time);
526ddf0e 157 }
158
159 return retEvent;
160}
161
a1a695e5 162//________________________________________________________________
fa97f7bd 163void AliToyMCEventGeneratorSimple::RunSimulation(const Int_t nevents/*=10*/, const Int_t ntracks/*=400*/, const Int_t rate/*=50*/)
a1a695e5 164{
165 //
166 // run simple simulation with equal event spacing
fa97f7bd 167 // rate in kHz
a1a695e5 168 //
169
170 if (!ConnectOutputFile()) return;
cd8ed0ac 171 //initialise the space charge. Should be done after the tree was set up
172 InitSpaceCharge();
223d9e38 173
174 // number of tracks to simulate per interaction
32438f4e 175 fNtracks=ntracks;
223d9e38 176 // within one simulation the track count should be unique for effeciency studies
a06336b6 177 // don't use the track ID 0 since the label -0 is not different from 0
178 fCurrentTrack=1;
223d9e38 179
a1a695e5 180 Double_t eventTime=0.;
fa97f7bd 181 const Double_t eventSpacing=1./rate/1e3; //50kHz equally spaced
32438f4e 182 TStopwatch s;
a1a695e5 183 for (Int_t ievent=0; ievent<nevents; ++ievent){
32438f4e 184 printf("Generating event %3d (%.3g)\n",ievent,eventTime);
a1a695e5 185 fEvent = Generate(eventTime);
186 FillTree();
187 delete fEvent;
188 fEvent=0x0;
189 eventTime+=eventSpacing;
190 }
32438f4e 191 s.Stop();
192 s.Print();
193
a1a695e5 194 CloseOutputFile();
195}
526ddf0e 196
f356075c 197//________________________________________________________________
0403120d 198void AliToyMCEventGeneratorSimple::RunSimulationLaser(const Int_t nevents/*=1*/)
199{
200 //
201 // run simple simulation with equal event spacing
202 //
203
204 if (!ConnectOutputFile()) return;
205 //initialise the space charge. Should be done after the tree was set up
206 InitSpaceCharge();
207
208 // within one simulation the track count should be unique for effeciency studies
a06336b6 209 fCurrentTrack=1;
0403120d 210
211 Double_t eventTime=0.;
212 const Double_t eventSpacing=1./10.; //laser is running at 10Hz equally spaced
213 TStopwatch s;
214 for (Int_t ievent=0; ievent<nevents; ++ievent){
215 printf("Generating event %3d (%.3g)\n",ievent,eventTime);
216 fEvent = GenerateLaser(eventTime);
217 FillTree();
218 delete fEvent;
219 fEvent=0x0;
220 eventTime+=eventSpacing;
221 }
222 s.Stop();
223 s.Print();
224
225 CloseOutputFile();
226}
227
228//________________________________________________________________
f356075c 229AliToyMCEvent* AliToyMCEventGeneratorSimple::GenerateESD(AliESDEvent &esdEvent, Double_t time) {
230
231
232 AliToyMCEvent *retEvent = new AliToyMCEvent();
233 retEvent->SetT0(time);
234 retEvent->SetX(esdEvent.GetPrimaryVertex()->GetX());
235 retEvent->SetY(esdEvent.GetPrimaryVertex()->GetY());
236 retEvent->SetZ(esdEvent.GetPrimaryVertex()->GetZ());
237
238
239
240 if(!fNtracks) fNtracks = esdEvent.GetNumberOfTracks();
241 Int_t nUsedTracks = 0;
242 for(Int_t iTrack=0; iTrack<esdEvent.GetNumberOfTracks(); iTrack++){
243
244 AliESDtrack *part = esdEvent.GetTrack(iTrack);
245 if(!part) continue;
246 if (!fESDCuts->AcceptTrack(/*(AliESDtrack*)*/part))continue;
609d6d39 247 if(part->Pt() < 0.3) continue; //avoid tracks that fail to propagate through entire volume
f356075c 248 Double_t pxyz[3];
249 pxyz[0]=part->Px();
250 pxyz[1]=part->Py();
251 pxyz[2]=part->Pz();
252 Double_t vertex[3]={retEvent->GetX(),retEvent->GetY(),retEvent->GetZ()};
253
254 Double_t cv[21]={0};
255 Int_t sign = part->Charge();
0403120d 256 AliToyMCTrack *tempTrack = retEvent->AddTrack(vertex,pxyz,cv,sign);
7d26e3ce 257 // use unique ID for track number
258 // this will be used in DistortTrack track to set the cluster label
223d9e38 259 // in one simulation the track id should be unique for performance studies
260 tempTrack->SetUniqueID(fCurrentTrack++);
0403120d 261 DistortTrack(*tempTrack, time);
262 nUsedTracks++;
f356075c 263
264 if(nUsedTracks >= fNtracks) break;
265 }
266
267 return retEvent;
268}
269//________________________________________________________________
270void AliToyMCEventGeneratorSimple::RunSimulationESD(const Int_t nevents/*=10*/, const Int_t ntracks/*=400*/)
271{
272 //
273 // run simulation using esd input with equal event spacing
274 //
275
276 if (!ConnectOutputFile()) return;
277 TFile f(Form("%s%s",fInputFileNameESD.Data(),"#AliESDs.root"));
278 if(!f.IsOpen()) {
279 std::cout << "file "<<fInputFileNameESD.Data() <<" could not be opened" << std::endl;
280 return;
281 }
282 TTree* esdTree = (TTree*) f.Get("esdTree");
283 if(!esdTree) {
284 std::cout <<"no esdTree in file " << std::endl;
285 return;
286 }
287 InitSpaceCharge();
288 AliESDEvent* esdEvent = new AliESDEvent();
29b7f480 289
f356075c 290
291 esdEvent->ReadFromTree(esdTree);
29b7f480 292
223d9e38 293 fNtracks = ntracks;
294 // within one simulation the track count should be unique for effeciency studies
a06336b6 295 fCurrentTrack=1;
223d9e38 296
f356075c 297 Double_t eventTime=0.;
298 const Double_t eventSpacing=1./50e3; //50kHz equally spaced
299 TStopwatch s;
29b7f480 300
f356075c 301 fESDCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2011(kTRUE,1);
1e62e876 302 gRandom->SetSeed();
f356075c 303 Int_t nUsedEvents = 0;
304 for (Int_t ievent=0; ievent<esdTree->GetEntries(); ++ievent){
609d6d39 305 printf("Generating event %3d (%.3g)\n",nUsedEvents,eventTime);
f356075c 306 esdTree->GetEvent(ievent);
307 if(esdEvent->GetNumberOfTracks()==0) {
308 std::cout << " tracks == 0" << std::endl;
309 continue;
310 }
311
312 fEvent = GenerateESD(*esdEvent, eventTime);
609d6d39 313 if(fEvent->GetNumberOfTracks() >=10) {
314 nUsedEvents++;
315 FillTree();
316 eventTime+=eventSpacing;
317 }
f356075c 318 delete fEvent;
319 fEvent=0x0;
609d6d39 320
f356075c 321 if(nUsedEvents>=nevents) break;
322 }
323 s.Stop();
324 s.Print();
325
326 CloseOutputFile();
327}
328
609d6d39 329//________________________________________________________________
609d6d39 330void AliToyMCEventGeneratorSimple::RunSimulationBunchTrain(const Int_t nevents/*=10*/, const Int_t ntracks/*=400*/)
331{
332 //
333 // run simple simulation with equal event spacing
334 //
335
336 //Parameters for bunc
337 const Double_t abortGap = 3e-6; //
338 const Double_t collFreq = 50e3;
339 const Double_t bSpacing = 50e-9; //bunch spacing
340 const Int_t nTrainBunches = 48;
341 const Int_t nTrains = 12;
342 const Double_t revFreq = 1.11e4; //revolution frequency
343 const Double_t collProb = collFreq/(nTrainBunches*nTrains*revFreq);
344 const Double_t trainLength = bSpacing*(nTrainBunches-1);
345 const Double_t totTrainLength = nTrains*trainLength;
346 const Double_t trainSpacing = (1./revFreq - abortGap - totTrainLength)/(nTrains-1);
347 Bool_t equalSpacing = kFALSE;
348
349 TRandom3 *rand = new TRandom3();
350 //rand->SetSeed();
351
352 if (!ConnectOutputFile()) return;
353 //initialise the space charge. Should be done after the tree was set up
354 InitSpaceCharge();
355
356 fNtracks=ntracks;
223d9e38 357 // within one simulation the track count should be unique for effeciency studies
a06336b6 358 fCurrentTrack=1;
223d9e38 359
609d6d39 360 Double_t eventTime=0.;
29b7f480 361 // const Double_t eventSpacing=1./50e3; //50kHz equally spaced
609d6d39 362 TStopwatch s;
363 Int_t nGeneratedEvents = 0;
364 Int_t bunchCounter = 0;
365 Int_t trainCounter = 0;
366
367 //for (Int_t ievent=0; ievent<nevents; ++ievent){
368 while (nGeneratedEvents<nevents){
369 // std::cout <<trainCounter << " " << bunchCounter << " "<< "eventTime " << eventTime << std::endl;
370
371 if(equalSpacing) {
372 printf("Generating event %3d (%.3g)\n",nGeneratedEvents,eventTime);
373 fEvent = Generate(eventTime);
374 nGeneratedEvents++;
375 FillTree();
376 delete fEvent;
377 fEvent=0x0;
378 eventTime+=1./collFreq;
379 }
380 else{
381 Int_t nCollsInCrossing = rand -> Poisson(collProb);
382 for(Int_t iColl = 0; iColl<nCollsInCrossing; iColl++){
383 printf("Generating event %3d (%.3g)\n",nGeneratedEvents,eventTime);
384 fEvent = Generate(eventTime);
385 nGeneratedEvents++;
386 FillTree();
387 delete fEvent;
388 fEvent=0x0;
389
390 }
391 bunchCounter++;
392
393 if(bunchCounter>=nTrainBunches){
394
395 trainCounter++;
396 if(trainCounter>=nTrains){
397 eventTime+=abortGap;
398 trainCounter=0;
399 }
400 else eventTime+=trainSpacing;
401
402 bunchCounter=0;
403 }
404 else eventTime+= bSpacing;
405
406 }
407
408
409 }
410 s.Stop();
411 s.Print();
412 delete rand;
413 CloseOutputFile();
414}
29b7f480 415
416
417
418
419
420//________________________________________________________________
421Int_t AliToyMCEventGeneratorSimple::OpenInputAndGetMaxEvents(const Int_t type, const Int_t nevents) {
422
423
424
425 if(type==0) return nevents;
426
427 if(type==1) {
428
429 fInputFile = new TFile(Form("%s%s",fInputFileNameESD.Data(),"#AliESDs.root"),"read");
430 if(!fInputFile->IsOpen()) {
431 std::cout << "file "<<fInputFileNameESD.Data() <<" could not be opened" << std::endl;
432 return 0;
433 }
434 fESDTree = (TTree*) fInputFile->Get("esdTree");
435 if(!fESDTree) {
436 std::cout <<"no esdTree in file " << std::endl;
437 return 0;
438 }
439 fESDEvent = new AliESDEvent();
440 fESDEvent->ReadFromTree(fESDTree);
441 fESDCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2011(kTRUE,1);
442
443 fInputIndex = 0;
444
445 return fESDTree->GetEntries();
1e62e876 446 gRandom->SetSeed();
29b7f480 447 }
448
449
450 std::cout << " no available input type (toymc, esd) specified" << std::endl;
451 return 0;
452
453
454}
455
456
457//________________________________________________________________
458void AliToyMCEventGeneratorSimple::RunSimulation2(const Bool_t equalspacing, const Int_t type, const Int_t nevents, const Int_t ntracks) {
459
460 //type==0 simple toy
461 //type==1 esd input
462
463 Int_t nMaxEvents = OpenInputAndGetMaxEvents(type,nevents);
464 if (!ConnectOutputFile()) return;
465
223d9e38 466 InitSpaceCharge();
29b7f480 467
468
469 fNtracks = ntracks;
223d9e38 470 // within one simulation the track count should be unique for effeciency studies
a06336b6 471 fCurrentTrack=1;
223d9e38 472
29b7f480 473 Double_t eventTime=0.;
474 TStopwatch s;
475 // const Double_t eventSpacing=1./50e3; //50kHz equally spaced
476
477 Int_t nGeneratedEvents = 0;
478
479 for (Int_t ievent=0; ievent<nMaxEvents; ++ievent){
480 printf("Generating event %3d (%.3g)\n",ievent,eventTime);
481
482 Int_t nColls = 0;
483 Double_t spacing = 0;
484 GetNGeneratedEventsAndSpacing(equalspacing, nColls, spacing);
485
486 for(Int_t iColl = 0; iColl<nColls; iColl++){
487 if(type==0) fEvent = Generate(eventTime);
488 else if(type==1) fEvent = GenerateESD2(eventTime);
489 if(fEvent) {
490 FillTree();
491 delete fEvent;
492 fEvent=0x0;
493 nGeneratedEvents++;
494 }
495 }
496 eventTime+=spacing;
497 if(nGeneratedEvents >= nevents) break;
498 }
499 s.Stop();
500 s.Print();
501
502 CloseOutputFile();
503 CloseInputFile();
504
505
506}
507//________________________________________________________________
508AliToyMCEvent* AliToyMCEventGeneratorSimple::GenerateESD2(Double_t time) {
509
510 //test that enough tracks will pass cuts (and that there is tracks at all)
511 Bool_t testEvent = kTRUE;
512 while(fESDTree->GetEvent(fInputIndex) && testEvent) {
513
514 Int_t nPassedCuts = 0;
515 for(Int_t iTrack = 0; iTrack<fESDEvent->GetNumberOfTracks(); iTrack++){
516 if(fESDCuts->AcceptTrack(fESDEvent->GetTrack(iTrack)) && (fESDEvent->GetTrack(iTrack)->Pt() < 0.3) ) nPassedCuts++;
517 }
518 if(nPassedCuts<fNtracks || (fNtracks>100 && nPassedCuts <100) ) fInputIndex++; //100 is a perhaps bit arbitrary, do we want the events were only a small number of tracks are accepted?
519 else testEvent = kFALSE;
520 }
521
522 if(fInputIndex>=fESDTree->GetEntries()) return 0;
523
524
525
526 AliToyMCEvent *retEvent = new AliToyMCEvent();
527 retEvent->SetT0(time);
528 retEvent->SetX(fESDEvent->GetPrimaryVertex()->GetX());
529 retEvent->SetY(fESDEvent->GetPrimaryVertex()->GetY());
530 retEvent->SetZ(fESDEvent->GetPrimaryVertex()->GetZ());
531
532
533
534 if(!fNtracks) fNtracks = fESDEvent->GetNumberOfTracks();
535 Int_t nUsedTracks = 0;
536 for(Int_t iTrack=0; iTrack<fESDEvent->GetNumberOfTracks(); iTrack++){
537
538 AliESDtrack *part = fESDEvent->GetTrack(iTrack);
539 if(!part) continue;
540 if (!fESDCuts->AcceptTrack(/*(AliESDtrack*)*/part))continue;
541 if(part->Pt() < 0.3) continue; //avoid tracks that fail to propagate through entire volume
542
543 Double_t pxyz[3];
544 pxyz[0]=part->Px();
545 pxyz[1]=part->Py();
546 pxyz[2]=part->Pz();
547 Double_t vertex[3]={retEvent->GetX(),retEvent->GetY(),retEvent->GetZ()};
548 Double_t cv[21]={0};
549 Int_t sign = part->Charge();
550
0403120d 551 AliToyMCTrack *tempTrack = retEvent->AddTrack(vertex,pxyz,cv,sign);
29b7f480 552 // use unique ID for track number
553 // this will be used in DistortTrack track to set the cluster label
223d9e38 554 // in one simulation the track id should be unique for performance studies
555 tempTrack->SetUniqueID(fCurrentTrack++);
0403120d 556 DistortTrack(*tempTrack, time);
557
29b7f480 558
559 if(nUsedTracks >= fNtracks) break;
560 }
561 fInputIndex++;
562
563 return retEvent;
564}
565
566//________________________________________________________________
0403120d 567AliToyMCEvent* AliToyMCEventGeneratorSimple::GenerateLaser(Double_t time)
568{
569 //
570 // Generate an Event with laser tracks
571 //
572
573 AliToyMCEvent *retEvent = new AliToyMCEvent();
574 retEvent->SetEventType(AliToyMCEvent::kLaser);
575
576 retEvent->SetT0(time);
577 retEvent->SetX(0);
578 retEvent->SetX(0);
579 retEvent->SetZ(0);
580
581 AliTPCLaserTrack::LoadTracks();
582 TObjArray *arr=AliTPCLaserTrack::GetTracks();
583
584 //since we have a laser track force no material budges
585 Bool_t materialBudget=GetUseMaterialBudget();
586 SetUseMaterialBudget(kFALSE);
cedbc66e 587
588 //the laser tracks have partially large inclination angles over the pads
589 // -> relax the propagation contraint and switch off error messages
0403120d 590 SetIsLaser(kTRUE);
cedbc66e 591 Int_t debug=AliLog::GetDebugLevel("","AliToyMCEventGeneratorSimple");
592 AliLog::SetClassDebugLevel("AliToyMCEventGeneratorSimple",-3);
0403120d 593
594 for (Int_t iTrack=0; iTrack<arr->GetEntriesFast(); ++iTrack){
595 AliExternalTrackParam *track=(AliExternalTrackParam*)arr->At(iTrack);
596 AliToyMCTrack *tempTrack = retEvent->AddTrack(AliToyMCTrack(*track));
597 // for laser only TPC clusters exist
598 tempTrack->SetUniqueID(fCurrentTrack++);
599 MakeTPCClusters(*tempTrack, time);
600 }
601
602 SetIsLaser(kFALSE);
cedbc66e 603 AliLog::SetClassDebugLevel("AliToyMCEventGeneratorSimple",debug);
0403120d 604 SetUseMaterialBudget(materialBudget);
605 return retEvent;
606}
607
608//________________________________________________________________
29b7f480 609void AliToyMCEventGeneratorSimple::GetNGeneratedEventsAndSpacing(const Bool_t equalSpacing, Int_t &ngen, Double_t &spacing)
610{
611
612 static Int_t bunchCounter = 0;
613 static Int_t trainCounter = 0;
614
615 if(equalSpacing) {
616 ngen =1;
617 spacing = 1./50e3; //50kHz equally spaced
618 return;
619 }
620 else if(!equalSpacing) {
621 //parameters for bunch train
622 const Double_t abortGap = 3e-6; //
623 const Double_t collFreq = 50e3;
624 const Double_t bSpacing = 50e-9; //bunch spacing
625 const Int_t nTrainBunches = 48;
626 const Int_t nTrains = 12;
627 const Double_t revFreq = 1.11e4; //revolution frequency
628 const Double_t collProb = collFreq/(nTrainBunches*nTrains*revFreq);
629 const Double_t trainLength = bSpacing*(nTrainBunches-1);
630 const Double_t totTrainLength = nTrains*trainLength;
631 const Double_t trainSpacing = (1./revFreq - abortGap - totTrainLength)/(nTrains-1);
632 Double_t time = 0;
633 Int_t nCollsInCrossing = 0;
634 while(nCollsInCrossing ==0){
635
636
637 bunchCounter++;
638
639 if(bunchCounter>=nTrainBunches){
640
641 trainCounter++;
642 if(trainCounter>=nTrains){
643 time+=abortGap;
644 trainCounter=0;
645 }
646 else time+=trainSpacing;
647
648 bunchCounter=0;
649 }
650 else time+= bSpacing;
651
652
653 nCollsInCrossing = gRandom -> Poisson(collProb);
654 //std::cout << " nCollsInCrossing " << nCollsInCrossing << std::endl;
655 }
656 ngen = nCollsInCrossing;
657 if(nCollsInCrossing > 1)std::cout << " nCollsInCrossing " << nCollsInCrossing << std::endl;
658 spacing = time;
659 return;
660
661 }
662
663}
664
665//________________________________________________________________
666Bool_t AliToyMCEventGeneratorSimple::CloseInputFile()
667{
668 //
669 // close the input file
670 //
671 if (!fInputFile) return kFALSE;
672 fInputFile->Close();
673 delete fInputFile;
674 fInputFile=0x0;
675
676 return kTRUE;
677}