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