Replace the PHOS version from IHEP to Run1 and fix the compilation, when possibletest...
[u/mrichter/AliRoot.git] / test / generators / herwig / Config.C
CommitLineData
97bd8e56 1// One can use the configuration macro in compiled mode by
2// root [0] gSystem->Load("libgeant321");
5889a25e 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()
97bd8e56 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>
5889a25e 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"
97bd8e56 52#endif
53
54enum PDC06Proc_t
55{
56 kPythia6, kPhojet, kHerwig, kRunMax
57};
58
59const char * pprRunName[] =
60{ "kPythia6", "kPhojet", "kHerwig" };
61
62enum Mag_t
63{
64 kNoField, k5kG, kFieldMax
65};
66
67const char * pprField[] =
68{ "kNoField", "k5kG" };
69
70enum PprTrigConf_t
71{
72 kDefaultPPTrig, kDefaultPbPbTrig
73};
74
75const char * pprTrigConfName[] =
76{ "p-p", "Pb-Pb" };
77
78static PprTrigConf_t strig = kDefaultPPTrig;// default PP trigger configuration
79
80//--- Functions ---
81
82class AliGenPythia;
83AliGenerator *MbPythia();
84AliGenerator *MbPhojet();
85AliGenerator *Herwig();
86void ProcessEnvironmentVars();
87
88// Geterator, field, beam energy
89static PDC06Proc_t proc = kHerwig;
90static Mag_t mag = k5kG;
91static Float_t energy = 10000; // energy in CMS
92//========================//
93// Set Random Number seed //
94//========================//
95TDatime dt;
96static UInt_t seed = dt.Get();
97
98// Comment line
99static TString comment;
100
101void 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
93133609 139 AliSimulation::Instance()->SetTriggerConfig(pprTrigConfName[strig]);
97bd8e56 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:
5889a25e 201 default:
97bd8e56 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
07787d1d 327 AliITS *ITS = new AliITSv11("ITS", "ITS v11");
97bd8e56 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
047922b1 356 AliZDC *ZDC = new AliZDCv4("ZDC", "normal ZDC");
97bd8e56 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
5889a25e 399 AliPHOS *PHOS = new AliPHOSv1("PHOS", "Run1");
97bd8e56 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
c76b9217 419 AliEMCAL *EMCAL = new AliEMCALv2("EMCAL", "EMCAL_COMPLETEV1");
97bd8e56 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
438AliGenerator* 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
454AliGenerator* MbPhojet()
455{
456 comment = comment.Append(" pp at 14 TeV: Phojet low-pt");
457
458#if defined(__CINT__)
257e1b2b 459 gSystem->Load("libDPMJET"); // Parton density functions
97bd8e56 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
475AliGenerator* Herwig()
476{
477 comment = comment.Append("pp at 14 TeV: Herwig");
478
479#if defined(__CINT__)
257e1b2b 480 gSystem->Load("libHERWIG"); // HERWIG library
97bd8e56 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
500void 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}