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