577c6af1bf1aaea30a6f107b10384cd953e5c494
[u/mrichter/AliRoot.git] / test / generators / therminator / Config.C
1 //
2 // Configuration for the first physics production 2008
3 //
4
5 // One can use the configuration macro in compiled mode by
6 // root [0] gSystem->Load("libgeant321");
7 // root [0] gSystem->SetIncludePath("-I$ROOTSYS/include -I$ALICE_ROOT/include\
8 //                   -I$ALICE_ROOT -I$ALICE/geant3/TGeant3");
9 // root [0] .x grun.C(1,"Config.C++")
10
11 #if !defined(__CINT__) || defined(__MAKECINT__)
12 #include <Riostream.h>
13 #include <TRandom.h>
14 #include <TDatime.h>
15 #include <TSystem.h>
16 #include <TVirtualMC.h>
17 #include <TGeant3TGeo.h>
18 #include "STEER/AliRunLoader.h"
19 #include "STEER/AliRun.h"
20 #include "STEER/AliConfig.h"
21 #include "PYTHIA6/AliDecayerPythia.h"
22 #include "PYTHIA6/AliGenPythia.h"
23 #include "TDPMjet/AliGenDPMjet.h"
24 #include "STEER/AliMagFCheb.h"
25 #include "STRUCT/AliBODY.h"
26 #include "STRUCT/AliMAG.h"
27 #include "STRUCT/AliABSOv3.h"
28 #include "STRUCT/AliDIPOv3.h"
29 #include "STRUCT/AliHALLv3.h"
30 #include "STRUCT/AliFRAMEv2.h"
31 #include "STRUCT/AliSHILv3.h"
32 #include "STRUCT/AliPIPEv3.h"
33 #include "ITS/AliITSv11.h"
34 #include "TPC/AliTPCv2.h"
35 #include "TOF/AliTOFv6T0.h"
36 #include "HMPID/AliHMPIDv3.h"
37 #include "ZDC/AliZDCv4.h"
38 #include "TRD/AliTRDv1.h"
39 #include "TRD/AliTRDgeometry.h"
40 #include "FMD/AliFMDv1.h"
41 #include "MUON/AliMUONv1.h"
42 #include "PHOS/AliPHOSv1.h"
43 #include "PHOS/AliPHOSSimParam.h"
44 #include "PMD/AliPMDv1.h"
45 #include "T0/AliT0v1.h"
46 #include "EMCAL/AliEMCALv2.h"
47 #include "ACORDE/AliACORDEv1.h"
48 #include "VZERO/AliVZEROv7.h"
49 #endif
50
51 enum PDC06Proc_t
52 {
53   kPythia6, kPhojet, kTherminator, kRunMax
54 };
55
56 const char  * pprRunName[] =
57 {
58   "kPythia6", "kPhojet", "kTherminator"
59 };
60
61 enum Mag_t
62 {
63   kNoField, k5kG, kFieldMax
64 };
65
66 const char * pprField[] =
67 {
68 "kNoField", "k5kG"
69 };
70
71 enum PprTrigConf_t
72 {
73 kDefaultPPTrig, kDefaultPbPbTrig
74 };
75
76 const char * pprTrigConfName[] =
77 {
78 "p-p","Pb-Pb"
79 };
80
81 static PprTrigConf_t strig = kDefaultPbPbTrig;// default PP trigger configuration
82
83 //--- Functions ---
84 class AliGenPythia;
85 AliGenerator *MbPythia();
86 AliGenerator *MbPhojet();
87 AliGenerator *Therminator();
88 void ProcessEnvironmentVars();
89
90 // Geterator, field, beam energy
91 static PDC06Proc_t proc = kTherminator;
92 static Mag_t mag = k5kG;
93 static Float_t energy = 10000; // energy in CMS
94 //======================= =//
95 // Set Random Number seed //
96 //========================//
97 TDatime dt;
98 static UInt_t seed = dt.Get();
99
100 // Comment line
101 static TString comment;
102
103 void Config()
104 {
105
106   // Get settings from environment variables
107   ProcessEnvironmentVars();
108
109   gRandom->SetSeed(seed);
110   cerr<<"Seed for random number generation= "<<seed<<endl;
111
112   // Libraries required by geant321
113 #if defined(__CINT__)
114   gSystem->Load("liblhapdf"); // Parton density functions
115   gSystem->Load("libEGPythia6"); // TGenerator interface
116   gSystem->Load("libpythia6"); // Pythia
117   gSystem->Load("libAliPythia6"); // ALICE specific implementations
118   gSystem->Load("libgeant321");
119 #endif
120
121   new TGeant3TGeo("C++ Interface to Geant3");
122
123   //=======================================================================
124   //  Create the output file
125
126
127   AliRunLoader* rl=0x0;
128
129   cout<<"Config.C: Creating Run Loader ..."<<endl;
130   rl = AliRunLoader::Open("galice.root",
131       AliConfig::GetDefaultEventFolderName(),
132       "recreate");
133   if (rl == 0x0)
134   {
135     gAlice->Fatal("Config.C","Can not instatiate the Run Loader");
136     return;
137   }
138   rl->SetCompressionLevel(2);
139   rl->SetNumberOfEventsPerFile(1000);
140   gAlice->SetRunLoader(rl);
141
142   // Set the trigger configuration: proton-proton
143   AliSimulation::Instance()->SetTriggerConfig(pprTrigConfName[strig]);
144   cout<<"Trigger configuration is set to  "<<pprTrigConfName[strig]<<endl;
145
146   //
147   //=======================================================================
148   // ************* STEERING parameters FOR ALICE SIMULATION **************
149   // --- Specify event type to be tracked through the ALICE setup
150   // --- All positions are in cm, angles in degrees, and P and E in GeV
151
152
153   gMC->SetProcess("DCAY",1);
154   gMC->SetProcess("PAIR",1);
155   gMC->SetProcess("COMP",1);
156   gMC->SetProcess("PHOT",1);
157   gMC->SetProcess("PFIS",0);
158   gMC->SetProcess("DRAY",0);
159   gMC->SetProcess("ANNI",1);
160   gMC->SetProcess("BREM",1);
161   gMC->SetProcess("MUNU",1);
162   gMC->SetProcess("CKOV",1);
163   gMC->SetProcess("HADR",1);
164   gMC->SetProcess("LOSS",2);
165   gMC->SetProcess("MULS",1);
166   gMC->SetProcess("RAYL",1);
167
168   Float_t cut = 1.e-3; // 1MeV cut by default
169   Float_t tofmax = 1.e10;
170
171   gMC->SetCut("CUTGAM", cut);
172   gMC->SetCut("CUTELE", cut);
173   gMC->SetCut("CUTNEU", cut);
174   gMC->SetCut("CUTHAD", cut);
175   gMC->SetCut("CUTMUO", cut);
176   gMC->SetCut("BCUTE", cut);
177   gMC->SetCut("BCUTM", cut);
178   gMC->SetCut("DCUTE", cut);
179   gMC->SetCut("DCUTM", cut);
180   gMC->SetCut("PPCUTM", cut);
181   gMC->SetCut("TOFMAX", tofmax);
182
183   //======================//
184   // Set External decayer //
185   //======================//
186   TVirtualMCDecayer* decayer = new AliDecayerPythia();
187   decayer->SetForceDecay(kAll);
188   decayer->Init();
189   gMC->SetExternalDecayer(decayer);
190
191   //=========================//
192   // Generator Configuration //
193   //=========================//
194   AliGenerator* gener = 0x0;
195
196   switch (proc)
197   {
198     case kPythia6:
199     gener = MbPythia();
200     break;
201     case kPhojet:
202     gener = MbPhojet();
203     break;
204     case kTherminator:
205     gener = Therminator();
206     break;
207   }
208
209   // PRIMARY VERTEX
210   //
211   gener->SetOrigin(0., 0., 0.); // vertex position
212   //
213   // Size of the interaction diamond
214   // Longitudinal
215   Float_t sigmaz = 5.4 / TMath::Sqrt(2.); // [cm]
216   if (energy == 900)
217   sigmaz = 10.5 / TMath::Sqrt(2.); // [cm]
218   //
219   // Transverse
220   Float_t betast = 10; // beta* [m]
221   Float_t eps = 3.75e-6; // emittance [m]
222   Float_t gamma = energy / 2.0 / 0.938272; // relativistic gamma [1]
223   Float_t sigmaxy = TMath::Sqrt(eps * betast / gamma) / TMath::Sqrt(2.) * 100.; // [cm]
224   printf("\n \n Diamond size x-y: %10.3e z: %10.3e\n \n", sigmaxy, sigmaz);
225
226   gener->SetSigma(sigmaxy, sigmaxy, sigmaz); // Sigma in (X,Y,Z) (cm) on IP position
227   gener->SetCutVertexZ(3.); // Truncate at 3 sigma
228   gener->SetVertexSmear(kPerEvent);
229
230   gener->Init();
231
232   // FIELD
233   AliMagF* field = 0x0;
234
235   if (mag == kNoField)
236   {
237     comment = comment.Append(" | L3 field 0.0 T");
238     field = new AliMagF("Maps","Maps", -1., -1., AliMagF::k2kG);
239   }
240   else if (mag == k5kG)
241   {
242     comment = comment.Append(" | L3 field 0.5 T");
243     field = new AliMagF("Maps","Maps", -1., -1., AliMagF::k5kG);
244   }
245   printf("\n \n Comment: %s \n \n", comment.Data());
246
247   TGeoGlobalMagField::Instance()->SetField(field);
248
249   rl->CdGAFile();
250
251   Int_t iABSO = 1;
252   Int_t iACORDE= 0;
253   Int_t iDIPO = 1;
254   Int_t iEMCAL = 0;
255   Int_t iFMD = 1;
256   Int_t iFRAME = 1;
257   Int_t iHALL = 1;
258   Int_t iITS = 1;
259   Int_t iMAG = 1;
260   Int_t iMUON = 1;
261   Int_t iPHOS = 1;
262   Int_t iPIPE = 1;
263   Int_t iPMD = 0;
264   Int_t iHMPID = 1;
265   Int_t iSHIL = 1;
266   Int_t iT0 = 1;
267   Int_t iTOF = 1;
268   Int_t iTPC = 1;
269   Int_t iTRD = 1;
270   Int_t iVZERO = 1;
271   Int_t iZDC = 1;
272
273   //=================== Alice BODY parameters =============================
274   AliBODY *BODY = new AliBODY("BODY", "Alice envelop");
275
276   if (iMAG)
277   {
278     //=================== MAG parameters ============================
279     // --- Start with Magnet since detector layouts may be depending ---
280     // --- on the selected Magnet dimensions ---
281     AliMAG *MAG = new AliMAG("MAG", "Magnet");
282   }
283
284   if (iABSO)
285   {
286     //=================== ABSO parameters ============================
287     AliABSO *ABSO = new AliABSOv3("ABSO", "Muon Absorber");
288   }
289
290   if (iDIPO)
291   {
292     //=================== DIPO parameters ============================
293
294     AliDIPO *DIPO = new AliDIPOv3("DIPO", "Dipole version 3");
295   }
296
297   if (iHALL)
298   {
299     //=================== HALL parameters ============================
300
301     AliHALL *HALL = new AliHALLv3("HALL", "Alice Hall");
302   }
303
304   if (iFRAME)
305   {
306     //=================== FRAME parameters ============================
307
308     AliFRAMEv2 *FRAME = new AliFRAMEv2("FRAME", "Space Frame");
309     FRAME->SetHoles(1);
310   }
311
312   if (iSHIL)
313   {
314     //=================== SHIL parameters ============================
315
316     AliSHIL *SHIL = new AliSHILv3("SHIL", "Shielding Version 3");
317   }
318
319   if (iPIPE)
320   {
321     //=================== PIPE parameters ============================
322
323     AliPIPE *PIPE = new AliPIPEv3("PIPE", "Beam Pipe");
324   }
325
326   if (iITS)
327   {
328     //=================== ITS parameters ============================
329
330     AliITS *ITS = new AliITSv11("ITS","ITS v11");
331   }
332
333   if (iTPC)
334   {
335     //============================ TPC parameters =====================
336
337     AliTPC *TPC = new AliTPCv2("TPC", "Default");
338   }
339
340   if (iTOF)
341   {
342     //=================== TOF parameters ============================
343
344     AliTOF *TOF = new AliTOFv6T0("TOF", "normal TOF");
345   }
346
347   if (iHMPID)
348   {
349     //=================== HMPID parameters ===========================
350
351     AliHMPID *HMPID = new AliHMPIDv3("HMPID", "normal HMPID");
352
353   }
354
355   if (iZDC)
356   {
357     //=================== ZDC parameters ============================
358
359     AliZDC *ZDC = new AliZDCv4("ZDC", "normal ZDC");
360   }
361
362   if (iTRD)
363   {
364     //=================== TRD parameters ============================
365
366     AliTRD *TRD = new AliTRDv1("TRD", "TRD slow simulator");
367     AliTRDgeometry *geoTRD = TRD->GetGeometry();
368     // Partial geometry: modules at 0,1,7,8,9,10,17
369     // starting at 3h in positive direction
370     geoTRD->SetSMstatus(2,0);
371     geoTRD->SetSMstatus(3,0);
372     geoTRD->SetSMstatus(4,0);
373     geoTRD->SetSMstatus(5,0);
374     geoTRD->SetSMstatus(6,0);
375     geoTRD->SetSMstatus(11,0);
376     geoTRD->SetSMstatus(12,0);
377     geoTRD->SetSMstatus(13,0);
378     geoTRD->SetSMstatus(14,0);
379     geoTRD->SetSMstatus(15,0);
380     geoTRD->SetSMstatus(16,0);
381   }
382
383   if (iFMD)
384   {
385     //=================== FMD parameters ============================
386
387     AliFMD *FMD = new AliFMDv1("FMD", "normal FMD");
388   }
389
390   if (iMUON)
391   {
392     //=================== MUON parameters ===========================
393     // New MUONv1 version (geometry defined via builders)
394
395     AliMUON *MUON = new AliMUONv1("MUON", "default");
396   }
397
398   if (iPHOS)
399   {
400     //=================== PHOS parameters ===========================
401
402     AliPHOS *PHOS = new AliPHOSv1("PHOS", "IHEP");
403   }
404
405   if (iPMD)
406   {
407     //=================== PMD parameters ============================
408
409     AliPMD *PMD = new AliPMDv1("PMD", "normal PMD");
410       }
411
412       if (iT0)
413 {
414   //=================== T0 parameters ============================
415 AliT0  *T0 = new AliT0v1("T0", "T0 Detector");
416 }
417
418 if (iEMCAL)
419 {
420   //=================== EMCAL parameters ============================
421
422   AliEMCAL *EMCAL = new AliEMCALv2("EMCAL", "EMCAL_COMPLETEV1");
423 }
424
425 if (iACORDE)
426 {
427   //  =================== ACORDE parameters ============================
428
429 AliACORDE *ACORDE = new AliACORDEv1("ACORDE", "normal ACORDE");
430 }
431
432 if (iVZERO)
433 {
434 //=================== ACORDE parameters ============================
435
436 AliVZERO *VZERO = new AliVZEROv7("VZERO", "normal VZERO");
437 }
438 }
439 //
440 //           PYTHIA
441 //
442
443 AliGenerator* MbPythia()
444 {
445 comment = comment.Append(" pp at 14 TeV: Pythia low-pt");
446 //
447 //    Pythia
448 AliGenPythia* pythia = new AliGenPythia(-1);
449 pythia->SetMomentumRange(0, 999999.);
450 pythia->SetThetaRange(0., 180.);
451 pythia->SetYRange(-12.,12.);
452 pythia->SetPtRange(0,1000.);
453 pythia->SetProcess(kPyMb);
454 pythia->SetEnergyCMS(energy);
455
456 return pythia;
457 }
458
459 AliGenerator* MbPhojet()
460 {
461 comment = comment.Append(" pp at 14 TeV: Phojet low-pt");
462 //
463 //    DPMJET
464 #if defined(__CINT__)
465 gSystem->Load("libDPMJET"); // Parton density functions
466 gSystem->Load("libTDPMjet"); // Parton density functions
467 #endif
468 AliGenDPMjet* dpmjet = new AliGenDPMjet(-1);
469 dpmjet->SetMomentumRange(0, 999999.);
470 dpmjet->SetThetaRange(0., 180.);
471 dpmjet->SetYRange(-12.,12.);
472 dpmjet->SetPtRange(0,1000.);
473 dpmjet->SetProcess(kDpmMb);
474 dpmjet->SetEnergyCMS(energy);
475
476 return dpmjet;
477 }
478
479 AliGenerator* Therminator()
480 {
481 comment = comment.Append(" pp at 14 TeV: Therminator");
482
483 #if defined(__CINT__)
484 gSystem->Load("libTTherminator"); // Therminator library
485 #endif
486 AliGenTherminator *genther = new AliGenTherminator();
487 genther->SetFileName("event.out");
488 genther->SetEventNumberInFile(1);
489 genther->SetTemperature(0.145);
490 genther->SetMiuI(-0.0009);
491 genther->SetMiuS(0.000);
492 genther->SetMiuB(0.0008);
493 genther->SetAlfaRange(8.0);
494 genther->SetRapRange(4.0);
495 genther->SetRhoMax(7.74);
496 genther->SetTau(9.74);
497 genther->SetModel("Lhyquid3D");
498 genther->SetLhyquidSet("LHC500C2030");
499
500 return genther;
501 }
502
503 void ProcessEnvironmentVars()
504 {
505 // Run type
506 if (gSystem->Getenv("CONFIG_RUN_TYPE"))
507 {
508   for (Int_t iRun = 0; iRun < kRunMax; iRun++)
509   {
510     if (strcmp(gSystem->Getenv("CONFIG_RUN_TYPE"), pprRunName[iRun])==0)
511     {
512       proc = (PDC06Proc_t)iRun;
513       cout<<"Run type set to "<<pprRunName[iRun]<<endl;
514     }
515   }
516 }
517
518 // Field
519 if (gSystem->Getenv("CONFIG_FIELD"))
520 {
521   for (Int_t iField = 0; iField < kFieldMax; iField++)
522   {
523     if (strcmp(gSystem->Getenv("CONFIG_FIELD"), pprField[iField])==0)
524     {
525       mag = (Mag_t)iField;
526       cout<<"Field set to "<<pprField[iField]<<endl;
527     }
528   }
529 }
530
531 // Energy
532 if (gSystem->Getenv("CONFIG_ENERGY"))
533 {
534   energy = atoi(gSystem->Getenv("CONFIG_ENERGY"));
535   cout<<"Energy set to "<<energy<<" GeV"<<endl;
536 }
537
538 // Random Number seed
539 if (gSystem->Getenv("CONFIG_SEED"))
540 {
541   seed = atoi(gSystem->Getenv("CONFIG_SEED"));
542 }
543 }