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