Extra header added to the list
[u/mrichter/AliRoot.git] / test / pileup / Config.C
1 // Pileup generation
2
3 enum PDC06Proc_t 
4 {
5   kPythia6, kPhojet, kRunMax
6 };
7
8 const char * pprRunName[] = {
9   "kPythia6", "kPhojet"
10 };
11
12 enum Mag_t
13 {
14   kNoField, k5kG, kFieldMax
15 };
16
17 const char * pprField[] = {
18   "kNoField", "k5kG"
19 };
20
21 //--- Functions ---
22 class AliGenPythia;
23 AliGenerator *MbPythia();
24 AliGenerator *MbPhojet();
25 void ProcessEnvironmentVars();
26
27 // Geterator, field, beam energy
28 static PDC06Proc_t   proc     = kPhojet;
29 static Mag_t         mag      = k5kG;
30 static Float_t       energy   = 10000; // energy in CMS
31 //========================//
32 // Set Random Number seed //
33 //========================//
34 TDatime dt;
35 static UInt_t seed    = dt.Get();
36
37 // Comment line
38 static TString comment;
39
40 void Config()
41 {
42     
43
44   // Get settings from environment variables
45   ProcessEnvironmentVars();
46
47   gRandom->SetSeed(seed);
48   cerr<<"Seed for random number generation= "<<seed<<endl; 
49
50   new TGeant3TGeo("C++ Interface to Geant3");
51
52   //=======================================================================
53   //  Create the output file
54
55    
56   AliRunLoader* rl=0x0;
57
58   cout<<"Config.C: Creating Run Loader ..."<<endl;
59   rl = AliRunLoader::Open("galice.root",
60                           AliConfig::GetDefaultEventFolderName(),
61                           "recreate");
62   if (rl == 0x0)
63     {
64       gAlice->Fatal("Config.C","Can not instatiate the Run Loader");
65       return;
66     }
67   rl->SetCompressionLevel(2);
68   rl->SetNumberOfEventsPerFile(1000);
69   gAlice->SetRunLoader(rl);
70   // gAlice->SetGeometryFromFile("geometry.root");
71   // gAlice->SetGeometryFromCDB();
72   
73   // Set the trigger configuration: proton-proton
74   AliSimulation::Instance()->SetTriggerConfig("p-p");
75
76   //
77   //=======================================================================
78   // ************* STEERING parameters FOR ALICE SIMULATION **************
79   // --- Specify event type to be tracked through the ALICE setup
80   // --- All positions are in cm, angles in degrees, and P and E in GeV
81
82   TVirtualMC * vmc = TVirtualMC::GetMC();
83
84     vmc->SetProcess("DCAY",1);
85     vmc->SetProcess("PAIR",1);
86     vmc->SetProcess("COMP",1);
87     vmc->SetProcess("PHOT",1);
88     vmc->SetProcess("PFIS",0);
89     vmc->SetProcess("DRAY",0);
90     vmc->SetProcess("ANNI",1);
91     vmc->SetProcess("BREM",1);
92     vmc->SetProcess("MUNU",1);
93     vmc->SetProcess("CKOV",1);
94     vmc->SetProcess("HADR",1);
95     vmc->SetProcess("LOSS",2);
96     vmc->SetProcess("MULS",1);
97     vmc->SetProcess("RAYL",1);
98
99     Float_t cut = 1.e-3;        // 1MeV cut by default
100     Float_t tofmax = 1.e10;
101
102     vmc->SetCut("CUTGAM", cut);
103     vmc->SetCut("CUTELE", cut);
104     vmc->SetCut("CUTNEU", cut);
105     vmc->SetCut("CUTHAD", cut);
106     vmc->SetCut("CUTMUO", cut);
107     vmc->SetCut("BCUTE",  cut); 
108     vmc->SetCut("BCUTM",  cut); 
109     vmc->SetCut("DCUTE",  cut); 
110     vmc->SetCut("DCUTM",  cut); 
111     vmc->SetCut("PPCUTM", cut);
112     vmc->SetCut("TOFMAX", tofmax); 
113
114
115
116
117   //======================//
118   // Set External decayer //
119   //======================//
120   TVirtualMCDecayer* decayer = new AliDecayerPythia();
121   decayer->SetForceDecay(kAll);
122   decayer->Init();
123   vmc->SetExternalDecayer(decayer);
124
125   //=========================//
126   // Generator Configuration //
127   //=========================//
128   // Create pileup generator
129   AliGenPileup *pileup = new AliGenPileup();
130
131   AliGenerator* gener = 0x0;
132   
133   if (proc == kPythia6) {
134       gener = MbPythia();
135   } else if (proc == kPhojet) {
136       gener = MbPhojet();
137   }
138   
139   // Set the pileup interaction generator
140   // The second argument is the pileup rate
141   // in terms of event rate per bunch crossing
142   pileup->SetGenerator(gener,0.01);
143   // Set the beam time structure
144   // Details on the syntax in STEER/AliTriggerBCMask
145   pileup->SetBCMask("72(1H1L)3420L");
146   // Examples of the pileup rate and beam structure settings
147   // Most of the information is taken from the LHC commissionning page
148   // rate from 0.01 (at 900GeV) to 0.76 (at 14TeV)
149   // 1 bunch/orbit     - bc-mask = "1H3563L"
150   // 43 bunches/orbit  - bc-mask = "43(1H80L)81L"
151   // 72 bunches/orbit  - bc-mask = "72(1H1L)3420L" (50ns mode)
152   // Please note that most of these setting should be cross-checked because
153   // for example the 43 bunches mode is taken at CMS IP and not the ALICE one.
154
155   // Generate the trigger interaction
156   pileup->GenerateTrigInteraction(kTRUE);
157
158   // PRIMARY VERTEX
159   //
160   pileup->SetOrigin(0., 0., 0.);    // vertex position
161   //
162   //
163   // Size of the interaction diamond
164   // Longitudinal
165   Float_t sigmaz  = 5.4 / TMath::Sqrt(2.); // [cm]
166   if (energy == 900)
167     sigmaz  = 10.5 / TMath::Sqrt(2.); // [cm]
168   //
169   // Transverse
170   Float_t betast  = 10;                 // beta* [m]
171   Float_t eps     = 3.75e-6;            // emittance [m]
172   Float_t gamma   = energy / 2.0 / 0.938272;  // relativistic gamma [1]
173   Float_t sigmaxy = TMath::Sqrt(eps * betast / gamma) / TMath::Sqrt(2.) * 100.;  // [cm]
174   printf("\n \n Diamond size x-y: %10.3e z: %10.3e\n \n", sigmaxy, sigmaz);
175     
176   pileup->SetSigma(sigmaxy, sigmaxy, sigmaz);      // Sigma in (X,Y,Z) (cm) on IP position
177   pileup->SetCutVertexZ(3.);        // Truncate at 3 sigma
178   pileup->SetVertexSmear(kPerEvent);
179
180   pileup->Init();
181
182   // FIELD
183   //
184   AliMagF* field = 0x0;
185   if (mag == kNoField) {
186     comment = comment.Append(" | L3 field 0.0 T");
187     field = new AliMagF("Maps","Maps", 0., 0., AliMagF::k5kGUniform);
188   } else if (mag == k5kG) {
189     comment = comment.Append(" | L3 field 0.5 T");
190     field = new AliMagF("Maps","Maps", -1., -1., AliMagF::k5kG);
191   }
192   printf("\n \n Comment: %s \n \n", comment.Data());
193   TGeoGlobalMagField::Instance()->SetField(field);
194     
195   rl->CdGAFile();
196
197   Int_t iABSO  = 1;
198   Int_t iACORDE= 0;
199   Int_t iDIPO  = 1;
200   Int_t iEMCAL = 0;
201   Int_t iFMD   = 1;
202   Int_t iFRAME = 1;
203   Int_t iHALL  = 1;
204   Int_t iITS   = 1;
205   Int_t iMAG   = 1;
206   Int_t iMUON  = 1;
207   Int_t iPHOS  = 1;
208   Int_t iPIPE  = 1;
209   Int_t iPMD   = 0;
210   Int_t iHMPID = 1;
211   Int_t iSHIL  = 1;
212   Int_t iT0    = 1;
213   Int_t iTOF   = 1;
214   Int_t iTPC   = 1;
215   Int_t iTRD   = 1;
216   Int_t iVZERO = 1;
217   Int_t iZDC   = 1;
218   
219
220     //=================== Alice BODY parameters =============================
221     AliBODY *BODY = new AliBODY("BODY", "Alice envelop");
222
223
224     if (iMAG)
225     {
226         //=================== MAG parameters ============================
227         // --- Start with Magnet since detector layouts may be depending ---
228         // --- on the selected Magnet dimensions ---
229         AliMAG *MAG = new AliMAG("MAG", "Magnet");
230     }
231
232
233     if (iABSO)
234     {
235         //=================== ABSO parameters ============================
236         AliABSO *ABSO = new AliABSOv3("ABSO", "Muon Absorber");
237     }
238
239     if (iDIPO)
240     {
241         //=================== DIPO parameters ============================
242
243         AliDIPO *DIPO = new AliDIPOv3("DIPO", "Dipole version 3");
244     }
245
246     if (iHALL)
247     {
248         //=================== HALL parameters ============================
249
250         AliHALL *HALL = new AliHALLv3("HALL", "Alice Hall");
251     }
252
253
254     if (iFRAME)
255     {
256         //=================== FRAME parameters ============================
257
258         AliFRAMEv2 *FRAME = new AliFRAMEv2("FRAME", "Space Frame");
259         FRAME->SetHoles(1);
260     }
261
262     if (iSHIL)
263     {
264         //=================== SHIL parameters ============================
265
266         AliSHIL *SHIL = new AliSHILv3("SHIL", "Shielding Version 3");
267     }
268
269
270     if (iPIPE)
271     {
272         //=================== PIPE parameters ============================
273
274         AliPIPE *PIPE = new AliPIPEv3("PIPE", "Beam Pipe");
275     }
276  
277     if (iITS)
278     {
279         //=================== ITS parameters ============================
280
281         AliITS *ITS  = new AliITSv11("ITS","ITS v11");
282     }
283
284     if (iTPC)
285     {
286       //============================ TPC parameters =====================
287
288         AliTPC *TPC = new AliTPCv2("TPC", "Default");
289     }
290
291
292     if (iTOF) {
293         //=================== TOF parameters ============================
294
295         AliTOF *TOF = new AliTOFv6T0("TOF", "normal TOF");
296     }
297
298
299     if (iHMPID)
300     {
301         //=================== HMPID parameters ===========================
302
303         AliHMPID *HMPID = new AliHMPIDv3("HMPID", "normal HMPID");
304
305     }
306
307
308     if (iZDC)
309     {
310         //=================== ZDC parameters ============================
311
312         AliZDC *ZDC = new AliZDCv4("ZDC", "normal ZDC");
313     }
314
315     if (iTRD)
316     {
317         //=================== TRD parameters ============================
318
319         AliTRD *TRD = new AliTRDv1("TRD", "TRD slow simulator");
320         AliTRDgeometry *geoTRD = TRD->GetGeometry();
321         // Partial geometry: modules at 0,1,7,8,9,10,17
322         // starting at 3h in positive direction
323         geoTRD->SetSMstatus(2,0);
324         geoTRD->SetSMstatus(3,0);
325         geoTRD->SetSMstatus(4,0);
326         geoTRD->SetSMstatus(5,0);
327         geoTRD->SetSMstatus(6,0);
328         geoTRD->SetSMstatus(11,0);
329         geoTRD->SetSMstatus(12,0);
330         geoTRD->SetSMstatus(13,0);
331         geoTRD->SetSMstatus(14,0);
332         geoTRD->SetSMstatus(15,0);
333         geoTRD->SetSMstatus(16,0);
334     }
335
336     if (iFMD)
337     {
338         //=================== FMD parameters ============================
339
340         AliFMD *FMD = new AliFMDv1("FMD", "normal FMD");
341    }
342
343     if (iMUON)
344     {
345         //=================== MUON parameters ===========================
346         // New MUONv1 version (geometry defined via builders)
347
348         AliMUON *MUON = new AliMUONv1("MUON", "default");
349     }
350
351     if (iPHOS)
352     {
353         //=================== PHOS parameters ===========================
354
355         AliPHOS *PHOS = new AliPHOSv1("PHOS", "noCPV_Modules123");
356     }
357
358
359     if (iPMD)
360     {
361         //=================== PMD parameters ============================
362
363         AliPMD *PMD = new AliPMDv1("PMD", "normal PMD");
364     }
365
366     if (iT0)
367     {
368         //=================== T0 parameters ============================
369         AliT0 *T0 = new AliT0v1("T0", "T0 Detector");
370     }
371
372     if (iEMCAL)
373     {
374         //=================== EMCAL parameters ============================
375
376         AliEMCAL *EMCAL = new AliEMCALv2("EMCAL", "EMCAL_FIRSTYEAR");
377     }
378
379      if (iACORDE)
380     {
381         //=================== ACORDE parameters ============================
382
383         AliACORDE *ACORDE = new AliACORDEv1("ACORDE", "normal ACORDE");
384     }
385
386      if (iVZERO)
387     {
388         //=================== ACORDE parameters ============================
389
390         AliVZERO *VZERO = new AliVZEROv7("VZERO", "normal VZERO");
391     }
392 }
393 //
394 //           PYTHIA
395 //
396
397 AliGenerator* MbPythia()
398 {
399       comment = comment.Append(" pp at 14 TeV: Pythia low-pt");
400 //
401 //    Pythia
402       AliGenPythia* pythia = new AliGenPythia(-1); 
403       pythia->SetMomentumRange(0, 999999.);
404       pythia->SetThetaRange(0., 180.);
405       pythia->SetYRange(-12.,12.);
406       pythia->SetPtRange(0,1000.);
407       pythia->SetProcess(kPyMb);
408       pythia->SetEnergyCMS(energy);
409       
410       return pythia;
411 }
412
413 AliGenerator* MbPhojet()
414 {
415       comment = comment.Append(" pp at 14 TeV: Pythia low-pt");
416 //
417 //    DPMJET
418       AliGenDPMjet* dpmjet = new AliGenDPMjet(-1); 
419       dpmjet->SetMomentumRange(0, 999999.);
420       dpmjet->SetThetaRange(0., 180.);
421       dpmjet->SetYRange(-12.,12.);
422       dpmjet->SetPtRange(0,1000.);
423       dpmjet->SetProcess(kDpmMb);
424       dpmjet->SetEnergyCMS(energy);
425
426       return dpmjet;
427 }
428
429 void ProcessEnvironmentVars()
430 {
431     // Run type
432     if (gSystem->Getenv("CONFIG_RUN_TYPE")) {
433       for (Int_t iRun = 0; iRun < kRunMax; iRun++) {
434         if (strcmp(gSystem->Getenv("CONFIG_RUN_TYPE"), pprRunName[iRun])==0) {
435           proc = (PDC06Proc_t)iRun;
436           cout<<"Run type set to "<<pprRunName[iRun]<<endl;
437         }
438       }
439     }
440
441     // Field
442     if (gSystem->Getenv("CONFIG_FIELD")) {
443       for (Int_t iField = 0; iField < kFieldMax; iField++) {
444         if (strcmp(gSystem->Getenv("CONFIG_FIELD"), pprField[iField])==0) {
445           mag = (Mag_t)iField;
446           cout<<"Field set to "<<pprField[iField]<<endl;
447         }
448       }
449     }
450
451     // Energy
452     if (gSystem->Getenv("CONFIG_ENERGY")) {
453       energy = atoi(gSystem->Getenv("CONFIG_ENERGY"));
454       cout<<"Energy set to "<<energy<<" GeV"<<endl;
455     }
456
457     // Random Number seed
458     if (gSystem->Getenv("CONFIG_SEED")) {
459       seed = atoi(gSystem->Getenv("CONFIG_SEED"));
460     }
461 }